CC BY
54
СЛОЖНАЯ ВОЛНОВАЯ ДИНАМИКА АНСАМБЛЯ НЕЙРОНОПОДОБНЫХ ЭЛЕМЕНТОВ СО СЛОЖНОПОРОГОВЫМ ВОЗБУЖДЕНИЕМ
В.И. Некоркин, Д.С. Щапин, А.С. Дмитричев
Проведено исследование пространственно-временной динамики системы, моделирующей коллективное поведение ансамбля электрически связанных нервных клеток. Моделью элемента ансамбля является уравнение ФитцХью - Нагумо со сложнопоро-говым возбуждением. Изучены гетероклинические траектории системы и ассоциирующиеся с ними волновые фронты. В системе для бегущих волн обнаружено существование гетероклинического контура, образованного сепаратрисными многообразиями двух седло-фокусов. Показано, что наличие такого контура свидетельствует о сложной пространственно-временной динамике ансамбля - в виде ромбоподобных и нестационарных колебательных волновых структур.
Введение
Нелинейные волновые процессы, играющие фундаментальную роль в физике и математике, имеют также важное значение для нейронных систем [1-3]. С ними связано выполнение основной задачи нейронных систем - межнейронной коммуникации и обработки информации [4-7]. В нейронных системах наиболее известны и изучены волновые явления в виде бегущих фронтов и импульсов. Как правило, характеристики этих волн зависят только от локальных свойств среды и в определенных пределах не зависят от начальных условий, а для взаимодействия волн между собой и с границами среды типично свойство аннигиляции. Такие волны принято называть автоволнами [8]. Наряду с автоволнами в нейронных системах существуют волны, рассмотренные в работах [9-11], свойства и характеристики которых принципиально отличаются от автоволн. Они обладают частицеподобными свойствами и при взаимодействии друг с другом и с границами среды могут «переключаться» в новое состояние и образовывать сложные, в том числе, фрактальные [12,13], пространственно-временные структуры активности и др. [14-16].
Эта работа посвящена исследованию нелинейных волновых структур в системе следующего вида:
( йу = /(и) - V] + й{щ-х - 2щ + и^+г), \ = е( 9(из) - vj -1),
] = 1,2,... И, йо(г) = йг(г), UN+l(í) = им(г),
где
f (u) = u — u3/3, g(u)
au, u < 0, eu, u > 0.
(2)
Заметим, что производная функции g(u) имеет точку разрыва u = 0, в которой происходит скачок на величину ||3 — a| и поэтому система (1) требует доопределения при uj = 0, j = 1, 2,... N. Будем считать, что g'(0) = |3. В этом случае траектории системы (1) будут непрерывны, а соответствующая производная при uj = 0 имеет конечный скачок.
Система (1) описывает ансамбль нейроноподобных элементов типа Фитц-Хью - Нагумо со сложнопороговым возбуждением - в модели элемента пороговые свойства определяются сепаратрисой седла, сложным образом локализованной на фазовой плоскости [17,18]. Исследование системы (1) представляет интерес в связи с двумя проблемами нейродинамики. Первая из них - выяснение механизмов коллективного возбуждения в виде кластеров активности в ансамбле электрически связанных нейронов. В этом случае переменная uj качественно описывает динамику мембранного потенциала j-го нейрона ансамбля, а Vj - совокупное действие всех ионных токов, проходящих через мембрану этого нейрона и отвечающих за восстановление покоя мембраны. Поэтому переменную Vj часто называют восстанавливающей. Параметр d моделирует действие электрического синапса и характеризует величину взаимодействия между нейронами. Параметр I контролирует уровень деполяризации мембраны, а е (е > 0) - скорость изменения ионных токов. Для этой проблемы важны возбудимые свойства нейронов и поэтому здесь параметр е выбирается достаточно малым. Параметры a и в (a, в > 0) описывают нелинейные свойства ионных токов. Вторая проблема, ассоциирующаяся с системой (1), это исследование процессов распространения волн активности в нейронных ансамблях. В этом случае переменные uj и Vj являются интегральными характеристиками некоторой группы нейронов (uj - усредненное по группе значение мембранного потенциала, Vj - ионных токов). Для данной задачи важны бистабильные свойства среды и переходы «покой - возбужденное состояние». Поэтому параметр е может быть здесь и не малым. В этой статье мы уделяем основное внимание второй проблеме и номер j элемента ансамбля будем трактовать как пространственную координату.
Статья структурирована следующим образом. В разделе 1 кратко излагается динамика элемента ансамбля. В разделе 2 вводится система для бегущих волн и исследуются ее основные свойства. Раздел 3 посвящен исследованию гетероклини-ческих траекторий системы для бегущих волн. В разделе 4 изучается гетерокли-нический контур и ассоциирующиеся с ним волновые структуры. В Заключении представлено краткое обсуждение полученных результатов.
Рассмотрим систему, описывающую динамику отдельного элемента ансамбля
1. Динамика элемента ансамбля
u = f (u) — v, V = е( g(u) — v — I).
Исследование динамики системы (3) проведено в [18]. Приведем здесь краткое описание необходимых в дальнейшем динамических режимов системы (3) и их бифуркаций.
Будем рассматривать систему (3) в области параметров, в которой она одновременно имеет три состояния равновесия. Из (2) и (3) следует, что это условие выполняется, если параметр
I е (-
где 1п
2(1 - в)3/2
1ш
2(1 - а)3/2
3
При выполнении (4) система (3) имеет следующие состояния равновесия:
01« А О2(и02,у°2) и ОзСи°^з°), где
и? = -2^1 - асов(-),
V0 = аи1,
и°2 =<
-2^1 - асов(Ф + -3), = аи°, I < 0, 33
-2VT—всо8(З + у), у°° = ви°, 1> 0,
и
3 = -2^1-3 сов( | + 2Л), ^3 = ви°,
сов ф
-31
2(1 - Р)3/2'
соэ —
-31
2(1 - а)3/2'
(5)
Состояние равновесия О2 является седлом, а состояния равновесия О1 и О3 являются узлами или фокусами, которые в зависимости от параметров могут быть устойчивыми (верхний индекс 8) или неустойчивыми (индекс и). Зафиксируем для определенности параметры
а = 0.8, в = 0.9, I = 0.024,
(6)
а параметр е будем рассматривать в качестве контрольного. Разбиение полупрямой е > 0 на области, соответствующие различным структурно-устойчивым фазовым портретам системы (3), осуществляется бифуркационными значениями е = ег (г = 1,..., 6) (рис. 1). Фазовые портреты, отвечающие этим областям параметров системы (3), представлены на рис. 2: (а) соответствует значениям 0 < е < е1, (б) е1 < е < е2, (в) е2 < е < е3, (г) е3 < е < е4, (д) е4 < е < е5, (е) е5 < е < е6, (ж) е > е6. В точках, разделяющих эти области, происходят следующие бифуркации:
Рис. 1. Однопараметрическая бифуркационная диаграмма системы (3) для а = 0.8; в = 0-9 I = 0.024, где ех = 0.5019...; £2 = 0.5762 ... £з = 0.5881...; £4 = 0.6293...; £5 = 0.6338... £6 = 0.6640...
• при е = £г состояние равновесия Ог изменяет свою устойчивость и в результате бифуркации Андронова - Хопфа рождается (в сторону е > ег) неустойчивый предельный цикл I";
• при е = е2 цикл I" «влипает» в «малую» петлю, образованную сепаратрисами W1U и Wl седла О2;
• при е = е3 происходит бифуркация Андронова - Хопфа, в результате которой неустойчивый фокус О3 становится устойчивым и на фазовой плоскости появляется (в сторону е > ез) неустойчивый предельный цикл 1";
• при е = е4 смыкаются сепаратрисы W2U и W1s, образуют так называемую «большую» петлю сепаратрис (гомоклиническую орбиту), из которой рождается (в сторону е > е4) «большой» неустойчивый предельный цикл Ьи, охватывающий все три состояния равновесия;
• при е = е5 предельные циклы Ьа и ЬБ сливаются, образуя двукратный предельный цикл, который исчезает при е > е5;
• при е = е6 образуется «малая» петля (гомоклиническая орбита), образо-
Ч и W|,
цикл I".
ванная сепаратрисами W2U и в которую «влипает» неустойчивый предельный
2. Система для бегущих волн
Рассмотрим бегущие волны системы (1) в длинноволновом приближении, когда характерные пространственные масштабы волны значительно превосходят собственные пространственные масштабы ансамбля (1). Будем искать решение системы (1) в следующем виде:
Uj (t) = u(4),Vj (t) = v(l), (7)
где ^ = t — jh (h > 0) - бегущая координата. Подставляя (7) в систему (1), получим
du
— = f (u) — v + d(u(l + h) — 2u£) + u(l — h)), ^ (8)
— = e(g(u) — v — I).
В рамках длинноволнового приближения из (8) имеем
du f ( ) + d2u
^н (9)
с^ = e(g(u) — v —1),
dsн
где = с^, а 0 < с < l/(hVd). Система (9) является системой для бегущих волн.
С помощью замены г = V - д(и) +1 система (9) преобразуется к следующему виду:
и = у,
у = су + ^(и) - I + г, сг = -ег - сд'(и)у,
(10)
где ^(и) = -/(и) + д(и); точкой обозначено дифференцирование по Система (10) имеет трехмерное фазовое пространство, которое обозначим через С. Нетривиальным ограниченным траекториям системы (10) отвечают бегущие со скоростью V = 1/Н волны в ансамбле (1), профили которых определяются свойствами этих траекторий [19]. Например, если система (10) имеет в С предельный цикл, то в (1) может распространяться периодическая волна; гомоклинической к состоянию равновесия траектории отвечает бегущий импульс; гетероклинической траектории - бегущий фронт (кинк) и т. д. В этой работе мы покажем, что в системе (10) возможно существование гетероклинического контура, в окрестности которого существуют нетривиальные гиперболические множества, с которыми ассоциируются сложные волновые структуры ансамбля.
2.1. Состояния равновесия системы (10) и их локальные бифуркации.
При выполнении (4) система (10) имеет три состояния равновесия О^и^, 0,0), 02(и2, 0, 0) и 03(и3, 0, 0), где и° определяется выражением (5). Как и в разделе 1, зафиксируем параметры а, в и I (см. условия (6)) и будем рассматривать параметрическую задачу на плоскости (е, с). На рис. 3 представлены результаты анализа
Рис. 3. Разбиение плоскости параметров (£,с) на области, соответствующие различному качественному расположению корней характеристического уравнения: а - для состояния равновесия Ох, б - для 02, в - для 03. Значения параметров а = 0.9; в = 0.8; I = 0.024
качественного расположения корней характеристического уравнения на комплексной плоскости для состояний равновесия О1, О2 и О3. Для каждого из состояний равновесия плоскость (е, с) разбивается на три области линиями А^ и Si (г = 1, 2, 3):
А = с = </ ^^ , 0 * е^К)
е(Р - е)
А2 = с Ч -щ-е,0 * е<в
Аз = с = ,/ , 0 * е<Г(4)
, №§) - е
Si = {-р\2 + 4p3ri + 4д3 - 18рдг + 27г2 = 0} ,
Р = - - c,qi = —е + ¡'(и0),П = -Е'(и0).
Из представленного на рис. 3 расположения корней характеристических уравнений для О1, О2 и Оз следует, что О1 и Оз являются либо неустойчивыми фокусами, либо седловыми состояниями равновесия (седло-фокусами или седлами) с одномерным устойчивым и двумерным неустойчивым многообразиями, а О2 - либо устойчивый фокус, либо седловое состояние равновесия (седло-фокус или седло) с одномерным устойчивым и двумерным неустойчивым многообразиями. Отметим, не останавливаясь подробно, что при пересечении линий Ai в системе (10) происходит бифуркация Андронова - Хопфа [20], приводящая к рождению предельных циклов в окрестности соответствующих состояний равновесия 0i. При этом около состояний равновесия О1 и О3 появляются неустойчивые, а около О2 - устойчивый предельные циклы.
Поскольку нас интересуют гетероклинические траектории, образованные многообразиями состояний равновесия О1 и О3, далее будем рассматривать систему (10) в области параметров, где О1 и Оз являются седловыми состояниями равновесия. Обозначим через ШБ(01) и Ши(О1) (Ш8(03) и Ши(03)) устойчивое и неустойчивое многообразия состояния равновесия О1 (соответственно О3). Многообразие Ши(О1) (Ши(О3)) состоит из точки О1 (О3) и двух выходящих траекторий - сепаратрис, которые обозначим через Ш1и(О1) и Ш21(О1) (Ш"(О3) и Ш2и(О3)). Двумерное многообразие ШБ(01) (Ш8(03)) представляет собой поверхность, проходящую через точку О1 (О3) и состоящую из траекторий, асимптотически стремящихся при Ь ^ к состоянию равновесия. Анализ собственных векторов матриц, линеаризованных в точках О1 и О3 систем, показывает, что неустойчивые многообразия Ши(О1) и Ши(О3) касаются в точках О1 и О3, соответственно, прямых
и - и0 = У = £ (11)
1 X(0i) к' К '
где A,(0i) - положительный корень характеристического уравнения для состояния равновесия О^ а
к = X2(0i) - ск(Ог) - ^(и0), г = 1,3.
Из (11) следует, что одна из сепаратрис Ж!и(О!) (соответственно (О3)) выходит в полупространство у > 0, а другая Ж"(О!) (Ж"(О3)) - в полупространство у < 0.
2.2. Поверхности без контакта. Проведем исследование нелокального поведения многообразий седловых состояний равновесия с помощью «конических» бесконтактных поверхностей [12,21,22]. Введем в рассмотрение функцию
Z2 У2 Г
V(u, y,z) = — - (2е - вг2 + (2е - в) J (F(n) - I)dn.
(12)
uo
Производная функции V(u, y, z) в силу системы (10) имеет вид
е
V = --z2 - [2е - в + g'(u)]zy - с(2е - в)у2. (13)
Нетрудно видеть, что при е > в будет выполняться неравенство
V(u,y,z) < 0. (14)
При этом = 0 лишь в точках прямой {y = z = 0}. Следовательно, при е > в на любой поверхности уровня V(u, y,z) = C, C = const (за исключением точек {V(u, y,z) = C} П {y = z = 0}) векторное поле системы (10) ориентировано внутрь области, границей которой является эта поверхность. Из (12) следует, что сечения поверхностей V(x, y,z) = C плоскостями u = const представляет собой либо гиперболы, либо две пересекающиеся в точке y = z = 0 прямые. С другой стороны, плоскости y = const в зависимости от значения этой константы задают в сечении либо одну, либо две замкнутые линии (в частном случае замкнутые линии могут вырождаться в точку). Граничному, разделяющему эти случаи значению y = const соответствует существование в сечении замкнутой линии с одной точкой самопересечения (типа «восьмерки»). Таким образом, поверхности V(u, y,z) = C являются замкнутыми. Траектории системы (10) в силу (14) пересекают эти поверхности в одном направлении, так что вдоль траекторий значения C, соответствующие поверхностям уровня, убывают. Поэтому после пересечения какой-либо поверхности уровня V(u, y, z) = C траектория системы (10) не может покинуть область фазового пространства G, ограниченную этой замкнутой поверхностью. Следовательно, при е > в система не имеет замкнутых траекторий.
Установим с помощью поверхностей V(u, y, z) = C поведение сепаратрис седловых состояний равновесия Oi, O3. Для этого рассмотрим поверхности уровня, проходящие через эти точки
V (Oi) = {V (u,y,z)=0},
V (O3) = {V (u,y,z)= Co}, u3
Co = (2е - в) J (F(n) - I)dn .
u
Рис. 4. а - качественный вид поверхности V + (О1) при I > I* (I* = 0, 0168), б - качественный вид поверхности V + (О3) при I > I*
Существует два принципиально различных вида этих поверхностей: для значений I > I* и I < I*, где I* определяется условием
"3
J (п) - Т)дц = 0.
0
1
Для значений параметров, удовлетворяющих (6), I* = 0.0168.
Пусть I > I*. В этом случае поверхность V(0\) при у > 0 1 имеет качественный вид, представленный на рис. 4, а, а поверхность V(03) - на рис. 4, б. V(03) состоит из Оз и двух компонент V +(0з) и V-(0з), расположенных соответственно при у > 0 и у < 0. В силу (14), сепаратриса ЖЦ(03) выходит из точки 03 внутрь области фазового пространства, ограниченной поверхностью V+(03). Поскольку область целиком расположена при у > 0, в силу первого уравнения системы (10), переменная и(Ь) вдоль сепаратрисы Ш^(03) непрерывно возрастает. Кроме того, согласно (14), сепаратриса ЖЦ(03) пересекает поверхности V(и, у, г) = С (где С < С0), расположенные внутри поверхности V + (03), в одном направлении. Следовательно, сепаратриса ЖЦ(03) при Ь ^ стремится в бесконечность {и(Ь),у(Ь) ^ оставаясь внутри области, ограниченной поверхностью V +(03). Совершенно аналогично можно показать, что сепаратриса ЖЦ (03) при Ь ^ стремится в бесконечность {и(Ь),у(Ь) ^ не выходя из области, ограниченной поверхностью V-(03). Таким образом, при I > I*, е > в обе сепаратрисы состояния равновесия 03 гомоклинических и гетероклинических орбит образовывать не могут.
При I < I*, е > в поверхность V(01) имеет лишь одну точку пересечения с плоскостью {у = 0} - состояние равновесия 01, а V(03) пересекает {у = 0} в точке 03 по замкнутой линии, охватывающей состояние равновесия 01. Принимая во внимание свойства поверхностей V(01) и V(03), устанавливаем, что в этом случае не могут образовывать гомоклинические и гетероклинические траектории сепаратрисы состояния равновесия 01.
1В силу инвариантности функции V(и, у, г) относительно замены у ^ —у, качественный вид поверхностей V(О1) и V(О3) при у < 0 будет точно такой же, как при у > 0.
Нетрудно видеть, что при I = I* поверхности V(Oi) и V(Оз) сливаются (в этом случае C = 0), образуя единую поверхность, пересекающую плоскость y = 0 в точках Oi и Оз. Поэтому при I = I*, е > в все сепаратрисы состояний равновесия Oi и Оз стремятся в бесконечность и система (10) гомоклинических и гетероклини-ческих траекторий не имеет.
2.3. Динамика системы (10) при c ^ 1. Представим для удобства систему (10) в виде одного уравнения
с v +(е - е2)й + c[f'(u) - е]и - e[F(и) - I] = 0. (15)
Введем в (15) новую бегущую переменную cZ = и представим в виде следующей системы:
du
dZ =Ун'
dz=*н, (16)
dz
ЦdZT = (! - + [е - f'(и)]Ун + e[F(и) - I],
где ц = 1 /с2. При с » 1 система (16) имеет малый параметр ц при старшей производной и ее динамика является релаксационной. Применяя к системе (16) методы [23,24] релаксационных колебаний, устанавливаем, что в фазовом пространстве G существует двумерное неустойчивое многообразие медленных движений, близкое к поверхности
гн + (е - f'(u)y + e(F(u) - I) = 0. (17)
На этом многообразии движения системы (16) определяются двумерной системой, близкой к системе
{й = ун,
(18)
У = f (и)Ун - е^(u) - I + Ун).
Вне малой окрестности неустойчивого многообразия движения системы (16) происходят по быстрым траекториям, близким к прямым {и = и0, Ун = У°}, где и0, У° - const. По этим траекториям изображающая точка удаляется от многообразия медленных движений и уходит в бесконечность (|гн| ^ ж). С помощью замены Ун = f (и) - v система (18) может быть приведена к виду (3). Следовательно, на двумерном неустойчивом многообразии системы (16) движения близки к движениям системы (3), представленным на рис. 2. Так что для получения разбиения фазового пространства G на траектории при с ^ 1 достаточно к соответствующим фазовым портретам системы (3) «подклеить» быстрые неустойчивые траектории, близкие к прямым. В частности, отсюда вытекает, что при с ^ 1 все сепаратрисы седловых состояний равновесия Oi, Оз уходят на бесконечность и гомоклинические и гетеро-клинические траектории они не образуют.
3. Гетероклинические траектории
Напомним, что для состояний равновесия траектория называется гетерокли-нической, если она при ^ ^ и ^ ^ -то стремится к различным состояниям равновесия. Покажем, что такие траектории в системе (10) существуют. Для определенности рассмотрим систему при I > I* в области параметров е > р. В этой
области параметров, как установлено в разделе 2, сепаратрисы Жи(01), ЖЦ(03) и ЖЦ1 (О3) уходят в бесконечность и гетероклинических траекторий не образуют. Рассмотрим теперь поведение сепаратрисы Ж1и(01). Если Ж^(0-\) образует гетеро-клиническую траекторию, то она принадлежит устойчивому многообразию Ж8(0з) состояния равновесия О3. Следовательно, существование такой гетероклинической траектории можно установить, исследуя взаимное расположение ЖЦ(0-\) и Ж№(03) для различных значений параметров. Изучим взаимное расположение этих многообразий с помощью специальным образом построенных двумерных систем сравнения [19,25,26].
3.1. Системы сравнения. Рассмотрим систему
и = у,
у = су + ¥ (и) — I ,
(19)
где I - параметр. При I е (1тт>1тах) система (19) имеет три состояния равновесия: 01(и1, 0), 02(и2, 0) и 03(и3, 0), где щ определяются выражениями (5), в которых I нужно заменить на I. На рис. 5 представлено (см. [12]) разбиение плоскости параметров (I, с) на области, соответствующие различным фазовым портретам системы (19). Разбиение осуществляется прямыми I = 1тах, I = — !т1П и двумя би-
-0.035 /т1П / -0.001 /
Рис. 5. Разбиение плоскости параметров и фазовые портреты системы (19) (а = 0.8; в = 0.9; С0 = 0.450...; с2 = 0.494 ...)
фуркационными линиями с = Л+(1) и с = Ь-(1). Прямая I = 1тах соответствует слиянию состояний равновесия 01 и 02 и образованию седло-узла 01>2(-V1 - а, 0) с неустойчивой узловой областью. При I = -1т;п седло-узел 02,3(-V1 - в, 0) с неустойчивой узловой областью образуется от слияния состояний равновесия О2 и О3. Линия с = К+ (I) соответствует существованию гетероклинической траектории, образованной сепаратрисами седел О1 и О3, расположенными в полуплоскости у > 0, а с = Н-(!) - гетероклинической траектории, образованной сепаратрисами тех же седел, но в полуплоскости У < 0.
Введем в рассмотрение две системы сравнения
и = у,у = су + ^(и) - I + го, (20)
и = у, у = су + ^(и) - I - го, (21)
где го > 0 - произвольный параметр. Каждая из этих систем представляет собой систему вида (19), у которой параметр I = I - г0 и I = I + г0, соответственно. Нетрудно видеть, что при выполнении условия
го < го*, (22)
где г0 = тш{(!тах - I), (I + !тт)}, системы (20), (21) имеют по три расположенных на прямой у = 0 состояния равновесия. При этом координаты этих состояний равновесия удовлетворяют неравенствам
и+ < и° < и- < и- < и2 < и+ < и- < и° < и+, (23)
где знак «+» относится к состояниям равновесия системы (20), а «-» - системы (21).
Потребуем, чтобы параметры систем (20), (21) одновременно принадлежали области 1 или области 2 (см. рис. 5). Используя разбиение плоскости параметров (I, с), представленное на рис. 5, получаем, что параметры этих систем лежат в области 1, если выполняется неравенство
сЖ^^ + го), (24)
и в области 2, если
с <11+^ - го), (25)
где го удовлетворяет (22). При выполнении (24) или (25) обе системы имеют качественно одинаковые фазовые портреты, однако их состояния равновесия на прямой у = 0 смещены относительно друг друга (см. (23)). Взаимное расположение седел систем (20), (21) и их сепаратрис для областей параметров (24) и (25) имеет вид, представленный на рис. 6, а и б, соответственно. При этих значениях параметров сепаратрисы седел систем (20), (21) образуют так называемые [24,25] направляющие русла, с помощью которых мы построим в фазовом пространстве С инвариантные области, содержащие многообразия седловых состояний равновесия О1 и О3 системы (10). Поскольку системы (20), (21) не зависят от переменной г, в фазовом пространстве С сепаратрисы седел этих систем образуют цилиндрические поверхности. Рассмотрим ориентацию векторного поля системы (10) на этих поверхностях. Из (10), (20) и (21) имеем
dy\ _ (dy\ _ z — zo du) (10) V du) (20) У
< 0, при z < z0, y > 0,
(26)
-(^ _ ^ > 0, при z > -zo,y> 0. du (10) du (21) У
В силу (26) траектории системы (10) при y > 0 пересекают цилиндрические поверхности, образованные сепаратрисами седел системы (20) в сторону убывания координаты y, а сепаратрисами седел системы (21) - в сторону возрастания координаты y (см. рис. 6). Обозначим через g1 (g2) область фазового пространства G, расположенную при y > 0 и ограниченную цилиндрическими поверхностями, образованными выходящими (соответственно входящими) сепаратрисами седел систем (20), (21) и плоскостями z _ zo и z _ — zo (см. рис. 6). Рассмотрим ориентацию векторного поля системы (10) на этих плоскостях. Из системы (10) имеем
i £z0 zlz=zo < 0, если y >---,
ср (27)
i £z0 v '
z|z=-zo > 0, если y < —.
Следовательно, при выполнении неравенства
ы < с0 (28)
траектории системы (10) пересекают плоскости z _ ±z0 внутрь областей gi и g2. Поскольку седла 01 и 03 расположены внутри соответственно g1 и g2 (см. (23)), из установленной выше ориентации векторного поля (10) на границе gi и g2 следует, что сепаратриса Ж1и(01) и многообразие Ws(0s) остаются локализованными внутри областей gi и g2 при выполнении (28).
3.2. Функция расщепления. Рассмотрим «следы» областей g1 и g2 на плоскости {u _ u0}, то есть области S1 _ g1 П {u _ u°} и S2 _ g2 П {u _ u°}. Граница области S1 (S2) состоит из двух прямых y _ y+(c,I — z0) и y _ y-(c,I + z0) (соответственно y _ y+(c, I — z0) и y _ y- (с, I + z0)), образованных выходящими
сепаратрисами седел систем соответственно (20) и (21), и двух прямых г = ±го. Потребуем, чтобы
тах{у+(с, I - го), у-(с^ + го)} < ■ (29)
При одновременном выполнении (22), (29) сепаратриса Ш"(01) и многообразие Ш8(03) остаются внутри областей д1 и д2 до первого пересечения с плоскостью {и = Ц°}. Следовательно, при выполнении условий (22), (29) сепаратриса Ш1и(01) пересекает плоскость {и = и°} в некоторой точке М^у", г"), а Ш8(03) - по некоторой кривой №(у = у8(г)), соединяющей плоскости г = ±го. Для характеристики взаимного расположения сепаратрисы Ш"(01) и многообразия Ш8(03) введем функцию «расщепления»2
р+(е,с) = у" - у8(г")■
Функция р+(е, с) определена и непрерывна, если одновременно выполнены условия (22) и (29), то есть в некоторой области параметров А, задаваемой неравенством (29), где го изменяется от нуля до го. Нетрудно видеть (см. рис. 6), что при выполнении (24) справедливо неравенство
у + (с^ - го) >у-(с^ + го), (30)
а при выполнении (25) - противоположное (30) неравенство. Поэтому на плоскости {и = ио} область Б 1 целиком (Б 1 П Б2 = 0) расположена правее (по координате у) области Б2 при выполнении (24), а при выполнении (25) - наоборот, Б2 расположена целиком правее 51. Поскольку М| £ 51 и N8 £ 52, из такого взаимного расположения 51 и Б2 вытекают следующие неравенства
р+ (е, с) > 0, если (е, с) е А+ , р+ (е, с) < 0, если (е, с) е А+ ,
где
А+ = {ОП+! + го), е > сву+(с, ^ - го) },
го
А+ = {с<К+^ - го),е>Фу-(с,:[ + го) },
2 го
0 < го < го ■
Заметим, что границы областей А+ и А+ определяются кривыми, заданными в параметрическом виде, в котором роль варьируемого параметра играет го. Причем при го ^ 0 эти кривые приближаются к значению с = ). Вид областей А+ и А+ на плоскости (е, с) представлен на рис. 7, а.
В силу непрерывности функции р+ (е, с), из (31) следует, что при изменении параметров с переходом из области, например А+, в область А+ функция р+(е, с) по крайней мере один раз обращается в нуль. Следовательно, существует бифуркационное множество Н + = {р+(е, с) = 0}, соответствующее в фазовом пространстве С
2Конечно, функция р+(е,с) зависит не только от е и с, но и от всех остальных параметров модели (10), однако, поскольку все остальные параметры были зафиксированы (см.(6)), эту зависимость мы не указываем.
(31)
(32)
0.0 4.0 8.0 12.0 6 0.0 4.0 8.0 12.0 е
а о
Рис. 7. а - области А+ и А+ (I* <1 = 0.4; I* = 0.0168); б - области А- и А- (I* > I = 0.1)
системы (10) гетероклинической кривой, образованной сепаратрисой Ш" (О1) и многообразием Ш Б(О3).
Совершенно аналогично, в случае I < I*, вводя функцию расщепления р- (е, с) для характеристики взаимного расположения Ш2 (О3) и Ш8(О1), можно показать существование бифуркационного множества Н- = {р-(е,с)}, отвечающего гетероклинической траектории, образованной этими многообразиями. На рис. 7, б представлены области А- и А-, для точек которых функция расщепления р-(е, с) имеет разные знаки.
4. Гетероклинический контур и сложные волновые структуры
Значения бифуркационных множеств Н+ и Н- были найдены нами путем численного построения функций расщепления р+ (е, с) и р-(е,с) и изучения их нулей. Для этого функции расщепления были определены и для значений параметров е < в (напомним, что в разделе 3 выполняется условие е > в). Например, функция р+(е, с) вводилась следующим образом. В фазовом пространстве С системы (10) строилась вспомогательная сфера 2 достаточно малого радиуса с центром в состоянии равновесия О3, на которой численно находился «след» сепаратрисы Ш1и(О1), то есть точка Ти = Ши(О{) П 2. В качестве функции р+ (е, с) выбиралось минимальное расстояние между точкой Ти и плоскостью Ш[ос(О33) - касательной плоскостью к многообразию Ш13 (О3). Аналогично, при е < в была введена функция расщепления р-(е, с).
Результаты численного исследования множеств Н + и Н- представлены на рис. 8, а. На плоскости (е, с) множества Н + и Н- представляют собой бифуркационные кривые, имеющие общую точку Н0(е0 ,с0). Точке Н0 отвечает в фазовом пространстве С гетероклинический контур С, образованный сепаратрисами Ши(О1) и Ш2 (О3) и многообразиями ШБ(О3) и ШБ(О1), соответственно (рис. 8, б). Бифуркационное множество системы (10) в окрестности Н0 не исчерпывается элементами Н + и НПокажем, что существуют бифуркационные значения параметров, соответствующие гомоклиническим траекториям (петлям сепаратрис), образованным многообразиями состояний равновесия О1 и О3. Будем варьировать параметры е и с в окрестности Н0, например, вдоль кривой Нтак чтобы выполнялось неравенство р+ (е, с) < 0. Такое изменение параметров приведет к тому, что сепаратриса Ш1и(О1)
Рис.8. а - бифуркационные множества Н+ , Н и бифуркационные множества Г+, Г в окрестности точки Н0(е0 = 0.60593, с0 = 1.0664); б - гетероклинический контур, соответствующий точке Н0 (а = 0, 9; в = 0, 8; I = 0, 024)
«покинет» многообразие Ш13 (Оз) и начнет двигаться в О в окрестности сепаратрисы Ш2и(Оз). С другой стороны, в пространстве параметров существуют области, для которых взаимное расположение сепаратрисы ЖЦ (Оз) и Ш8(О\) характеризуется тем, что р"(е, с) > 0 и р"(е, с) < 0. Отсюда, поскольку сепаратриса Ш1(О1) локализована в окрестности сепаратрисы ЖЦ (Оз), вытекает существование областей параметров, для которых сепаратриса Ш!и(О!) расположена по разные стороны от многообразия ЖБ(О!). Очевидно, что между этими областями существует бифуркационное множество Г+, отвечающее гомоклиническим траекториям системы (10), образованным сепаратрисой Ш1и(О!), возвращающейся в О! по многообразию Ж8(О!). Совершенно аналогично устанавливается существование в окрестности Н0 множества Г", соответствующего гомоклиническим траекториям, образованным сепаратрисой Ш2Ц(Оз). Результаты численного построения Г+ и Г" представлены на рис. 8, а. При выполнении (6) в точке Н0 и ее окрестности состояния равновесия О1 и Оз являются седло-фокусами (см. раздел 2.1), для которых так называемые седловые величины
о(Ог) = 1з(Ог) + ИеХ^Ог), I = 1,3,
где А,з(Ог) и Х12(Ог) - корни характеристического уравнения для Ог, являются положительными. В соответствии с теоремой Шильникова [27] такие гомоклинические орбиты ассоциируются с хаотической динамикой, так как в окрестности каждой из петель Г+ и Г" (как в момент существования, так и при разрушении в обе стороны) существует нетривиальное гиперболическое множество, включающее бесконечное множество периодических седловых траекторий и др.
Таким образом, для значений параметров из окрестности Н0 система (10) для бегущих волн демонстрирует чрезвычайно сложную динамику и можно ожидать, что и пространственно-временное поведение всего ансамбля (1) в этом случае будет нетривиальным. Проведенное нами численное моделирование динамики ансамбля (1) подтвердило справедливость этого утверждения. Был рассмотрен ансамбль (1) из шестисот элементов (М = 600), для которого численно решалась задача с начальными условиями для различных значений параметра е. Во всех случаях начальные условия были одни и те же и имели для (0) и Vj (0) вид, представлен-
ный на рис. 9, а. В случаях, когда параметр е был выбран «вдали» от точки Н 0, в ансамбле (1) распространялись лишь простые волновые фронты, которые при столкновении друг с другом аннигилировали и переводили ансамбль в пространственно-однородное состояние, соответствующее либо состоянию равновесия О3 (рис. 9, б) либо - О1 (рис. 9, в). Совершенно другое поведение демонстрирует ансамбль (1), когда значение е выбрано в окрестности точки е = е0 .В этом случае наличие нетривиального гиперболического множества и, следовательно, бесконечного множества неустойчивых периодических волн и других нетривиальных волновых движений в ансамбле (1) приводили к тому, что волновые фронты принципиально изменяли свою динамику. Теперь волновые фронты «проходят» друг сквозь
Рис. 9. а - начальные условия для системы (1) (переменным Пу (0) соответствуют значения, отмеченные сплошной линией, а переменным Уу(0) - пунктирной, ^ = 1,..,^, N = 600); б - динамика в случае е =1 (градация светлого тона соответствует значению иу = и0); в - динамика волновых фронтов в случае е = 0.71 (градация темного тона соответствует значению и у = и0); г - фракталоподобная волновая структура при е = 0.590; д - нестационарная колебательная волновая структура при е = 0, 575. Значения параметров: а = 0.9; в = 0.8; I = 0.024; й = 1.0
друга и отражаются от границ подобно классическим солитонам. В результате такого поведения в ансамбле (1) образуются сложные волновые структуры. Два примера таких структур представлены на рис. 9, г, д. Волновая структура на рис.9, г является ромбоподобной, а на рис. 9, д - нестационарной колебательно-волновой, которая устанавливается в результате сложного взаимодействия волновых фронтов.
Заключение
Проведено исследование пространственно-временной динамики ансамбля (1), моделирующего коллективное поведение электрически связанных нейронов, обладающих колебательными свойствами. Обнаружено, что коллективная активность ансамбля существенно зависит от значения параметра е, который управляет скоростью изменения ионных токов через мембрану нейрона, то есть от скорости изменения восстанавливающей переменной v. При достаточно медленном или достаточно быстром изменении v ансамбль (1) ведет себя как типичная система «реакция - диффузия» с автоволновой динамикой. В этом случае ведущую роль играют волновые фронты, взаимодействие которых приводит систему в пространственно-однородное состояние. Совершенно другую динамику демонстрирует ансамбль (1), если параметр е лежит в некоторой средней области. Даже очень простые начальные условия вызывают в ансамбле сложные пространственно-временные структуры. Для этих значений е система (1) перестает быть автоволновой и волны начинают взаимодействовать подобно классическим солитонам.
Мы установили, что сложная волновая динамика системы (1) ассоциируется с наличием у системы для бегущих волн гетероклинического контура, образованного многообразиями двух седло-фокусов, имеющих превалирующее неустойчивое направление. Другими словами, существование такого контура можно рассматривать как некоторый признак возможной сложной волновой динамики системы.
Работа поддержана грантами РФФИ (проекты 06-02-16137, 05-02-17441).
Библиографический список
1. Murray J.D. Mathematical Biology, Second Corrected Edition. Springer-Verlag, Berlin, 1993.
2. Winfree A.T. The geometry ofBiological Time. Springer-Verlag, New-York, 1980.
3. Scott A. Neuroscience: a mathematical premier. Springer-Verlag, Berlin, 2002.
4. Koch C. Biophysics of computation: information processing in single neurons. Oxford University Press, 1998.
5. Казанцев В.Б., Некоркин В.И. Динамика колебательных нейронов. Иформа-ционные аспекты // Нелинейные волны - 2002 / Ред. А.В. Гапонов-Грехов, В.И. Некоркин. Нижний Новгород: ИПФ РАН, 2003.
6. Борисюк Г.Н., Борисюк Р.М., Казанович Я.Б., Иваницкий Г.Р. Модели динамики нейронной активности при обработке информации мозгом - итоги «десятилетия» // УФН. 2002. Т. 172, № 10. С. 1189.
7. Потапов А.Б., Али М.К. Нелинейная динамика обработки информации в нейронных сетях // Новое в синергетике: Взгляд в третье тысячелетие / Ред. Г.Г. Ма-линецкий, С.П. Курдюмов. М.: Наука, 2002.
8. Васильев В.А., Романовский Ю.М., Яхно В.Г. Автоволновые процессы. М.: Наука, 1987.
9. Морнев О.А., Асланиди О.В., Алиев Р.Р., Чайлахян. Солитонный режим в уравнении ФитцХью - Нагумо: отражение сталкивающихся импульсов возбуждения//ДАН. 1996. Т. 347. С. 123.
10. Асланиди О.В., Морнев О.А. Об отражении бегущих импульсов возбуждения // Биофизика. 1996. Т. 41. С. 953.
11. Асланиди О.В., Морнев О.А. Могут ли стакивающиеся нервные импульсы отражаться? // Письма в ЖЭТФ. 1997. Т. 65. С. 553.
12. Некоркин В.И. Бегущие импульсы в двухкомпонентной активной среде с диффузией // Изв. вузов. Радиофизика. 1988. Т. 31, №1. С. 41.
13. Nekorkin V.I., Kazantsev V.B. Autowaves and solitons in three component reaction-diffusion system // Int. J. Bifurcation and Chaos. 2002. Vol. 12, № 11. P. 2421.
14. Hayase Y. Collision and self-replication of pulses in a reaction diffusion system // J. of the Physical Society of Japan. 1997. Vol. 66, № 9. P. 2584.
15. Hayase Y., Ohta T. Self replicating pulses and Sierpinski gaskets in exitable media // Phys. Rev. E. 2000. Vol. 62, № 5. P. 5998.
16. Kazantsev V.B., Nekorkin V. I., Binczak S., Bilbault J.M. Spiking patterns emerging from wave instabilities in one-dimensional neural lattice // Phys. Rev. E. 2003. Vol. 68. P. 017201.
17. Kazantsev V.B. Selective communication and information processing by exitable systems // Phys. Rev. E. 2001. Vol. 64. P. 056210.
18. Некоркин В.И., Дмитричев А. С., Щапин Д.С., Казанцев В.Б. Динамика модели нейрона со сложнопороговым возбуждением // Математическое моделирование. 2005. Т. 17, № 6. C. 75.
19. Nekorkin V.I., Velarde M.G. Sinergetic phenomena in active lattices. Springer-Verlag, 2002, 357 p.
20. Марсден Дж., Мак Кракен М. Бифуркация рождения цикла и ее приложения. М.: Мир, 1980. 368 с.
21. Максимов А.Г., Некоркин В.И. Гетероклинические траектории и фронты сложной формы модели ФитцХью - Нагумо // Математическое моделирование. 1990. Т. 2, № 2. C. 129.
22. Nekorkin V.I., Chua L.O. Spatial disorder and wave fronts in a chain of coupled Chua's circuits // Int. J. Bifurcation and Chaos. 1993. Vol. 3. P. 1281.
23. Митропольский Ю.А., Лыкова О.Б. Интегральные многообразия в нелинейной механике. М.: Наука, 1973.
24. Андронов А.А., Витт А.А., Хайкин С.Э. Теория колебаний. М.: Физматгиз, 1959.
25. Белых В.Н., Некоркин В.И. О качественном исследовании многомерной фазовой системы // Сибирский матем. журнал. 1977. Т. 18, № 4. С. 723.
26. Belykh V.N. Homoclinic and heteroclinic linkages in concrete systems: nonlocal analysis and model maps // Amer. Math. Soc. Transl. (2) 2000. Vol. 200. P. 51.
27. Шильников Л.П. К вопросу о структуре расширенной окрестности грубого состояния равновесия типа седло-фокус // Матем. сб. 1971. 81(123), № 1. С. 92.
Институт прикладной физики РАН, Поступила в редакцию 10.04.2006
Н.Новгород После доработки 06.11.2006
COMPLEX WAVE DYNAMICS OF ENSEMBLE OF NEURON-LIKE ELEMENTS WITH COMPLEX THRESHOLD EXCITATION
V.I. Nekorkin, D.S. Shapin, A.S. Dmitrichev
We present the analysis of spatiotemporal dynamics in the system modeling collective behaviour of ensemble of electrically coupled neuronal cells. The dynamics of local element is described by the FitzHugh - Nagumo system with complex threshold excitation. Heteroclinic orbits and corresponding wave fronts are investigated. We show that in the phase space of system for traveling waves there exist heteroclinic cycle formed by separatrix manifolds of two saddle-foci. It is shown that the existence of such cycle leads to complex spatiotemporal dynamics of ensemble including rhomb-like and nonstationary oscillating wave structures.
Некоркин Владимир Исаакович - родился на Украине (1948), окончил радиофизический факультет Нижегородского государственного университета (1971). Защитил диссертацию на соискание ученой степени кандидата физико-математических наук (1981, ННГУ) и доктора физико-математических наук (1992, СГУ). С 1971 года работал в Нижегородском институте прикладной математики, с 1983 по 1986 год в Нижегородском политехническом институте, с 1986 года и по настоящее время работает на радиофизическом факультете ННГУ, профессор. Область научных интересов - теория колебаний и волн, структурообразование и нелинейные волны в непрерывных и дискретных средах, синхронизация и пространственно-временной хаос. Имеет более 100 научных публикаций в указанных направлениях, в том числе две монографии (в соавторстве). В качестве приглашенного профессора читал лекции в Мадридском и Стэнфордском университетах. Работал в Калифорнийском университете. E-mail: [email protected]
Щапин Дмитрий Сергеевич - родился в Нижнем Новгороде (1981), закончил Нижегородский государственный университет им. Н.И. Лобачевского (2003). Учится в аспирантуре Института прикладной физики РАН. Область научных интересов - теория нелинейных колебаний и волн, нейродинамика.
Дмитричев Алексей Сергеевич - родился в с. Слободское Кстовского района Нижегородской области (1980), окончил Нижегородский государственный университет им. Н.И. Лобачевского (2002). После окончания ННГУ работает в ИПФ РАН.
