Содержание к диссертации
Введение
ГЛАВА 1. О численных методах, применяемых для моделирования сложных многокомпонентных биологических систем 37
1.1. Об устойчивости нелинейных разностных схем 38
1.2. О свойствах разностных схем для решения модельного квазилинейного уравнения теплопроводности 42
1.3. О численном решении систем уравнений типа реакция- диффузия 52
Заключение к Главе 1 61
ГЛАВА 2. Популяционные модели и квазилинейные уравнения параболического типа 62
2.1 Модели одной популяции. Автомодельные решения 63
2.2 Модель конкурирующих популяций 76
2.3 Модель формирования зон роста у растений 84
Заключение к Главе 2 90
ГЛАВА 3. Исследование феноменологической математической модели свертывания крови in vitro 91
3.1 Двухавтоволновая модель динамики свертывания крови 94
3.2 Результаты численных исследований в ID случае 99
3.3 Двумерные структуры в модели свертывания крови 110
3.4 Резонансные явления в системе типа реакция - диффузия 120
3.5 Качественное исследование начального этапа формирования структур в модели свертывания крови 132
3.6 Моделирование роста оторвавшегося тромба в пристеночном потоке 142
Заключение к Главе 3 154
ГЛАВА 4. Исследование математической модели свертывания крови с учетом гипотезы о переключении активности тромбина 155
4.1 Модель динамики свертывания крови с учетом гипотезы о переключении активности тромбина 155
4.2. О существовании пространственно-неоднородных стационарных решений в модели динамики свертывания крови с учетом гипотезы о переключении активности тромбина 163
4.3. Об устойчивости пространственно-неоднородных стационарных решений в модели динамики свертывания крови с учетом гипотезы о переключении активности тромбина 167
Заключение к Главе 4 175
ГЛАВА 4. Математические модели химических реакций на клеточной мембране с учетом электрических полей 176
4.1 Математическая модель пространственно-временных структур в системе «реакция-диффузия» при приложенном внешнем электрическом поле 178
4.2 Параметрический резонанс и формирование диссипативных структур в двухкомпонентной системе при воздействии периодического электрического поля 193
4.3 О границах применимости моделей 205
Заключение к Главе 5 208
Заключение. Основные выводы 210
Литература 213
- О свойствах разностных схем для решения модельного квазилинейного уравнения теплопроводности
- Модель формирования зон роста у растений
- Качественное исследование начального этапа формирования структур в модели свертывания крови
- Об устойчивости пространственно-неоднородных стационарных решений в модели динамики свертывания крови с учетом гипотезы о переключении активности тромбина
О свойствах разностных схем для решения модельного квазилинейного уравнения теплопроводности
Математические модели роста стеноза предложены в [162-164] и активно используются в расчетах. Так, модель стеноза [163] была использована в [165] для расчета осесимметричных стационарных течений идеальной несжимаемой жидкости в области локального сужения сосуда. В работах [166, 167] модели Янга развития стеноза были соединены с моделями движения в сосудах жидкости со сложными реологическими характеристиками: жидкости Кэссона [168] и двухвяз-костой моделью [169], представляющей собой вариант модели жидкости Шведова — Бингама [170]. Также гемодинамические модели со сложными реологическими соотношениями рассматривались в [169, 170]. Все перечисленные работы про описании течений крови в сосудах во-первых предполагают, что реологические свойства крови постоянны в течение всего процесса (что, вообще говоря, неверно [164]). Во-вторых, гемодинамические факторы во всех работах рассматриваются в отрыве от описания процессов свертывания. В свою очередь, при рассмотрении моделей свертывания крови гидродинамические потоки, играющие ключевую роль, в расчет не принимаются.
В работах [171, 172] с участием автора диссертации строятся минимальные (по числу учитываемых метаболитов свертывания и химических реакций) модели, учитывающие как гемодинамические, так й коагулогические свойства крови. При этом работа [171] посвящена исследованию влияния сдвиговых течений на свойства свертываемости крови, а в [172] рассмотрена система с учетом простейшей модели тромболизиса. Течения в сосуде (который моделируется плоским каналом) учитывают как изменение формы за счет развития стеноза (модели Янга [163]), так и изменения формы расчетной области за счет быстрого формирования тромба. При формулировке данных математических моделей и проведении численных расчетов существенную роль оказало использование результатов, вошедших в Главу 2 диссертационной работы.
При численном моделировании задач популяционной динамики, формирования структур в модели свертывания крови in vitro и с учетом гемо динамических течений использовались конечно-разностные методы [173, 174, 175, 176, 177]. Во многих расчетах зарубежных авторов для гемодинамических задач используются различные варианты метода конечных элементов [166, 167, 169]. Рассмотрение их не входит в число основных целей данной работы, их приложения к рассматриваемым задачам не обсуждаются.
Затронутые в данной работе темы продолжают оставаться в центре внимания. Так, различные биологические системы продолжают давать обширный материал для создания новых нелинейных математических моделей, а математические модели служат основой для создания новых математических методов. Обзор современного состояния математического аппарата для изучения нелинейных систем дан, например, в [7, 178].
Отметим, что модели сложных нелинейных явлений не обязательно должны включать в себя большое количество переменных и управляющих параметров. В соответствии с подходами, развивающимися в рамках современной нелинейной динамики и синергетики [7, 50-52, 179], в системе можно выделить сравнительно небольшое число переменных и параметров (либо входящих в систему, либо агрегированных), называемых параметрами порядка. При этом параметры порядка определяют качественно поведение сложной системы и ее эволюцию на больших временах. Существует также ряд «достаточно простых» математических моделей, называемых иногда «базовыми математическими моделями» [12, 37]. Набор базовых систем описывает набор эффектов, наблюдаемых, например, в реакционно-диффузионных системах. Близкие к традиционным подходам синергетики и прикладной нелинейной динамики использовались при изучении режимов с обострением в задачах для квазилинейных уравнений параболического типа [6, 109-111]. Роль базовых моделей при этом играют автомодельные решения (или приближенные автомодельные решения) и собственные функции горения нелинейной среды [109, ПО]. О роли параметров порядка в квазилинейных задачах подробно рассказано в [6, 179].
Подходы, развиваемые автором в данной диссертации, по своим идеям близки как к традиционным подходам синергетики, так и подходам, применяемым при исследовании режимов с обострением. Первые две главы посвящены рассмотрению довольно простых моделей (по уровню сложности соизмеримых с «базовыми»). Такие системы уравнений позволяют качественно объяснить некоторые явления в моделях биологических систем, описываемых полулинейными или квазилинейными уравнениями параболического типа. Для исследования свойств таких математических моделей применяются как качественные, так и численные методы. Вторая глава посвящена рассмотрению двух простейших моделей свертывания крови. При этом даже если биохимическое значение некоторых переменных и констант моделей остается не вполне ясным, их можно трактовать как параметры порядка системы. Тем более, что результаты численных экспериментов находятся в качественном соответствии с лабораторными данными.
В четвертой главе, напротив, при исследовании явлений в окрестности мембраны использована одна из базовых моделей — брюсселятор.
Основные публикации. По теме диссертации автором опубликованы 22 статьи и подготовлено 16 докладов на конференциях разного уровня. Это статьи [68-71,98,99, 125, 171, 172, 184, 185,202,203,211 218,232,234,236,237,248, 252, 253]. Тезисы докладов включены в список литературы в конце работы под номерами [204, 235, 257-270].
Основные результаты работы также докладывались на научных семинарах в ИПМ РАН им. М. В. Келдыша (руководители семинара — д.ф.-м.н. Г. Г. Малинецкий и чл.-корр. РАН С. П. Курдюмов), ВЦ РАН (рабочая группа по математическому моделированию в экологии, руководитель — д.ф.-м.н. Ю. М. Свирежев), лаборатории ИФА РАН (руководитель семинара — д.ф.-м.н. Д. О. Логофет), семинаре «Нелинейный мир» на биологическом факультете МГУ (руководители семинара — д.ф.-м.н. Г. Ю. Ризниченко и чл.-корр. РАН А. Б. Рубин), семинаре «Нелинейная физика» (руководитель семинара — академик РАН Е. П. Велихов) и других.
Автор считает своим долгом выразить благодарность всем, с кем посчастливилось сотрудничать при постановке и решении задач, включенных в данную диссертационную работу. Это д. б. н., профессор Ф. И. Атауллоханов, к. ф. -м. н. Н. В. Белотелов, к. ф. -м. н. Е. В. Гельфанд, к. ф. -м. н. А. П. Гузеватых, к. ф. -м. н. Г. Т. Гурия, к. ф. -м. н. О. В. Демин, к. ф. -м. н. В. И. Зарницына, А. И. Лаврова, к. ф. -м. н. О. Л. Морозова, к. ф. -м. н., доц. Т. Ю. Плюснина, д. ф. -м. н., профессор Г. Ю. Ризниченко, чл.-корр. РАН А. Б. Рубин, к. ф. -м. н. Т. К. Старожилова, д. ф. -м. н., профессор А. П. Черняев.
Автор благодарен Российскому фонду фундаментальных исследований за поддержку работ, гранты №96-01-01306,99-01-01145. Значительная часть результатов, полученных при проведении работ по этим проектам, включена в главы 2 и 3 настоящей работы.
Модель формирования зон роста у растений
В первой главе сформулирована и доказана теорема об устойчивости разностных схем по линейному приближению. Отметим, что полученный результат носит несколько ограниченный характер, так как использует информацию о свойствах точного решения дифференциальной задачи, которое,--как правило, не извест-но априори. Полностью исследовать такую схему на устойчивость по линейному приближению возможно лишь для некоторых тестовых задач, примеры которых рассматриваются в тексте параграфа. Более того, если для дифференциальной задачи не установлено существование решения или такое решение не единственно, то приведенная схема исследования по линейному приближению становится неприемлемой.
Для систем уравнений типа «реакция-диффузия» достаточно лишь предположения о существовании и единственности решения (установить которое достаточно трудно, но в случае моделирования физических (биологических) процессов существование и единственность поставленной задачи обычно подразумевается) и информации о некоторых его свойствах. Так, для уравнений типа «реакция-диффузия» концентрации реагентов (или численность популяций для популяци-онньгх задач) не могут принимать отрицательные значения, что существенно облегчает исследование задач на устойчивость. В таких системах вначале исследуется дифференцируемость кинетической части системы как функции многих переменных в фазовом пространстве задачи, а затем исследуется устойчивость по линейному приближению.
Рассмотрено также приложение метода неопределенных коэффициентов к исследованию свойств одномерных уравнений типа «реакция-диффузия». На основе проведенных исследований в одномерном случае осуществляется выбор конкретных разностных схем метода расщепления по физическим процессам для решения задач, рассматриваемых в следующих главах диссертации.
Рассмотрение математических моделей биологических систем начнем классических математических моделей. Ка к уже было отмечено во введении, математические модели популяционного типа исторически были первыми математическими моделями биологических процессов. Хотя простейшие популяционные модели не описывают, по-видимому, никаких реальных процессов, они тем не менее позволяют прояснить и понять многие качественные механизмы природных явлений. Так как кинетическая часть уравнений Лотки - Вольтерра была фактически заимствована из химической кинетики (она тождественна описанию химических реакций второго порядка), то провести четкую грань между популяционными моделями и другими моделями типа «реакция-диффузия» затруднительно. К тому же часто популяционные модели включаются в качестве составной части в модели более сложных систем биологической природы. В качестве примера отметим модель агрегации тромбоцитов [78], где система Лотки - Вольтерра используется для описания взаимодействия между агрегатами.
В данной главе рассмотрены модели с достаточно простыми источниковыми (реакционными) частями, характерными для задач популяционной экологии. Новым элементом является применение квазилинейных уравнений для описания пространственно-временной динамики как изолированных, так и взаимодействующих популяций. Хотя модели такого рода встречаются в современной литературе [32], они, с одной стороны, не получили широкого распространения ввиду недостаточной обоснованности уравнений, с другой стороны, они пока недостаточно изучены. В цитированной выше работе Марри [32] применяются модификации методов фазовой плоскости, в то время как более естественным представляется использование других методов, в частности, построение автомодельных решений системы других типов, не ограничиваясь рассмотрением автомодельно-сти типа бегущей волны.
Применительно к моделированию задач химической кинетики, системы уравнений похожего вида рассмотрены в [113]. В заключительном параграфе данной главы рассматривается модель формирования зон роста у растений. Так как в системе рассмотрены реакции наработки ауксина и его взаимодействия с ингибитором (второго порядка) и модель транспорта ауксина с участием переносчиков, то с математической точки зрения она не отличается от рассматриваемых популяционных моделей.
Проблема математического описания пространственно-временной динамики одной или нескольких взаимодействующих популяций имеет длительную историю. Центральное место в ней занимает задача анализа взаимодействия миграционных популяционных процессов с демографическими. Под демографическими понимаются процессы рождения и смерти. Популяционные модели, учитывающие только их, широко известны и достаточно хорошо развиты [193, 194]. Это так называемые точечные модели, основным допущением при их построении является предположение о бесконечно быстром перемешивании особей на рассматриваемом ареале. Оно справедливо, если рассматриваемый ареал достаточно мал по сравнению с радиусом индивидуальной активности особей [4]. При нарушении этого условия в популяционных моделях необходимо учитывать миграцию. Самым простым и широко используемым в настоящее время положением является гипотеза о случайных блужданиях особей по пространству. Это предположение позволяет обосновать использование в качестве инструмента моделирования уравнения типа «реакция - диффузия», где в качестве реакционной части используются правые части точечных моделей, а коэффициенты диффузии (подвижности особей) считаются постоянными [1, 4, 30-32].
В рамках этих моделей удается объяснить волны распространения числен-ности популяции при заселении ареала, а также существование сложной пространственно- временной динамики численности. Тем не менее, центральная гипотеза таких моделей о случайном характере перемещения особей по ареалу обитания не выдерживает критических замечаний биологов.
Наблюдения за поведением хвоели сто грызущих насекомых [195] показывают, что миграции являются неотъемлемой чертой экологии этих видов. Эти миграции, в свою очередь, имеет разную интенсивность и направленность в зависимости от локальной плотности популяции. Так, специфика пространственного перемещения шелкопряда в лесных массивах позволяет выделить ряд функциональных типов миграции, отличающихся интенсивностью и направленностью потоков. Аналогичные явления наблюдаются в популяциях мышей и леммингов [196]. Для учета подобных феноменологических наблюдений за пространственно - временной динамикой численности популяций необходимо введение миграционных членов, зависящих от численности. Одним из способов учета этого явления — введение в уравнение типа реакция - диффузия градиентных членов, учитывающих различные таксисы. Существует ряд работ, в которых проводится исследование такого рода моделей [32, 43—45, 197].
Качественное исследование начального этапа формирования структур в модели свертывания крови
Малые по амплитуде или/и ширине возмущения стремятся к нулевому стационарному состоянию и не способны инициировать автокаталитические процессы. Существование порога по амплитуде возмущения следует из (3.4). Наличие порога по пространственному размеру начального возмущения объясняется влиянием диффузионных потоков на начальном этапе эволюции. Они уменьшают амплитуду активатора до значения ниже порогового, стремясь сгладить неоднородности его распределения. Подробнее проблема роли параметров начального возмущения на эволюцию системы рассматривается ниже.
В случае достаточно большого начального возмущения в системе наряду с диффузией происходит автокаталитическое образование активатора. Наличие диффузии и автокаталитического роста приводит к распространению автоволнового фронта активатора. Его рост останавливает ингибитор, образование которого происходит с некоторой задержкой по времени (лагом). Это обусловлено различием констант скоростей реакций а и (3 на три порядка. Как следует из (3.2), в области, где концентрация активатора велика, возможно распространение автоволны ингибитора.
В результате взаимодействия реагентов первоначально локализованная концентрационная структура активатора распадается на две страты, симметричные относительно центра деления (точки х = 0), см. рис. 3.1а. Концентрация реагентов постепенно уменьшается, а резкие фронты сглаживаются за счет диффузии. Дальнейшая эволюция страт зависит от параметров модели. Имеется два различных пути эволюции системы.
В первом амплитуда страт становится подпороговой, и система приходит в стационарное состояние (9 = О, ф = 0). Вследствие этого режима формируется локализованная полимерная (фибриновая) структура из одного структурного элемента (полосы или кольца, в зависимости от рассматриваемой задачи).
Во втором случае количество активатора в стратах достаточно для того, чтобы в них начался самоускоряющийся процесс его производства и последующего взаимодействия с ингибитором. Каждая страта делится на две новые (см. рис. 3.16). При этом они могут уже не быть симметричными относительно центра деления из-за наличия ингибитора, наработанного на предыдущих этапах развития. В отличие от первого деления, далее возможны три различных типа поведения системы. a) Все страты затем приходят к нулевому стационарному состоянию. Формируется периодическая полимерная структура, состоящая из трех структурных элементов. Следует отметить, что финальное число структурных элементов фибрина может быть больше трех. Например, когда решение представляет собой комбинацию нескольких этапов Ь) и этапа а). Важно, что сформировавшаяся полимерная структура занимает не всю рассматриваемую область вплоть до правой границы, а только часть ее. b) К нулевому стационарному состоянию на данном этапе приходят страты, расположенные ближе к центру возбуждения системы см. рис. \в Оставшиеся вновь участвуют в процессах роста и деления с тремя представленными здесь альтернативами дальнейшего развития. с) Амплитуда концентрации активатора в стратах, расположенных ближе к центру возбуждения, остается надпороговой. В этом случае имеет место так называемый эффект кинетического эха, см. рис. 3.1г. Рост и распад внутренних страт, как правило, происходит в том же месте, где развивались их предшественники. В результате соответствующий элемент полимерной структуры становится более контрастным по сравнению с остальными. При этом страты, расположенные дальше от центра возбуждения, могут, как остаться надпороговыми и внести вклад в дальнейшую эволюцию системы, так и исчезнуть. Если в системе наблюдается более одной автокаталитической вспышки активатора, то ее эволюция может быть описана в виде комбинации нескольких этапов а)-с). При этом из-за наличия ингибитора от предыдущих актов деления, время обострения для разных страт активатора может различаться. Как следствие, различается и амплитуда страт. Процесс их развития в случае многократного повторения сценариев Ь)-с) будет выглядеть как автоволновой (типа бегущего импульса или бегущей волны, генерирующей эховолны). При этом будут наблюдаться пульсации амплитуды и связанные с ними пульсации скорости волнового фронта. Все режимы структурооб-разования можно расклассифицировать следующим образом.
Режим с образованием одной локализованной полимерной структуры. Основные стадии этого режима были описаны выше. Характерное время образования фибриновой полосы составляет 4-10 мин и зависит от параметров модели.
Режим бегущей пульсирующей волны. Этот режим представляет собой многократное повторение сценария Ь) до тех пор, пока волна не достигнет границ рассматриваемой области. В результате формируется периодическая полимерная структура в виде последовательности фибриновых полос или системы концентрических колец, занимающая всю расчетную область. Режим с образованием эховолны. В этом случае кроме волны активатора, распространяющейся от места возбуждения среды к границам рассматриваемой области, образуется волна, распространяющаяся в обратную сторону. Она, в свою очередь, может породить вторую волну, распространяющуюся вперед вслед за первой. Более того, одна волна может порождать эховолны несколько раз. Формируется сложная периодическая полимерная структура, занимающая всю рассматриваемую область. Ее элементы могут иметь различную амплитуду, а в некоторых случаях и ширину.
Режим с образованием полимерной структуры из конечного числа элементов. Такая структура получается в результате комбинации сценариев Ь) или Ь) и с), заканчивающейся сценарием а, когда волны активатора еще не достигли границ рассматриваемой области. В ID плоском случае получается структура из нескольких фибриновых полос, в цилиндрически-симметричном — структура типа мишени. Интересно, что режимы, подобные описанным, наблюдались в экспериментах по изучению свертывания крови in vitro [64-66,217].
Об устойчивости пространственно-неоднородных стационарных решений в модели динамики свертывания крови с учетом гипотезы о переключении активности тромбина
В данной главе рассмотрены модели, описывающие протекание химических реакций в пространственно-распределенных системах. При этом на конечный ответ существенное влияние может оказать не только учет стохастического переноса заряженных частиц, но и их направленное движение в электрическом поле (электрический ток). Модели с учетом электрических эффектов необходимо строить как для описания движения крови (эритроциты и тромбоциты несут на себе электрический заряд, а внутренняя поверхность сосуда, как правило, также имеет поверхностно-распределенный заряд), так и для описания других процессов на уровне клеточных мембран и ансамблей клеток.
С другой стороны, такие исследования актуальны при проведении математического моделирования внешнего слабого электровоздействия на объекты, когда напряженность приложенного электрического поля достаточно мала и неспособна существенно изменить температуру объекта (воздействие не является тепловым).
В настоящей главе рассмотрение процессов ведется, в основном, на модельном уровне. Связано это с тем, что, несмотря на накопленный экспериментальный материал, в настоящее время не построены реалистические модели трансмембранных токов, пригодные для использования в качестве составной части более сложных математических моделей.
Неравномерное распределение потенциала по поверхности мембран и связанные с ним токи вдоль мембран — распространенное явление среди биологических объектов. В качестве примера укажем на корневую систему растений или клетки зиготы морской водоросли Fucus. Биологической моделью для изучения более сложных систем такого типа служит одноклеточная нитчатая водоросль Chara. При описании даже самых простых объектов, на мембране которых имеется неравномерное распределение потенциалов, используются либо аппроксимации экспериментальных зависимостей, либо функции, физический смысл которых не до конца ясен. Часто функциональные зависимости не расшифровываются. Как правило, при рассмотрении электрических эффектов в моделях реакционно-диффузионного типа (устоялось также другое название — модели типа реакция - электродиффузия) условие нейтральности среды приобретает вид условия локальной электронейтральности, т. е. алгебраическая сумма зарядов в любой точке сплошной среды должна быть равна нулю. Условия такого типа справедливы, если дебаевский радиус мал, например, для растворов с высокой ионной силой. В противном случае, при соблюдении условий нейтральности «в целом» локальные нарушения нейтральности могут иметь место. При этом перестают быть справедливыми уравнения модели, применяемые «по аналогии» с химическими системами. Отдельно встает проблема адекватного описания слабых электрических и электромагнитных воздействий на биологические объекты. В таких случаях необходим вывод уравнений математической модели на основе фундаментальных физических законов и обобщения экспериментальных данных. Попытка такого рассмотрения процессов вблизи клеточной мембраны предпринимается в настоящей главе.
Отметим, что рассматриваемые в данной главе задачи имеют практический интерес не только для описания процессов в окрестности клеточных мембран. Как уже было упомянуто выше, электрические взаимодействия существенны при моделировании агрегации тромбоцитов [78]. Кроме того, как отмечалось выше, во многих биологических задачах существенны модели переноса при наличии тонкого слоя, в котором процесс переноса существенно отличен от основной массы жидкости. В качестве примера приведем недавнюю работу [247], в которой при описании гемостаза введены «химический пограничный слой». В последнем процесс переноса метаболитов существенно отличен от процессов переноса в основном объеме сосуда. Кроме этого, в [247] вводится и «тромбоцитный пограничный слой» (большей толщины), в котором существенно отличен от основного объема перенос тромбоцитов. В данной главе строятся простые математические модели, которые могут найти свое приложение и при математическом моделировании процесса свертывания крови.
В первом параграфе настоящей главы рассмотрена модель процессов электродиффузии в примембранном слое [248]. В модели учитывались происходящие на мембране химические реакции, кулоновские взаимодействия, случайные перемещения частиц (диффузия) и влияние внешнего электрического поля. На основе линейного анализа системы сделан вывод, что под влиянием внешнего электрического поля в системе возможно образование диссипативной структуры (ДС), которая отличается от ДС в случае отсутствия внешнего поля следующими характеристиками. Во-первых, она медленно движется (дрейфует) в одном направлении, при этом меняется число структурных элементов. Во-вторых, в дисперсионное соотношение входят и нечетные степени волнового числа к, что может приводить к появлению изолированных солитонообразных структур. Образование ДС может быть обусловлено проявлением не только диффузионной неустойчивости типа Тьюринга, но и дисперсионной неустойчивостью под воздействием внешнего электрического поля.
Начиная с классической работы Тьюринга [5] и работ школы Пригожина -Гленсдорфа - Николиса стало понятным, что изменение параметров распределенной системы типа «реакция-диффузия» при определенных соотношениях коэффициентов диффузии может приводить к возникновению структур в системах с первоначально однородным стационарным состоянием [76,77, 134]. Бифуркации в таких системах представляют собой изменение типа пространственно-временного режима. Например, возникновение диссипативной структуры из гомогенного стационарного состояния, переход стационарной структуры к бегущим или стоячим волнам.
