2011
ВЕСТНИК ПЕРМСКОГО УНИВЕРСИТЕТА
Математика. Информатика. Механика Вып. 2(6)
УДК 519.6
Численно-аналитические схемы анализа динамических систем с последействием
И. Е. Полосков
Пермский государственный университет, 614990, Пермь, ул.Букирева, 15 polosk@psu.ru
На основе схемы расширения фазового пространства рассматриваются различные типы динамических систем с последействием. Эта схема применяется для исследования линейных и нелинейных систем дифференциальных уравнений (ДУ) с одиночными и кратными постоянными запаздываниями, ДУ с переменными запаздываниями, ДУ нейтрального типа, специальных форм интегро-дифференциальных уравнений и др. Соответствующие алгоритмы реализованы в виде программ на входном языке пакета Ма'ЪЬета'Ыса. Представлен ряд примеров анализа систем различных классов.
Ключевые слова: динамическая система; последействие; фазовое пространство; расширение; метод-, моделирование; численное интегрирование
Введение
Начиная с середины прошлого века значительный интерес как с теоретической, так и практической точки зрения вызывают проблемы, связанные с анализом явлений, описываемых функционально-дифференциальными уравнениями (ФДУ) и их частными формами, такими как ДУ с запаздыванием (ДУсЗ), ДУ нейтрального типа (ДУНТ), интегро-дифференциальные уравнения (ИДУ) и др. [1-4]. Сначала это были задачи управления, а затем и биологии, механики, физики, химии, медицины, экономики, атомной энергии, теории информации и т.д.
Анализ систем, описываемых ФДУ, встречает значительные трудности. Как правило, обычные схемы исследования нелинейных систем (линейный анализ устойчивости, изучение бифуркаций, стандартные методы возмущений) непригодны для анализа ФДУ, решения которых обладают новыми свойствами (например, возможна автосинхронизация с временным запаздыванием, даже в системах первого порядка
Работа выполнена при финансовой поддержке Министерства образования и науки Российской федерации (проект № 1.25.11).
©И. Е. Полосков, 2011
допустимы колебательные режимы за счет наличия периодических запаздываний). Присутствие последействия может приводить к неустойчивой и неудовлетворительной работе технических объектов, а также их плохой управляемости.
В настоящее время существует только несколько областей исследования ФДУ, где ведутся интенсивные научных разработки. Среди них качественный анализ существования и устойчивости решений ФДУ; аналитические результаты для линейных систем, полученные методом шагов; использование процедур усреднения ДУ с учетом малости запаздывания; различные численные интеграторы и другие методики для получения частных решений задач [4-6] и др.
Наша схема анализа таких систем, как и стохастических дифференциальных уравнений с запаздыванием [7-10], базируется на комбинации метода шагов и расширении фазового пространства (МШРФП) и позволяет рассматривать различные формы ФДУ с одной точки зрения. При этом в большинстве вариантов схемы ошибка метода отсутствует. Кроме того, исчезают многие проблемы, возникающие при реализации процедур прямого численного интегрирования ДУ с запаздыванием.
В данной работе мы представляем детали вариантов схемы для анализа ДУ с кратными постоянными запаздываниями (разд. 1), ДУ нейтрального типа (разд. 2), ДУ с кусочно-
постоянными запаздываниями (разд. 3) и ДУ (ИДУ) с сосредоточенными и распределенными запаздываниями (разд. 4). Применение этой схемы демонстрируется на примерах. Инструментом вычислений являлся хорошо известный, мощный и удобный пакет компьютерной алгебры (ПКА) Mathematica [11].
1. ДУ с кратными постоянными запаздываниями
Рассмотрим систему нелинейных обыкновенных дифференциальных уравнений с кратными постоянными запаздываниями следующего вида:
x'(t) = f V (x(t), x(t -т), x(t - 2т),
x(t - VT),t), t > tv = to + VT. (1)
Здесь x G Rn - фазовый век тор, t - время, т = const > 0 ' - символ дифференцирования по t V G N f v = {fiv}T : Rn(v+1) x [to, to) ^ Rn t
x
функция, T - символ транспонирования.
x(t)
полуинтервалах (to,ti], (ti,t2], ..., (tv—i,tv] удовлетворяет системам ДУ вида
x'(t) = f o(x(t),t)
x'(t) = f i(x(t), x(t - T),t),
x'(t) = f v— 1 (x(t), x(t - T),...,
x(t - (v - 1)т),t),
(2)
где f q = {fqi}T, q = 1, 2,..., v - 1 - заданные
s Є [0, т], tq = to + q • т,
0,1, 2,...,
Aq
[tq, tq+1], sq s + tq, xq(s) x(sq
уд = жд (0) = ж^- 1 (т), ^0 = Жо,
21 =Со1( Ж1, Ж0) , £2 =Со1( Ж2, Ж1, Ж0)
Со1(Жд, Ж№-1, ..., Ж0) = {хд 1, 2, ...,
Хд-1,1, Хд-1,2, ..., Хд-1,п, ..., Х01, Хо2, ..., Х0п}Т.
Используя эти обозначения, мы построим цепочку систем ОДУ для векторов г0, г15 г2, ..., гд,
..., принадлежащих семейству вложенных фазовых пространств М" С М2" С М3" С ... С
Кп(№+1) с ...
Рассмотрим последовательность отрезков
Ді, * > 0.
00. Начнем с сегмепта До. Вектор хо(в), определенный на До, является решением системы ОДУ вида
X о(в) = f о( хо(в),во), хо(0) = Уо. (3)
Здесь и далее точка над векторами хі означает производную по переменной в.
10. Рассмотрим отрезки До апсі Д1 и представим систему ОДУ для расширенного вектора со1(х1, хо) так:
X о(в) = f о(хо(в),в^, хо(0) = Уо,
хі(в) = /і(хі(в), хо(в),вО, хі(0) = уі. (4)
V0. Теперь обратимся к временным сегментам До, ДІ5 ..., Д^ и построим систему ОДУ для вектора следующим образом:
х о(в) = / о(хо(в),в^, хо(0) = уо,
х і(в) = / ^хі(в), хо(в),ві), хі(0) = уі,
xv(s) = f v(xv(s),..., Xl(s),
xo(s),s^, Xv (0) = у v.
(5)
вектор-функции, удовлетворяющие тем же условиям, что и f а ж(£) равен у^и £ = £0.
Чтобы исследовать поведение вектора ж(£) для £ > ^0, здесь и далее используется наша схема анализа различных ДУ с запаздываниями. Эта схема основана на известной идее (классического метода шагов) преобразования ДУсЗ в обыкновенные ДУ (ОДУ), но в специальной форме. Для этой цели мы расширяем фазовое пространство рассматриваемой системы ДУ и вводим следующие обозначения:
N°. На этом этапе рассматриваются отрезки Д0, Д15 ..., Дд и строится система ОДУ для вектора гм:
Ж 0(й) = f 0( Ж0(в),в0), Ж0(0) = у0,
Ж 1(й) = f ^Ж1(й), Ж0(в),в^, Ж1(0) = у1,
ЖV(в) = ^(Ж„ («), ..., Ж1(в),
Ж0(й),й^ , XV (0) = у V,
Ж N (в) = f Д Жд (в), ...,
Жд^+1^), Жд-V (в),вд) , Жд (0) = Уд.
(6)
Несложно увидеть, что на всех этапах расчеты можно выполнять с помощью любого стандартного численного интегратора и только для в € [0, т], а следовательно, нет необходимости сохранять историю вычислений, что обычно требуется в специализированных алгоритмах решения ДУсЗ. Отсюда программная реализация рассмотренной схемы на каждом этапе с помощью любого ПКА будет включать такие шаги, как доопределение начальных условий и формирование системы ОДУ в символьной форме с
последующим интегрированием построенной системы.
Рассмотренный алгоритм был реализован в виде функции (процедуры), написанной на входном языке ПКА Mathematica, вызов которой имеет стандартный вид для этого ПКА:
AfterEffect[ eqs_, vars_, tq_, tauq_,
initVal^, t0_, tk_, delta_, NstepsPr^, maxSt_, meth_]
где eqs - список, содержащий векторы f 0, f l5 f v в стандартной для ПКА Mathematica форме, vars - вектор неизвестных, tq - имя временной переменной, tauq - имя переменной, определяющей запаздывание, initVal - вектор начальных условий для первой системы ДУ, tO — начальное время, tk - конечное время, delta - величина запаздывания, stepsPr - величина шага вывода результата в рамках отрезка длиной в запаздывание, maxSt - максимально допустимое число шагов интегрирования (параметр функции NDSolve ПКА Mathematica), meth - имя метода численного интегрирования ОДУ (параметр функции NDSolve).
Пример 1. Рассмотрим уравнение Ван дер Поля-Дюффинга:
x"(t) + 2а [ж2 (t - 2т) - 1] •x/(t)+
+ w2x(t) + Yx3(t - т) = 0, t > О,
x/(t) = 0, x//(t) =0, t G (-2т, 0],
ж(-2т) = 0,02, ж/(-2т) = 0,05,
(7)
где а, ш, 7 - постоянные. Поведение системы демонстрируется па рис.1 для а = 0,5, ш =1, 7 = 0,1, (а) т = 1, (Ь) т = 0,5.
2. ДУ нейтрального типа
Теперь обратимся к системе ДУ в форме
ж'(^) +2 ^к(^) ж'(4 - кт) =
к=1
= f V (ж(^), Ж (4 - т), ж (4 - 2т),..., ^
ж(4 — ^т),4), 4 > ^ = ^0 + ^т,
где основные обозначения соответствуют введенным выше, а матрицы-функции (4) непре-
рывны по переменной 4.
Как и в разд. 1, здесь предполагается, что фазовый вектор ж(£) на полуинтервалах (40,41], (£1^2], •••, (^^ъ^] удовлетворяет следующим системам ДУ: ж'(4) =
ж'(4) + Фи(4) ж'(4 - т) =
= f 1(ж(4), ж(4 - т),4),
(9)
V-1
ж'(4) + Е ^-1,к(4) ж'(4 - кт) =
к=1
= f V-1 (ж(4), ж(4 - т), ...,
ж(4 - (V - 1)т),4)
с начальным условием ж(^) = У0-
Если воспользоваться схемой из разд. 1, то мы можем получить цепочку систем ОДУ, которая подобна последовательности (3)-(6). Кроме того, очевидно, что правые части уравнений в обеих цепочках эквивалентны, а левые части уравнений, соответствующих (9), содержат суммы слагаемых вида
Гт
Жт(в) + Е ^гт,к(вт) Жт-к (в), (ю)
где
fc=l
І V, V + 1 ^ т ^ N.
Теперь мы можем исключить производные
хі
стей уравнений и получить следующую цепочку уравнений.
00 До
тор хо(в), определенный на До, удовлетворяет системе ОДУ
х о(в) = / +(хо(в),во), хо(0) = уо, (И)
где
/О (хо(в) во) = /о(хо(в) во) .
10 До
Ді
со1(хі, хо) можно представить так:
х о(в) = / +(хо(в),во), хо(0) = уо,
хі(в) = / +(хі(в), хо(в), (12)
вьво), хі(0) = уі
с
/ +(хі(в), хо(в), ві, во) =
= /1 (хі(в), хо(в),ві) - ^іі(ві) /+(хо(в),во).
V0. Здесь мы обозреваем временные отрезки До, Дь..., Д^ и получаем систему ОДУ для вектора в следующем виде:
х о(в) = / + ( хо(в),во), хо(0) = уо,
Xl(s) = f +(Xl(s), xo(s),
sbso), Xl(0) = У1,
xv(s) = f +(Xv (s), ..., x 1 (s), Xo(s),sv,..., s 1, so), Xv(0) = yv
Рис.1
где
/ и (х^(в), ..., x1(в), x0(в), в^, ..., в1, во)
= / v( х^ (в),..., хі(в), хо(в), вД -
V
— Ху Qvfc (^ ) / V— к (XV — к(в), ..., х1(в),
к=1
x0(в), ^ —..., в1, во) .
1\10. Теперь рассматриваются сегменты До,
Д15 ..., Д^ и строится система ОДУ для вектора
х о(в) = / +(хо(в),в^, хо(0) = уо,
хі(в) = / +(хі(в), хо(в),вьво), хі(0) = уі,
х V (в) = / ^ х„ (в),..., хі(в),
хо(в)^ , ...,вьво), ху (0) = у V,
Жд(в) = f ЛЖд (в),..., Ж1(в),
Ж0(в),вд, ...,«1,«0) , Жд (0) = уд
(14)
с
f д (жд (в), ..., Ж1(в), Ж0(в), вд,..., «1, =
= f ДЖд (в), ..., Жд^+1^), Жд-V(в), вд) -д
- Ху Qvk(sN ) f д - к( Жд-к (в), ..., Ж1(в),
к=1
Ж0(в), вд-к, ..., в1, .
Итак, цепочка построена, и ее форма пригодна для прямого численного интегрирования с помощью встроенных функций ПК А. Например, специальный вариант рассмотренного алгоритма был реализован в виде программы ЫетгЬгОЕд
на входном языке системы Ма'ЪЬета'Ыса, которая может быть вызвана также как функция КИетЕИесЬ (разд. 1).
Пример 2. Рассмотрим обобщенное уравнение Дюффинга:
х''(4) + с0х ''(4 - т) + 2а1х'(4)+
+ 2«2х'(4 - т) + ш2х(4) + ш|х(4 - т)+
+ 71Ж3(4) + 72Ж3(4 - т) = 0, 4> 0, (15)
х'(4) = 0, х''(4) =0, 4 € (-т, 0],
х(-т) = 5, х' (-т) = 0,
где е0, о*, ш*, 7* (г =1,2)- постоянные. На рис.2 приведены графики х(4) и х '(4) для системы (15) с т = 0,5, С0 = 0,75, «1 = 0,2, «2 =0 Ш1 = 2,
Ш2 = 0, 71 = 0,2, 72 = 0.
3. ДУ с кусочно-постоянными запаздываниями
Рассмотрим задачу решения нелинейных дифференциальных систем со специальной формой переменного запаздывания - кусочнопостоянной. Важность исследования систем с подобным типом запаздывания следует из возможности аппроксимации на замкнутом множестве с наперед заданной точностью любой непостоянной функции (в частности произвольного непрерывного переменного запаздывания) с помощью линейной комбинации кусочнопостоянных функций.
В данном случае система ДУсЗ принимает форму
Ж^) = f (ж(4), ж(4 - тД ж(4 - )2),...,
ж(4 - тДг), г>^0, (16)
где основные обозначения те же, что и в разделе
1, а
тг (4) = ^гд т ^д < 4 ^ ^д+1, ^д = ^0 + Ч ' т,
Рис.2
т = const У 0, Nrq Є {0,1, 2,...},
г = 1 2,...Х q = 0, 1 2,...
Например, [t - MJ, M = 0,1, 2,... - одна из регулярных форм для функций t - тг (t), где [zJ
- целая часть z [12].
Предположим, что для t Є (t0,t0),
to = min (t - Nrq т) = to - Noт, No > 0,
t^io ,r,q
x(t)
ОДУ:
x
(t) = f o(x(t),t), x(to) = У0. (17)
Задача состоит в вычислении значений вектора ж для любо го 4 ^ ^о.
Как и выше, мы можем просмотреть после-
Дд
временной точкой будет ^о, а 4д = + Ч • т, Ч =
0, 1, V-, N ...
0° - (N0 - 1)°. На всех этих этапах мы имеем сходную ситуацию, связанную с идентичностью ОДУ на сегментах Д0, Д15 ..., Дк, к = 0, 1, ..., N0 - 1:
x0(s) fo(xo(s), so),
Xl(s) = f o(Xl(s), si),
xo(0) = У0, Xl(0) = У1,
(18)
X k (s) = f o( Xk (s),Sk), Xk (0)
yk.
(N0)° - (N)°. Системы ОДУ для полуотрез-ков Дд^ Ддо+ъ •••, Дь к = N0, N0+1, ..N также имеют сходство, но их внутренние структуры, вследствие наличия запаздываний, отличаются:
Xo(s) = fo(xo(s), so), Xl(s) = fo{Xl(s), S1),
xo(0) = У0, Xl(0) = У1,
(19)
Ждо-1(в) = f о(^о-^^ вдо-0 ,
Ждо-1(0) = Удо-Ъ Ж до (в) = f (Ждо (в) Ждо-д11 (в)
Ждо-д21 ..., Ждо-да1 (в) ,
Ждо (0) = Удо,
Ждо + 1(в) = ^^о + ^в^ Ждо + 1-д12 (в)
Ждо + 1-д22 (в) ..., Ждо + 1-да2 (в^ вдо + 0
Ждо + 1(0) = Удо + Ъ
Жек(в) Жк-д12 (в)
Жк-д22 (в) ..., Жк-дй2 вк) ,
Жк(0) = Ук.
Очевидно, что проблема алгоритмически решена, т.к. методика без труда может быть реализована в форме программы на входном языке какого-либо ПК А.
Пример 3. Рассмотрим линейное дифференциальное уравнение второго порядка с двумя переменными запаздываниями:
х ''(4) + 2а1х'(4) + 2а2х '(4 - т1(4^ +
+ 71х(4) + 72х^ - т2(4)) = 0,
4 > 0, (20)
х '(4) = 0, х''(4) = 0, 4 € (?о, 0],
х(^о) = 5, х ' (?о) = 0,
где а*, 7* (г = 1,2) - постоянные. Рпс.З демонстрирует поведение х(4) и х '(4), описываемых уравнением (20) с т = 0,125, «1 = 0,05, «2 = -0,04, 71 = 4, 72 = -3 п(г) = |_10вт22^т, т2(4) = [Ювт2^ • т.
Рис.З
4. ДУ с сосредоточенными и распределенными запаздываниями
Рассмотрим систему ДУ (ИДУ) следующего вида:
х(і) = /(х(і), х(і - т), Рт[х(і)],і),
і > і і = іо + Т,
(21)
где основные обозначения те же, что и в разд. 1, а
г
Рт[ж(г)] = У д(ж(0),0) ^0
г-т
(д(-, •) - заданная, непрерывная по всем аргументам функция). Как и выше, предположим, что для 4 € (40,41] вектор ж(4) удовлетворяет следующей системе ОДУ:
х
(і) = / о(x(і),і), х(іо) = уо. (22)
Ясно, что главной причиной, препятствующей применению нашей схемы для уравнений (21), является присутствие векторного функционала Рт[х(і)] в правых частях этих уравнений. Но существует возможность обойти эту проблему. Действительно, мы можем определить новые векторы х^(в) с помощью соотношения
£к + 5
^(в) = У д(хй(0 - ій),0) ^0, в є [0, т].
£к
Это приводит к тому, что для і є Дк, в є [0, т] и 0к = 0 - і к мы получим
£
Рт[х(і)] = у д(х(0),0) ^0 =
£к
У д(хй—і(0й—1),0) й0-
£к-£к — 1+5
^хй—і(0й—і), 0) ^0+
£к —1
+^(в) = «й—і(т) - і(в) + ^(в),
т.е.
(24)
Рт[ж(4)] = ^к(в) + «к-1(т) - «к-1(в). (23)
Теперь мы в рамках применимости нашей схемы и можем начинать просмотр последовательности отрезков Д*.
0°. Сегмент До будет первым. Векторы жо(в) и ^о(в), определенные на этом сегменте, удовлетворяют системе
Ж о(«) = f о( жо(в),во), жо(0) = Уо,
V о(«) = д(жо(в),во) , «о(0)=0.
1 ° Д0 Д1
Ж о(«) = f о( жо(в),во), жо(0) = Уо,
V о(«) = д(жо(в),во) , «о(0)=0,
Ж 1(й) = f (ж1(й), жо(в),
«1(в) + ^о(т) - «о(в),в1), Ж1(0) = у 1,
V 1(в) = ^Ж1(в),в^, «1(0)=0,
°
N°. На этом этапе рассматриваются отрезки Д0, Д1; ..., Дд и система ОДУ для вектора
хо(в) = /о(хо(в), во),
Vо(в) = ^хо(в), во) , хі(в) = /(хі(в), хо(в),
хо(0) = уо,
«о(0) =0,
£к £к + 5
У д(хк—і(0к—1),0) Й0 +1 д(хк(0к),0) Й0 = Мв) + г>о(т) - Ыв),в0, хі(0) = уі
£к
Vі(в) = д(хі(в), ві),
•иі(0) = 0,
£—т
£—т
Рис.4
... ... ... ... ... ... ...
Жд (в) = ^Жд (в), Жд-1 (в), «д (в) + «д -1 (т )-- «д-1(в),вд) , Жд (0) = уд,
« д (в) = д(жд («),«д), «д (0)=0.
И мы снова получили желаемый результат.
Пример 4. Рассмотрим следующую форму уравнения Ван дер Поля:
х//(t) + 2a[x2(t - т) - 1 •х/(t)+ t
+ w2x(t) + y j х(в) de = 0, t> 0,
t-т
х/(t) = 0, x "(t) =0, t Є (-т, 0], х(-т) = 0,02, х/ (-т) = 0,05,
(27)
где а ш, 7 — постоянные величины. Некоторые из результатов вычислений отображены на рис.4 для а = 0,05, ш = 2, 7 = 5, (а) т = 1, (Ь) т = 1,2.
Заключение
Алгоритмы, описанные в данной работе, могут быть эффективно реализованы на основе любого современного ПКА, такого как Мар1е или Ма1;1аЪ, и использованы для изучения многих типов систем с последействием. Кроме того, различные схемы МШРФП, подобные изложенным выше, применялись для исследования поведения фазовых векторов, которые являлись решениями стохастических ДУ со специализированными формами запаздывания и ИДУ.
Список литературы
1. Азбелев Н.В., Максимов В.П., Рахматул-
лина Л. Ф. Введение в теорию функциональнодифференциальных уравнений. М.: Наука, 1991. 280 с.
2. Беллман Р., Кук К. Дифференциальноразностные уравнения. М.: Мир, 1967. 548 с.
3. Хейл Дж. Теория функционально-дифференциальных уравнений. М.: Мир, 1984. 421 с.
4. Эльсгольц Л.Э., Норкин С.Б. Введение в теорию дифференциальных уравнений с отклоняющим аргументом. М.: Наука, 1971. 296 с.
5. Bellen A., Zennaro М. Numerical Methods For Delay Differential Equations. Oxford: University Press, 2003. 416 p.
6. Shampine L.F., Gladwell I., Thompson S. Solving ODEs with Malab. Cambridge: University Press, 2003. 272 p.
7. Полосков II. 1%. Расширение фазового пространства в задачах анализа дифференциальноразностных систем со случайным входом / / Автоматика и телемеханика. 2002. № 9. С.58-73.
8. Полосков И.Е. Компьютерное моделирование динамики загрязнения бассейна реки с учетом запаздывания и случайных факторов / / Вычислительные технологии. 2005. Т. 10, № 1. С.103-115.
9. Poloskov I.E. Symbolic-numeric algorithms for analysis of stochastic systems with different forms of aftereffect // Proc. in Applied Mathematics and Mechanics (PAMM). 2007. Vol.7, Is.l. P.2080011-2080012.
10. Malanin V.V., Poloskov I.E. About some schemes of study for systems with different forms of time aftereffect // Proc. of the IUTAM Symp. on Nonlinear Stochastic Dynamics and Control
(Hangzhou, China) / W.Q.Zhu, Y.K.Lin, G.C.Cai (Eds.): IUTAM Bookseries, Vol. 29. Dordrecht: Springer, 2011. I’.55 01.
11. Wolfram S. The Mathematica Book: 5th ed. Champaign, II: Wolfram Media, 2003. 1488 p.
12. Gopalsamy K., Gydry I., Ladas G. Oscillations of a class of delay equations with continuous and piecewise constant arguments / / Funkcialaj Ekvacioj. 1989. Vol.32. P.395-406.
Symbolic and numeric schemes for study of dynamic systems with aftereffect
I. E. Poloskov
Perm State University, 614990, Perm, Bukirev st.,15 polosk@psu.ru
Different types of dynamic retarded systems are under consideration on the base of our scheme of phase space expansion. This scheme is applied for investigation of linear and nonlinear differential difference equations with single and multiple constant delays, differential equations with variable delays, neutral delay differential equations, special types of integro-differential equations, and some others. Corresponding algorithms were implemented with the help of Mathematica-code programs. Some examples of analysis of different systems are presented.
Keywords: dynamic system,; aftereffect; phase space; extension; phase vector, method; modeling-, numeric integration