CC BY
11МОДЕЛИРОВАНИЕ УСТОЙЧИВОСТИ НА ОСНОВЕ ОТНОШЕНИЯ ВОЗМУЩЕНИЯ РЕШЕНИЯ СИСТЕМЫ ОБЫКНОВЕННЫХ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ К ВОЗМУЩЕНИЮ НАЧАЛЬНЫХ ДАННЫХ
Ромм Я.Е.
Таганрогский институт имени А.П. Чехова (филиал) ФГБОУ ВО "РГЭУ (РИНХ)", д-р техн. наук, профессор
SIMULATION OF THE STABILITY BASED ON THE RATIO OF THE PERTURBATION SOLUTIONS OF SYSTEMS OF ORDINARY DIFFERENTIAL EQUATIONS TO THE PERTURBATION OF THE INITIAL DATA
Romm Ya.E.
A.P. Chekhov Taganrog Institute (branch) Rostov State University of Economics (RSUE), Doctor of Engineering Science, Professor
АННОТАЦИЯ
Строятся численные признаки устойчивости по Ляпунову решений систем обыкновенных дифференциальных уравнений. Их основой являются необходимые и достаточные условия устойчивости, выраженные в асимптотических границах отношения возмущения решения к варьируемому возмущению начальных значений. Численное моделирование реализуется дискретизацией отношения по ходу решения системы. Границы роста выходной последовательности определяют признаки устойчивости в сравнительно общих условиях. Описаны эксперименты для линейных и нелинейных систем с программным формированием признаков при помощи разностного решения по методу Эйлера.
ABSTRACT
We construct the numerical signs of the Lyapunov stability of solutions of systems of ordinary differential equations. Their basis is the necessary and sufficient conditions for stability, expressed in asymptotic boundaries of the ratio perturbation solutions to the perturbation of variable initial values. Numerical modeling is implemented by sampling the ratio in the process of solving the system. The boundaries of increase of the output sequence determine the stability in relatively general terms. Experiments for linear and nonlinear systems with formation of signs using a finite difference solution for the Euler method.
Ключевые слова: обыкновенные дифференциальные уравнения, устойчивость по Ляпунову, численное моделирование.
Keywords: ordinary differential equations, Lyapunov stability, numerical simulation.
Постановка вопроса. Анализ устойчивости по Ляпунову (ниже устойчивости) решений обыкновенных дифференциальных уравнений (ОДУ) рассматривается в аспекте численного моделирования. Для приложений устойчивость можно моделировать с помощью приближенных решений ОДУ, применяя, в частности, рекуррентные [1 - 3] и аддитивные преобразования [4] разностных методов. Помимо необходимости учитывать погрешность приближений, в случае нелинейных систем возникает ограничение, связанное с построением таких преобразований, - разность возмущенного и невозмущенного решений не должна обращаться в ноль на полуоси, это сужает область приложений результатов анализа. Целью сообщения является увеличить общность метода, сняв данное ограничение, и моделировать устойчивость без преобразования приближенных решений. С этой целью оценивается рост отношения возмущения решения к варьируемому возмущению начальных значений, устойчивость определяется равномерной ограниченностью
отношения на полуоси в окрестности невозмущенных начальных данных. В компьютерной реализации числовая граница дискретно выводится по ходу решения. Это проще моделирования равномерной непрерывности решения на полуоси и непрерывной зависимости от начальных данных на основе определения устойчивости. Ниже излагаются необходимые и достаточные условия устойчивости, численные признаки на их основе, ориентированные на компьютерную реализацию, и результаты численного эксперимента.
Исходные предположения и предварительные преобразования. Рассматривается задача Коши для системы
dY = F(t,Y), Y(О=У0, (1) dt
где F(tJ) = (f,(t,Y), f2(t,Y),..., f (t,Y)), Y = (y,(f), y2(t),..., y„ (t)), Y0 = (y01, y02,..., y0„). Предполагается, что существует 50 > 0, такое, что
на полупрямой [i0, да) выполнены все условия су- нию I Y(t) || = max | ук (t) | , в численных эксперимен-ществования и единственности для невозмущен- 1iki"
ного решения и для каждого его возмущения тах используется эвклвдова норма. Предполага-
У = У (/), У(¿0) =У0, с начальным вектором в границах || У0 — У11 < 80. Применяются канонические согласованные нормы матрицы и вектора. По умолча-
ется,
что
области
R0: { t0 < t <да; Y (t),V7(t) : || f0 - Y0\\ < 50 } функция
F(t,Y) определена, непрерывна и непрерывно дифференцируема по t (в i0 - справа), компоненты функции удовлетворяют неравенству
|/,' (t,Y)- fk(t,Y)|<L\yk -~ L = const Vt e [tB,да), V(t, Y),(t,f )e R0,Vke \n . (2)
Из (2) следует условие Липшица:
||F(x,Y)-F(x,f)||<Z|| Y-Y I, ~ = const, V(x,Y), (x:
. В этих условиях определение устойчивости [5] упрощено: решение Y = Y (t) устойчиво, если Ve> 0 найдется А, 0<А<50, такое, что
J ~ - Y01 < А влечет || ~(t) -Y(t) || < e Vt e [t0,да). Решение асимптотически устойчиво, если оно устойчиво и найдется Д, 0<Д <Д, такое, что
I ~ - Y01| <А0 влечет Y(t) -Y(t)^0 при t^да.
t=const, t=ti+1, h=
(3)
Метод Эйлера решения задачи (1)
У)е Я0 У+1 = У + РЬ У)к,
включая запись с остаточным членом, на произвольном отрезке [/„,/] рассматривается в предположении, что значение t е р0, да) остается фиксированным, при этом индекс I неограниченно возрастает с убывающим равномерным шагом:
t+1 ~ to " i + 1
i=0,1, ..., tt+1 =tt + h, 0 < k < i.
(4)
В форме с остаточным членом метод (3) можно записать в виде
У+1 =у + У )к+0, а = (^,^,...,Чп1) , (5) где 1, к из (4), - остаточный член формулы Тейлора для к -й компоненты приближения: лагается, что
волной: ~i = ~ + F(ti,Y )h+Q, Q = (~ii,~,-,qj,
~k= 1 /k'(Iki,~(~ki))h2,ti <ti+1.
На время предварительного описания предпо-
1
4ki = 2 Жki,Yß,))h2, ti <^и <ti+1, ke 1,n .
Для
~(t)—л(t)*0 Vt е [t0,да), Уке1,и . (6)
Рассматриваются следующие преобразования возмущенного решения аналоги (5) отмечаются возмущения:
У,
k (i +1) J k (i +1)
у, (i+1) = ~ki - ytl +
L (ti ,Yi)-/k (t, ,Yi)
~ki- y,i
(~~ki- y,i) h+W
или
~ /"1 i n(k , . \ I n(k ) — fk (ti ,Yi) fk (ti ,Yi)
Ук(i+1)-Ук(i+1) = (1+D 'h)(yw -y,i)+wti-, ^i =---
y,i- Ук>
wki = 4ki - 4ki .
(7)
(8)
С учетом (2) при ограничении (6) выполняется неравенство
|/k (t,Y) - /k (t,f)'
—рт-^— <L, L = const Vt e[t0,да), V(t,Y),(t,J~)e , Vke1,n . (9)
у, (t)-(t)|
Отсюда
Из (8)
У
D^Z < L, L = const, Vi=0,1, 0,i, Vke1,n .
к (i+1) - Ук (i+1) =П (1+ D^ h )(~k 0 - Ук 0) + R 0k', R 0k' = Z П (1+ D™ h ) Wk (i-,-1) + w,
к(i-r-1) 1 ,yki '
(10)
(11)
где h из (4), w из (8), к e1, n. Для дальней-
Лемма 1. В рассматриваемых условиях
шего потребуются приводимые непосредственно I'm R0i =0 Vte[to,w),Vke1,n .
ниже утверждения и элементы доказательств из [2
Доказательство. В силу непрерывной диффе-
- 4]. Пока не оговорено иное, будут предполагаться , „
J ^ ■> г ^ ренцируемости Ри,У) в области ^ для
выполненными все рассматриваемые условия, 0
включая (6). ^ е[to,да)лV~е[^t] ЗСо = Со^) \
r
£ = 0
r = 0 £=0
\FXt ,Y)||<Q, \\F;(t ,Y)||<C0.
Поэтому \R(k) I <O.
^|<C0h2, г = 0,1,...,i; i=1,2,...; ¿=1,2,...,л . Из
(11) с учетом (10) \R (0*? |<2 (1+ Lh) 'C0h2, или, ' i ' i=0
Л!
|R(,*? |<С ((1+ Lh)i+1 -l)h . Отсюда и из (4)
L
t-t 0 h
R ^ < C^(eL(t-t0)-1)
(1+ Lh) h -1 | h . t-L
Следовательно,
i+1
и R(k) ^ 0 при i ^ да
V? е[?0,да), Уке 1,и . Лемма доказана. Согласно (11) имеет место Следствие 1. В тех же условиях
Ук (t) - y (t) = lim П (1+D™ h)(y 0 - y о) V t e[ to, да), V k e1, n.
(12)
Из доказательства леммы ясно, что на любом Лемма 2. В рассматриваемых условиях для
отрезке [?0, ?] сходимость (12) является равномер- устойчивости решения задачи (1) необходимо и до-
ной, условия леммы обеспечивают равномерную статочно, чтобы существовало 4,0 <4 < 5о, та-
сходимость метода Эйлера. кое, что для всех возмущенных решений
Из (12) и определения устойчивости вытекает
lim П (1 + А-к>)
Y (t), Y (t0) =Y , при ограничении верно неравенство
<C(1), C(1) = const, V?e[?0,да), Vкe!~n.
Y-Y <А
-*0 -*0 —
(13)
Для асимптотической устойчивости необходимо и достаточно, чтобы выполнялось предыдущее утверждение и существовало 4<4, такое,
На [t0, t] частичное произведение П (1 + Dt™h)
что Y - Y <А влечет
lim
lim П (1+D-h
при изменении i меняет одновременно все сомножители, при этом h в каждом из них изменяется согласно (4).
Выражение условий устойчивости через отношение возмущения решения к возмущению начальных значений. Из (12)
У* (t) - У* (t) = lim p (1+d (*h) Vte[t0, да), V *е1П . (15)
У* о - У* о i "да '=0
=0 Vke 1,n . (14)
Отсюда и из леммы 2 вытекает
Следствие 2. При ф ук0 формулировка и утверждение леммы 2 дословно сохраняются при замене (13) на соотношение
Ук (t) - Ук (0
Ук 0- Ук,
<C(1), C(1) = const, Vte[t0,да), Vke 1,n ,
(16)
и (14) - на соотношение
У* (t) - У* (t)
lim
t ^ да
= 0 V к e 1, n .
_ (17)
Ук 0 Ук 0
Соотношение (16) означает равномерную ограниченность дроби на полуоси в некоторой окрестности возмущения начальных значений. Дробь выражает отношение возмущения решения именно к вызвавшему его возмущению начальных значений, при всех вариациях которого в границах
| ук0 "Ук„|<41 выполняется
^к ? ~ ^к(0 = О (1) УГе[Г0, да), V к еЩ.
Ук 0 _ Ук 0
Расширенные условия устойчивости. Выше предполагалось выполненным (6). Требуется показать, что обе леммы и следствия сохраняются при нарушении этого ограничения. В дальнейшем предполагаются все исходные условия, за исключением (6).
Лемма 3. Пусть для невозмущенного и возмущенного решений Y(?), 7 (?) задачи (1) рассматри-
ваются * е 1, п, L из (2), t из (4) и любые компоненты у (t), y (t), включая случай tt ^ t , где У к ( t) = У к (t), причем в некоторой окрестности t выполнено У1с (t) Ф Ук (t) Vt Ф t . В этих условиях существует конечный предел
lim (ti ,Y) - (ti ,Y) = В. , или, h У* (ti) - У* (ti) *
lim f* (t-Y - f*(= ß , Iß I < l , V *е 1л . (18) - У* (t) - У* (t) P* , |P*' , , ( )
Доказательство. Согласно (4) оба соотношения эквивалентны, рассматривается (18). Если знаменатель под знаком предела не обращается в нуль, то с учетом (9) конечный предел существует по непрерывности числителя и знаменателя, и выполняется (18). Если знаменатель равен нулю в точке t , то в достаточно малой ее окрестности выполняется (2), где по условию Ук (t) ФУк(t) Vt Ф t , следовательно, выполняется (9). По условиям существования и единственности в области R 0 функции
F(x,Y), F(x,Y) и решения Y(t), Y(t) задачи (1)
е=0
t=0
е=0
t ^ да
t = 0
определены и непрерывны на [t0, да). Поэтому в (9) существуют пределы: lim (У (t) - у (t) )= cojt и
lim (f, (t,Y) - f, (t,Y))= c1k, а разности под знаками
пределов достигают конечных наибольших значений на любом отрезке [t0, T], включающем t вместе с рассматриваемой окрестностью. Отсюда сами пределы конечны: kj < да и вследствие (2)
й r f (t,Y) - ft (t,Y) Cli
предел дроби: lim —--= , при этом
Р Р tm У (t) -У, (t) coi ' P
< L. Для любой подпоследовательности
t ^t, t ^ t , / ^ да , существует тот же предел
ßk = ^ , |ßk| < L.
kJ < L\cn¡
В рассматриваемой окрестности
точки ~ , где у ^) * у (^ Vt * t , отношение непрерывных функций ограничено константой Ь, при t ^ t существуют конечные пределы числителя и знаменателя, следовательно, существует
lim f t ,Y ) - f ( t Y ) = „ '■m y (t, ) - y (t ) Pk '
Лемма доказана.
Лемма 4. Условия и утверждение леммы 3 распространяются на D ® из (11) V k е 1, n , We 1,i , i = 1,2, ... - в условиях этой леммы всегда существует конечный предел
lim fk (ti-< , Y-i ) fk ( ti-i , Y-< ) = ß
ßj < L, V k el, n .
(19)
Доказательство. Соотношение (19) эквива-
где h и i
лентно lim fk (fi-t ,Y-<) fk (ti-< ,Y-<) = ß .
*m уЛи-уЛи
из (4), как в (18). Для (19) с точностью до обозначений сохраняются рассуждения и оценки доказательства леммы 3, согласно которому конечный предел существует для любой подпоследовательности t, , ^ 7, t, 7, i ^ да :
]пт~ Л ,у) ла,у) = р, |р>к | < Ь . Лемма до-
^ уЛ)—УЛО 1 1
казана.
В условиях леммы 3 исключены неустранимые особенности в Дк) и Д® .
Следствие 3. Соотношение (12) будет сохраняться на произвольном отрезке [^, ^ , если условия леммы 1 расширить до условий леммы 3 и кроме того допускать любое конечное число точек, в которых на этом отрезке нарушается ограничение (6).
Доказательство. При однократном нарушении (6) ни одно из соотношений доказательства леммы 1 в расширенных условиях не изменится в силу того, что сохраняются (18) и (19), поэтому, если некоторое ^ е [^, ^ оказывается нулевой предельной точкой знаменателя в (18) или (19), соотношение (12) будет сохраняться. Если на [^, ^ количество точек £ , где у (^) = у (£), конечно, то
они взаимно отделены. В достаточно малой окрестности каждой из них нарушение ограничения (6) однократно, следовательно, как и в остальных точках [^, ^, сохраняется (12). Следствие доказано.
Из следствия 3, соотношения (12) и леммы 2 вытекает
Теорема 1. В исходных предположениях, с тем изменением, что на произвольном отрезке [^, ^
допускается любое конечное число нарушений ограничения (6), для устойчивости, а также для асимптотической устойчивости решения задачи (1)
Ук Й—<)—Ук (ti—<)
сохраняются необходимые и достаточные условия леммы 2, включая соотношения (13) и (14). Согласно (15) имеет место Следствие 4. При у0 *у0 Vkе1,п утверждения следствия 2, включая соотношения (16), (17), сохраняются в условиях теоремы 1.
Условия следствия 4 можно расширить. Пусть решение системы (1) устойчиво, произвольно зафиксировано к, у0 * у0, , да) взят отрезок Р0, t], пусть выполнены исходные условия, за исключением произвольного множества точек t е [^, ^, в которых нарушается (6): у (t) = у ( t). Согласно следствию 1 для разности у (/)—у ) в каждой точке ^ е ^^, t], взятой вне данного множества, сохраняется (12). В силу устойчивости в каждой такой точке с необходимостью выполняется (13) и (16), где в Д1 -окрестности невозмущенных
начальных значений константа С (1) не меняется. В точках t е [^, t], в которых у ( t) = у ^ ), функция
Ук а)—Ук а)
< с(1), с(1) = const, Vt e
V {t ). В
ре-
равна нулю. По непрерывности этой
Ук о — Ук о
функции каждое t можно заключить в столь малый интервал V ^ ), что
Ук (К)—Ук (К)
ук о — Ук,
зультате все точки отрезка р0, ^ можно покрыть системой интервалов и V^, ^), на каждом из которых выполнено (16), где С (1) не меняется. По лемме Бореля из системы |JV^, £) можно выбрать конечную подсистему, покрывающую р0, t]. На каждом интервале конечного открытого покрытия, следовательно, и во всех точках р0, t], для выбранных к и t выполнено (16). Рассуждения воспроизводятся для V к е 1, п с одним и тем же значением
C
(1)
c
C
С учетом произвольности выбора ? отсюда следует, что в случае устойчивости соотношение (16) с необходимостью выполнено
Vkе\п , У?е[?0,да).
Вывод опирается на (12), что априори предполагает непрерывную дифференцируемость Г (г,У) из (1).
С другой стороны, выполнение (16) является достаточным условием устойчивости без этого предположения. В самом деле, из (16)
|Ук (О-Ук (?)|< С 1 Ук 0 - Ук„| , С(1) = ССП81,0 < С(1)
, где ук0 - у0 ф 0. Отсюда Уе > 0 будет выполнено | У (?)-У (?) < е Vке 1,п, V?е[?0,да), лишь только
< C, C = const, Vte[t0,да) . (22)
C (1)| Ук 0 - Ук,
1< 8. I ~ (О-Y (о I <А
Пусть
А1 =-
Если
то
С(1)
||Г (?) - У (?)| < е У?е[ ?0, да).
По определению выполнение (16) - достаточное условие устойчивости без отмеченного ограничения.
Теорема 2. В исходных условиях, включая непрерывную дифференцируемость Г(?,У ) из (1), с тем изменением, что на произвольном отрезке [?0, ?], Vке 1,п , допускается любое множество точек ? , в которых у (?) = у ( ?), для устойчивости решения задачи (1) необходимо, чтобы существовало 4,0 < 4 < ^о, такое, что для всех возмущенных решений У (?), У (?0) =У0
7 -У0 <41 > ук 0 Ф Ук 0 У к е1, п
условия ее устойчивости в случае Y(t)=0 Vt е [t0,да), Y =0, при этом (20) переходит в соотношение
У* (t)
У* 0
Необходимые и достаточные условия асимптотической устойчивости точки покоя получаются из этой теоремы при переходе (21) в соотношение
- lim i М)
Vte[?0,да), V*е 1,п ™ | у*0
Численное моделирование (20) - (23) в качестве признаков устойчивости выполнимо по ходу приближенного решения задачи (1). Ниже будут даны программные примеры.
Из теоремы 2 и неравенства (2) вытекает Следствие 6. В условиях теоремы 2 для устойчивости решения задачи (1) необходимо, чтобы су-
= 0. (23)
при ограничении верно неравенство
Ук (t) - Ук (t)
<C(l), C(1) = const, Vte[t0,да) . (20)
yt 0 - yt о '
Эти же условия являются достаточными условиями устойчивости без предположения о непрерывной дифференцируемости правой части (1).
Для асимптотической устойчивости необходимо и достаточно, чтобы в соответственных условиях устойчивости существовало Д2 , такое, что
Y -Y
-*0 -*0
<А влечет:
lim
t да
= 0 .
(21)
| Ук (t) - У* (t)
[ У* 0 - У* 0 J *=1
Доказательство. Условия устойчивости при выполнении (20) обоснованы выше. Остается рассмотреть (21). Поскольку Уы -уы ф0 V*e 1,п , (21) выполняется тогда и только тогда, когда lim (у (t) - У (t)) = 0 V*е1, п , или,
Ьп| | { Ук(t) - Ук(t) L|| = 0. Иными нение (21) при
словами, выпол-
Y0 -Y0 <А2 эквивалентно опреде-
лению асимптотической устойчивости. Теорема доказана.
Следствие 5. Пусть система (1) имеет точку покоя. Теорема 2 дает необходимые и достаточные
ществовало
А,
0 < А< 50,
такое,
что
VY(?), Y(?0)=Y0, при ограничении 0 < Y -Ц <А7 верно неравенство
fk(t,~) - fk(t,Y)
Ук 0 - Ук 0
< Cf, Cf = const, Vte[t0,да)
. Для асимптотической устойчивости необходимо, чтобы решение было устойчиво и существовало
А f <А f =
такое,
что
lim
t ^ да
fk (t,Y) - fk (t,Y)
Ук 0 - Ук 0
Y-УЛ <Аf
=0 .
влечет
Теорема 2 и следствия применимы к нелиней ным и линейным системам.
Случай линейной системы. Можно ограничиться однородной системой
ау
— = A(t)Y, Y(t0)=Y0 dt
(24)
где элементы матрицы А (?), п х п, функции а^ = (?) V/, j'е1,п определены, непрерывны и, если это оговорено отдельно, непрерывно дифференцируемы V?е[?0, да).
Следствие 7. Следствие 5 дает необходимые и достаточные условия устойчивости и асимптотической устойчивости системы (24).
Доказательство. Все решения системы (24) одновременно устойчивы или неустойчивы [6], можно ограничиться точкой покоя, для которой следствие 5 дает необходимые и достаточные условия устойчивости. Все решения линейной системы асимптотически устойчивы тогда и только тогда, когда асимптотически устойчива точка покоя [6], необходимые и достаточные условия для этого сформулированы в этом же следствии 5. Следствие доказано.
Следствие 8. Для устойчивости системы (24) необходимо и достаточно выполнение (22) для ка-
8
к=1
n
к=1
кого-либо одного решения Y(t), Y(t0) =Y0, независимо от выбора ^, Уы ^0, е1,п . Аналогично, для асимптотической устойчивости необходимо и достаточно выполнение (23) для какого-либо одного решения Y(t).
В случае линейной системы выбор Д -окрестности излишен [6]. Поэтому в следствиях 7, 8 речь может идти лишь о выполнении соотношений (20), (21) и (22), (23).
Численный эксперимент. Используются эвклидовы нормы матрицы и вектора. I. В системе
dyi гг;—Г
d" = - y Wy.+ y 2 =
dy222
—^Г = ywy 1 + y 2 нУле-
вое решение устойчиво, все остальные неустойчивы [5]. С помощью метода Эйлера устойчивость нулевого решения иллюстрирует программа на основе (20), (22):
program Stability;{$APPTYPE CONSOLE} uses SysUtils; const h = 0.000001; tt=10000000;
var X,y1,y2,y10,y20,z1,z2,eps1,eps2: extended; k: longint; function f1 (x,y 1 ,y2 :extended): extended; begin f1:=-y2*sqrt(sqr(y1)+sqr(y2)); end; function f2(x,y 1 ,y2 :extended): extended; begin f2:=y1* sqrt(sqr(y 1 )+sqr(y2)); end;
begin eps1:=0.000003; eps2:=0.000003; y1:=0.0; y10:= y1+eps1; y2:=0.0; y20:= y2+eps2; k:=0; x:=0; while x <=1500 do begin y1:= y1 + h * f1(x,y1,y2); y2:= y2 + h * f2(x,y1,y2); y10:= y10 + h * f1(x,y10,y20); y20:= y20 + h * f2(x,y10,y20); z1 :=abs(y 1-y10)/eps 1; z2:=abs(y2-y20)/eps2; k:=k+1; if k = tt then begin writeln('x=',x:4,' ');
writeln('z1=',z1:4,' ','z2=',z2:4,' ','norma=',sqrt(sqr(z1)+sqr(z2))); k:=0 end; x:=x+h; end; readln end.
Признак устойчивости нулевого решения получается из значений (20) при у (^ = 0, у2 (0 = 0 , у (0) = 0, у2 (0) = 0, выводимых на отрезке [0,1500] :
IIг11=^/г2 + г2 : 1.4, 1.4, ... , 1.4.
Возмущение начальных значений eps1, eps2 можно варьировать без изменения границ || г ||: в качестве правой части (20) или (22) можно взять С(1) = 1.42 при 0 <||71| <Д,, где 4 = 0.000003 . При сужении Д -окрестности, например, eps1 = 4 /10, eps2 = Д /10 или
eps1 = Д /100, eps2 = Д /100 сохраняются те же
значения г , аналогично, в случае eps1 = Д /3000, eps2 = Д /3000 .
Признак неустойчивости нетривиальных решений программа даст при начальных значениях у1:= 0.0001; у2:= 0.0001 . На отрезке [0,1500] появится незначительный рост нормы
Щ: 1.41, 1.41, ... , 1.42, ..., 1.43, ..., 1.44.
Значения || г || сохраняются при увеличении 4 в < 3000 раз и при его уменьшении в обратной пропорции. Существенное влияние на рост нормы оказывает рост начальных значений. При у1 := 0.001; у2 := 0.001 получится
||Z||: 1.4 , 1.5 , ... , 1.6 , ... , 1.9 , ... , 2 , ..., 2.1 ,..., 3 , ..., 3.3. При y1 := 0.01; y2 := 0.01 рост нормы усилится -|Z||: 1.4 , 1.5 , 1.6 , ..., 2 , ..., 3 , ..., 4 , ..., 15 ,..., 16 ,..., 20 ,..., 24 ,..., 29 , ..., 30. При y1:= 0.1; y2:= 0.1 получится ||Z||: 2 , ... , 4 , ... , 12 , ... , 89 , ... , 91 , ... , 180 , ... , 251 , ... , 301 , ..., 312.
Признаки неустойчивости приведены для возмущений eps1:= 0.000001/3; eps2:= 0.000001/3 .
dv, „ 1 2
ddf=р^"+7 vkv2
Z at' Xbrf
II. В [4] показано, что для системы
^=U(t,V), r=1,2,
.., n, r Ф к,
(25)
f=0
i=0
V (t0)=0, ke1,n, где 0 < p,0<, и 0 < q, 0<b, i=0,1,...,q, r -е уравнения могут быть
произвольными в условиях существования и единственности, точка покоя неустойчива при значениях ^ = 0.5, а = -0.5, р=0.5; р = q = 3;
a = b =1; n = 2; u2 (t,v) =-1. Это подтверждает программа Stability при соответственных изменениях:
function f1(x,y 1,y2:extended): extended;
begin f1 :=1/2*y1/sqrt(x)+(1/sqr(x))*y1*sqr(y2)*(1+x+sqr(x)+x*sqr(x))*
(1+sqr(y2)+sqr(sqr(y2))+sqr(y2)*sqr(sqr(y2)))end;
function f2(x,y 1 ,y2 :extended): extended;
begin f2:=-x end;
Рост нормы можно видеть уже на [0.5, 2.5] (x:=0.5; while x <= 2.5 do) при tt = 100000 -||Z|| =Jz2 + z22 : 1.46 , 1.51 , ... , 1.44X101 , 1.31x102 , ... , 1.66x104 , ... , 5.60x1095 ,
Далее наступит переполнение. Ввиду нарушения (22) - это признак неустойчивости. Аналогично предыдущему, он сохраняется при вариациях возмущений начальных значений.
III. Точка покоя системы (25) неустойчива при Va >-1 [4]. Если в той же программе положить а = -1 (f1:=1/2*y1/x+...), образуется переполнение. При tt = 100000 на отрезке [0.5, 2.5] виден рост нормы
||Z||: 1.00 , 1.00 , ... , 1.002 , 1.00079 , ... , 9.9973х10-1
1.31х101853 . ... , 4.83х10207 ,
II г II: 5.23 , 1.49Х101 , 1.34х102 , . .
... , 1.22х101853.
IV. Точка покоя системы (25) устойчива, если а < -1, а перед вторым слагаемым правой части вместо «+» взят «-» [4]. Если в рассматриваемой программе взять а = -1.5 , а перед вторым слагаемым «-» (й:=1/2*у1Лэд11(х)*х)-(1/sqr(x))*y1*sqr(y2)*...), то
9.976х10":
9.89х10-1
9.89х10-:
Эти значения получены для [4]: а > -1, перед первым и вторым слагаемым пер-h = 0.000001; й =10000000; eps1:= 0.00003; ефй := атО^о уравнения в правой части знак «-», кроме того,
. Результат не изменится, если eps1, eps2 умножить или разделить на у, 1 < у < 1000 . В качестве константы (22) можно взять С = 1.003 при
0 < 11 ~ < 4, где 4 = 0.000003 .
V. Точка покоя системы (25) будет асимптотически устойчива при следующих видоизменениях
II г II: 6.93Х10"1 , ... , 1.08Х10-5 , ... , 1.55Х10-7, ... , 1.31х10"10
ы2 (/, V)=- уггв , а < р . Если в рассматриваемой программе снова взять а = -0.5 и изменить знак первого слагаемого на противоположный (Л:=-1/2*у1^1!(х)-(...), соответственно изменить правую часть второго уравнения (£2:=-у2*1^г1(х);), то значения нормы примут вид
, 7.15х10-14 ,
, 3.56х10-17 ,
, 3.09х10-1
Согласно (23) это признак асимптотической устойчивости. В качестве 4 < 4 для (23) можно взять 4 < 0.0000029, 4 = 0.000003 .
VI. Случай линейной системы. Пусть в (24) в качестве Л(() взяты матрицы
( -1 - 0.09 0.077 - 0.01 ^ Г 1 -0.09 0.077 - 0.01"
0.087 - 0.9 0.005 0.019 А — 0.087 0.9 0.005 0.019
- 0.034 0.034 -0.2 - 0.06 , A1 - - 0.034 0.034 0.2 - 0.06
V 0 - 0.022 0.092 -14 , , 0 - 0.022 0.092 14 ,
A -
Если в программе Stability удвоить размерность, дописать константы, переменные и функции,
(26)
и выводить норму отношения возмущения к возмущению начальных значений для четырех компонентов, то для матрицы Л получится
V 1-1
z -Л Iz2 : 1.02х100,_, 5Л6х10-17,..., 3.72х10-29,..., 8.07х10-86,..., 2.94х10-1
1.20х10-135.
Согласно (23) это признак асимптотической В матрице Л1 знак диагональных элементов
устойчивости, что соответствует диагональному положителен, рост1 нормы соответствует неустой-
преобладанию отрицательных значений. чивости -
||г II: 5.75x10°,..., 1.01Х1020,..., 3.68х10100,..., 1.85х10491,..., 4.30х10656, ..., 8.62х10772,..., 4.18х10908.
Признак не изменится при вариациях окрестности начальных значений, но для линейной си-
стемы это излишне, поскольку все определяет ограниченность, неограниченность или стремление к нулю любого отдельно взятого решения [6].
VII. Чтобы оценить степень достоверности численных признаков, целесообразно рассмотреть пограничные случаи. Для этого можно взять матрицу Фробениуса, строка которой соответствует характеристическому полиному с кратными (или
A 2 =
В них последняя строка - коэффициенты характеристических полиномов, за вычетом коэффициентов перед старшей степенью, взятые с обратным знаком. Для А 2 характеристический полином
напротив) комплексными корнями с нулевой действительной частью. Пусть в (24) в качестве A(t) взяты матрицы
f 0 1 0 0 0 0 1 f 0 1 0 0 0 0 1
0 0 1 0 0 0 0 0 1 0 0 0
0 0 0 1 0 0 0 0 0 1 0 0
, A =
0 0 0 0 1 0 0 0 0 0 1 0
0 0 0 0 0 1 0 0 0 0 0 1
,4 0 7 0 2 0 ; -4 -16 - 25 - 20 -10 - 4 ;
(27)
или,
имеет
кратные
мнимые
корни:
function f1 (x,y 1 ,y2,y3 ,y4,y5,y6:extended):extended; begin f1:=0*y1+1 *y2+0 *y3 +0 *y4+0 *y5+0 *y6 ; end; function f2(x,y 1 ,y2,y3 ,y4,y5,y6:extended):extended; begin f2:=0*y1+0 *y2+1*y3+0 *y4+0 *y5+0 *y6 ; end;
р (х) = хб - 2х4 - 7х2 - 4 ,
Р (х) = (х2 + 1)2(х2 - 4). В этом случае система не устойчива [7].
При соответственном увеличении размерно сти, задании переменных и функций
begin f5 :=0*y1 +0 *y2+0 *y3 +0 *y4+0 *y5+1*y6; end; function f6(x,y 1 ,y2,y3 ,y4,y5,y6:extended):extended; begin f6:= 4*y1+0*y2+7*y3+0*y4+2*y5+0*y6; end; модифицированная программа Stability выдаст значения
V t=1
z\= J2: 3.35x10',..., 3.37X1067,..., 3.71x10'46,..., 1.76x1024',..., 1.04x101102 ,..., 4.58x101302.
Вследствие нарушения (22) это численный признак неустойчивости.
Для А характеристический полином
Р (х) = хб + 4х5 +10 х4 + 20 х3 + 25 х2 +16 х + 4 , или,
Р (х) = (х +1)4 (х2 + 4), среди мнимых корней не имеет кратных, а действительные части остальных корней отрицательны, что соответствует неасимптотической устойчивости [7]. При изменениях
function f6(x,y 1 ,y2,y3 ,y4,y5,y6:extended):extended; begin f6:= -4*y1-16*y2-25*y3-20*y4-10*y5-4*y6 end;
модифицированная программа Stability выдаст значения
||Z|| : 1.25x10',...,2.25x10' ,..., 2.36x10' ,..., 2.21x10' ,..., 1.79x10' ,
1.55x10
, 2.31x10
, 2.22x10
Согласно (22) это признак неасимптотической коэффициентов задается в разделе констант про-
устойчивости. граммы [1, 4]. Его применение для (26), (27) дает
VIII. Программировать функции строк мат- результаты аналогичные приведенным выше.
рицы не вполне удобно при большом п. Можно Именно, в [2, 3] показано, что система (24) устой-
применить способ анализа, при котором матрица чива тогда и только тогда, когда
lim П (E + hA(tt _, ))
для асимптотической устойчивости необходимо и достаточно, чтобы выполнялось (28) и
< C , C = const, Vt e[t0 ,œ),
lim
lim П (E + hA(ti_t ))
= 0
Эквивалентная запись (28), (29):
lim П (E+hA(t, ))
< C, lim
lim П (E + hA(t, _, ))
= 0.
Для постоянной матрицы (28) и (29) примут
вид:
lim (E+hA)
< C, lim
0 t ^œ
lim (E+hA)
= 0.
(28) (29)
В этих соотношениях матрицу (E+hA)1 можно заменить на (E+hA )2t [1]:
t = 0
t ^ œ
t = 0
t ^ œ
t = 0
t = 0
lim (E+HA)1
< C, C - const;
lim
t ^ ад
lim (E+HA)2
- 0.
(30)
Следующая программа приближенно реали- program FrobeniusStability;{$APPTYPE
зует численную проверку (30) с равномерным ша- CONSOLE} uses SysUtils; гом:
label 60,111;
const n=6; h=1.1e-14; A1: array [1..n,1..n] of extended= ((0, 1, 0, 0, 0, 0), (0, 0, 1, 0, 0, 0), (0, 0, 0, 1, 0, 0), (0, 0, 0, 0, 1, 0), (0, 0, 0, 0, 0, 1), (-4, -16, -25, -20, -10, -4)); type matr=array[1..n,1..n] of extended; var A,C: matr; s,s0,x: extended; i,j,l,k,k0: integer; begin FOR i:= 1 TO n do for j := 1 to n do begin A[i,j]:= A1[i,j]*h; if i=j then A[i,j]:= A[i,j]+1 end; k0 :=0; 60: for I:= 1 TO n do for j:= 1 to n do begin s:= 0; for l:= 1 to n do s:= s + A[i,l]*A[lj]; C[i,j]:= s end; k:= k+1; x:= h*exp((k+1)*ln(2)); if abs(x) >= 1e5 then goto 111; for i:= 1 to n do for j:= 1 to n do A[ij]:= C[i,j]; s0:= 0; for i:= 1 to n do for j:= 1 to n do s0:= s0 + sqr(A[i,j]); s0:= sqrt(s0); if k0 >= 1 then write ('norma-: 2, s0: 2 ,' '); k0:=k0+1; goto 60; 111: writeln; writeln(' ':2,' ',h,' '); writeln; writeln(' ':2,' ','k=',k,' '); writeln(' ':2, 'x=',x:2,' '); readln end.
Непосредственно в данной записи моделиру- Если Л заменить матрицей Л из (26), то про-ется устойчивость системы с матрицей Л3 из (27). грамма выдаст
(E+HA)2
: 2.0х100
2.0х100
7.8х10-2, ... , 1.5х10-18, ... , 5.3х10-72, ..., 9.1х10-2282.
Согласно (30) это признак асимптотической устойчивости, как и результат предыдущей программы на основе (23).
(E+HA)2
: 2.0х100
1.0х10],
3.2х107
Если A заменить матрицей A1 из (26), программа FrobeniusStability выдаст переполнение. Рост нормы виден на отрезке [0,1000] (IF abs(x) >= 1e3 THEN goto 111;). , 1.0x1015, ... , 1.0x1060, ..., 1.2x10240.
Согласно (30) это признак неустойчивости, как и результат на основе (22).
Если перейти к матрице Л2 из (27), то программа выдаст переполнение, на отрезке [0,1000] можно видеть значения
(E+HA)2
: 2.4х100,..., 4.5х100,., 4.6х10\ ., 4.8х105, ., 6.6х1021, ., 2.3х1086, ..., 3.5х10344.
Согласно (30) это признак неустойчивости, как и результат на основе (22).
Если вернуться к матрице A3 из (27), та же программа на [0,1000] даст значения
(E+HA)2
: 2.4х100, ... , 7.3х100,
1.0х10
1.2х10
, 9.4х10
, 6.7х100,..., 1.1х10!.....9.6х100.
Согласно (30) это соответствует устойчивости, как и результат программы Stability на основе (22).
При интерпретации результатов численного моделирования учитывается, что в исходных ограничениях, при условии устойчивости, разностные методы накапливают погрешность не более чем линейно по длине промежутка (для метода Эйлера доказательство дано в [8, 9], для метода Эйлера-Коши - в [10]). Эксперименты по выбору шага показывают, что данные численные признаки устойчивости для линейных систем сохраняют достоверность в диапазоне 10 15 < h < 10 3 [11 - 13], для нелинейных систем - в диапазоне 10-6 < h < 10 3 [4, 14].
Промежуток, на котором признак достоверен, сокращается с ростом шага и с ростом размерности системы. Для конкретной размерности нужно отладить программу на промежутке достоверности для эталонных систем (устойчивых, неустойчивых и асимптотически устойчивых), после чего на том же промежутке выполнять компьютерный анализ устойчивости исследуемой системы той же размерности. Для формирования признаков применимы разностные методы высших порядков, а также другие приближенные методы, включая кусочно-интерполяционные приближения решений с итерационным уточнением [15, 16].
Расширенные условия устойчивости для аддитивных преобразований разностных решений. В исходных условиях (12) эквивалентно соотношению
lim £ln(1 + Di-t h) _
У (t)-yt(t)=(Ук„ -y„) V/ e[t„,«>), Vкe 1,n . Поэтому лемма 2 дословно сохраняется при замене соотношения (13) на соотношение
lim S ln (1+ D^h) < cu, cn = const, Vte[t0,<»), V к el, n,
(31)
(32)
соотношения (14) - на соотношение
lim lim £ln(1+D(k)h) = -œ Vkein . (33)
t ^ œ i ^œ /
В условиях теоремы 1 при достаточно малом h под знаком логарифма исключены неустранимые особенности - согласно (18), (19) lim D(k) =ßt < L = const, lim D™ = ßei < L Vke 1n
, и эквивалентность (12) и (31) сохраняется в условиях этой теоремы. Соответственно, (12) - (14) эквивалентны соотношениям (31) - (33). Таким образом, условия устойчивости (32), (33) сохраняются
lim 2Д® h < cn, cu = const, Vt e[t0,»), V к e1, n
как в ограничении (6), так и при его нарушении на любом конечном множестве точек произвольного
ln(1+D\l]h)
отрезка [t0, t]. С учетом lim-—-=1 то же
Di - 1h
относится к следствию из [4]:
Следствие 9. Формулировка и утверждения леммы 2 сохраняются в условиях теоремы 1 при замене соотношения (13) на соотношение вида
i = o
(34)
£ = 0
и (14) - на соотношение
lim lim £D(k)eh = -œ Vke 1,П . (35)
t ^ œ i ^ œ «
элементами
являются
дискретные значения функции
Dk (t) =
f (t ,Y ) - f (t ,Y )
Ä (7 ) - Jk (7 )
(36)
dV —
—=U(t,V), V(t0)=O, dt
(39)
Ограничение (6) снимается также с интегральных аналогов [4] (34), (35) где lim ¿DC-<h - предел
i ^^¿=0
интегральной суммы на [t0, t]. Согласно (8), (11) ее
D® = fk(ti-£ ,~i-£} fk(ti-£ ,Yi-£}, £ = 0,1, ...,i, как
yki-£ - yki-£
В условиях теоремы 1 функция (36), за исключением устранимых особенностей, определена и непрерывна V7e[t0, t], Vk из (1). Путем замены
lim ZD^h на интеграл [4] получится
i ^ϣ=0
Теорема 3. Утверждения следствия 9 сохраняются в условиях теоремы 1 при замене (34) на соотношение
J Dk (t ) dt < cn, cu = const, Vte[t0, œ), V k e 1, n , (37)
t 0
и (35) - на соотношение
lim Jd (7) dit =-œ Vk ein . (38)
t ^ œ J 10
Здесь и ниже интегралы могут оказаться несобственными с устранимыми особенностями подынтегральной функции. В дальнейшем это дополнительно не оговаривается.
Случай точки покоя. Пусть система (1) преобразована к виду, при котором анализ устойчивости решения сводится [6] к анализу устойчивости точки покоя:
где U(t,V)=(ui(t,V),u2(t,V),...,щ(t,V)), V = (Vl(t),V2(t),...,Vn(t)), O = (0,0,...,0)T . В случае (39) вычитаемые в числителе и знаменателе (36) становятся нулями. Ниже для простоты возмущение нулевого решения (39) не будет отмечаться волной. Аналогом теоремы из [4] является
Теорема 4. Теорема 3 сохраняется для нулевого решения системы (39) при замене (37) на соотношение
JUk(t,V) dt < с(11), c(11) = const, Vte[t0,да), Vkel^n , 10 Vk (t)
и (38) - на соотношение
щ(t,V) , w, —
J kK' ' dt = -да Vkel,n . t0 Vk(t)
Здесь и ниже подразумевается, что в точках устранимых особенностей t , где у (~) = 0, значе-Uk (t,V)
доопределяются:
. Ограничение (6) снимается
ния
V, (')
V,(~) V,(')
с утверждений [4], непосредственно вытекающих из теорем 3, 4. Два из них, теоремы 5, 6, в скорректированном виде даны ниже.
Теорема 5. Пусть 0 < '0, и для системы (39) выполнены условия теоремы 1. Если 3 Д > 0, такое что для всех решений V('), 0<|| V('0) || < Д, верно не-
равенство
< f (t) Vt e[t0, œ), и при этом
J f (t) dt < Ci, cu = const, Vkel, n, то точка покоя
V
k
системы устойчива. В частности, утверждение верно в случае f (t) = р/', '<-1, в=const, р=const . Если же V Д >0 3 V(t), 0< 11 V(t0) || <4, такое что
для некоторого k е1, n выполнено
u
—> gk(t) Vte[t0,да) , и при этом Jg4(t)dt=да , то Vk 10 точка покоя системы (39) неустойчива. В частности, утверждение верно в случае, если g (t) = р tа, -1<а, 0 <р, а=const, р=const .
= -ty\ + vkcos2 (vv) е
sin31 v|vr21 |а/ zbv'r
' ¿ = 0 ¿ = 0
Теорема 6. Если в условиях устойчивости теоремы 5 3 4 > 0, 4 < 4, такое что для всех решений системы (39) У(?) , 0<|| V?)||<4, Vке1,П
выполнено — < ^ (?) Vt е[?0, да), и Г /к (?) dt=-да,
то точка покоя системы асимптотически устойчива. В частности, утверждение верно в случае /к (?) = -р?в, где -1 < Р ,0 < р, Р=const, р=const.
Пример 1. Условиям теоремы 6 удовлетворяет система
=-r (1 + sin6(v3 + v3)) v/', ke1,n,re1,n, r Ф k, (40)
dt
где 0 < ?0 (?0)=0 Vke1,п , 0<р, -1 <р< 0; 0 < р , 0 < ц; а, Ь - вещественные числа р
, Vze0,ц. В самом деле, для (40) выполнено
J fk (t) dt = -да .
Кроме
того,
1 sin31 vkvr21 Za,t' Y,b,v'r
7 vR , 1 2/ \ V I ¿=0 1=0
- kрt' +—cos (v^) e 1 v f (t) = -V + -2 e .
< fk (t),
где
t2
J fk (t) dt = ^р
— lim t '+1--—
' + 1/-да ' + 1
С 1 +
Далее,
Ч 1
- lim - +—
t ^да t t v 1 l0
-r (1 + sin6(v3 + v3)) t" < f (t), где f (t) = r/ , r Ф k, и Jf (t) dt =-да. Следовательно, точка покоя си-
10
стемы асимптотически устойчива.
На основе теоремы 2 и следствия 5 для n = 2, t0 =0.5, ' = -0.5, р=0.5, p = 3, q = 3,
a = 1, b =1, ¿ = 0,1,2, это же утверждение иллюстрирует программа Stability при соответственных изменениях:
function f1 (x,y 1 ,y2 :extended): extended;
begin f1:=-1/2*y 1/sqrt(x)+1/sqr(x)*y 1*sqr(cos(y1*y2))*
exp(sin(sqr(y 1)*y1* sqr(y2)*( 1 +x+sqr(x)+x* sqr(x))*(1 +y2+sqr(y2)+sqr(y2) *y2)) * sqr(sin(sqr(y 1)*y1 *sqr(y2)*(1 +x+sqr(x)+x*sqr(x))*(1 +y2+sqr(y2)+sqr(y2) *y2))))end; function f2(x,y 1 ,y2 :extended): extended;
begin f2:=- 2/sqrt(x)*y2*(1+ sqr(sqr(sin(sqr(y1)*y1+sqr(y2)*y2)))* sqr(sin(sqr(y 1)*y1+sqr(y2) *y2)))end;
||Z|| = J z2 + z2 : 5.33x10-1,., 6.17x10-3,., 9.66x10-10,..., 9.24x10-14,..., 4.92x10-16,..., 2.56x10-16.
Если в (40) множитель — заменить на -1,
то
, откуда следует асимптотическая устойчивость ?2 точки покоя. При сохранении последнего набора
полученная система снова будет удовлетворять параметров соответственно измененная программа условиям теоремы 6, можно положить /к (?) = -кр?в БйЫШу в этом случае выдаст признак
Ш: 3.99Х10-5,., 7.65Х10-12,., 8.26х10-21,..., 2.31х10-38,..., 1.91х10-54,..., 1.42х10-63,..., 1.45х10-66.
Напротив, если множитель -1 заменить на +1
, то точка покоя системы окажется неустойчивой. Правая часть к -го уравнения системы с этими изменениями примет вид ык (? ) = ук (?) ^ (?), где
.in3 I vly"2 \TPatt1 JEbv^ , ,
cos2 (vtv ) e 1 U=0 '=0 JJ>(1 - e4 - e) e-. Поскольку ' < 0, то, начиная с некоторого / > t0, вы-
полняется
-kрtв| < (e4 + e) е-1 We[/,да). С
точно-
"31 vly] \ £я,/ Zbvf l.f=0 ¿=0
стью
uk(t)
до
устранимых
v,
(t)
=gk(t) Vte[t0, да),
особенностей где
g (t) > (1 - 2e4 - 2e) e"1 V t e [/, да). Априори
e
gk (t) = -kрt^ + cos2 (vkvr) e \ v'=0 ¿=0 . Если предположить обратное, - что нулевое решение устойчиво, то Ve> 0 3Д>0: 0< V(t0<А влечет
(t) |< e Vke1,n, и
К(t)vr(t)|< e2 Vke1,n, Vre1,n, Vt e[^,да). С учетом cos2 (vtvr )=1-sin2 (vtvr), в предположении, что e достаточно мaло, cos2(vkvr)> 1 -e4 -e. Отсюда J^(/)dt =да. Присоединение к промежутку инте-
можно выбрать так, что e4 + e < ^, тогда
,. 1 t 1 t
gk(t)>1 e-1 Vte[t1, да), J gk (t) dt > - J e~ldt, и
2 t1 2 t1
и
грирования [?0, ? ], где ^ (?) ограничена, не изменит значения | g;. (?) dt = да. По теореме 5 нулевое
? 0
решение неустойчиво вопреки предположению. г II: 1.75Х103,..., 6.05Х1022,., 8.43х1040,..., 1.
Следовательно, предположение неверно, и точка покоя неустойчива.
Видоизмененная для данного случая программа Stability при сохранении предыдущих параметров непосредственно даст численный признак неустойчивости
x1050,..., 1.17x1066,., 4.96x1074,..., 8.42x1077.
Целесообразно заметить, что можно параллельно приближать конечное число возмущенных решений в малой окрестности невозмущенных начальных данных с вычислением норм из (20), (21). Одновременный вывод норм по ходу решения позволит с большей достоверностью выявить наличие или отсутствие границ, определяющих устойчивость, по нескольким признаковым последовательностям.
От известных аналитических методов [7, 8, 17] изложенные отличаются построением численных признаков устойчивости на основе отношения возмущения решения к возмущению начальных значений и их компьютерной реализацией в сравнительно общем случае. Отличия относятся к актуальным исследованиям [18 - 20], компьютерному анализу устойчивости [21, 22] и вычислению функций Ляпунова [23 - 25].
Заключение. Предложен компьютерно-ориентированный метод анализа устойчивости по Ляпунову на основе приближенного решения задачи Коши для системы ОДУ. Необходимым и достаточным условием устойчивости решения является равномерная ограниченность на полуоси отношения возмущения к вызвавшему его возмущению начальных значений, варьируемому в некоторой окрестности, при этом для асимптотической устойчивости необходимо и достаточно стремление отношения к нулю с ростом независимой переменной. Численное моделирование устойчивости реализуется дискретизацией этого отношения по ходу решения системы. Результаты моделирования показывает достоверность численных признаков устойчивости для систем ОДУ сравнительно общего вида. Аналогичное исследование устойчивости при возмущении параметров может указать средства применения метода для мониторинга автоматизированных систем управления в реальном времени.
Работа поддержана грантом РФФИ: проект 16-07-00100а «Компьютерные методы варьируемого кусочно-полиномиального решения дифференциальных уравнений и анализа устойчивости».
Литература
1. Ромм Я.Е. Параллельные итерационные схемы линейной алгебры с приложением к анализу устойчивости решений систем линейных дифференциальных уравнений // Кибернетика и системный анализ, 2004, № 4, с. 119-142.
2. Ромм Я.Е. Мультипликативные критерии устойчивости на основе разностных решений обыкновенных дифференциальных уравнений // Кибернетика и системный анализ, 2006, № 1, с. 127-142.
3. Ромм Я.Е. Моделирование устойчивости по Ляпунову на основе преобразований разностных схем решений обыкновенных дифференциальных уравнений // Известия РАН. Математическое моделирование, 2008, Т. 20, №12, с. 105-118.
4. Ромм Я.Е. Компьютерно-ориентированный анализ устойчивости на основе рекуррентных преобразований разностных решений обыкновенных дифференциальных уравнений // Кибернетика и системный анализ, 2015, Т. 51, № 3, с. 107-124.
5. Чезари Л. Асимптотическое поведение и устойчивость решений обыкновенных дифференциальных уравнений. - М.: Мир, 1964, 478 с.
6. Демидович Б.П. Лекции по математической теории устойчивости. - СПб.: Изд-во «Лань», 2008, 480 с.
7. Гантмахер Ф.Р. Теория матриц. - М.: Наука, 1988, 552 с.
8. Ромм Я.Е. Программируемые критерии устойчивости по Ляпунову. I / ДЕП в ВИНИТИ 21.06.2005, № 879- В2005, 24 с.
9. Ромм Я.Е. Программируемые критерии устойчивости по Ляпунову. II / ДЕП в ВИНИТИ 21.06.2005, № 880- В2005, 26 с.
10. Катрич С.А. Разработка и исследование программного моделирования устойчивости решений нелинейных дифференциальных уравнений на основе разностных методов. Автореферат диссертации на соискание ученой степени канд. техн. наук. Изд-во ТРТУ, Таганрог, 2006, 20 с.
11. Ромм Я.Е., Заярный В.В. Численный эксперимент по компьютерному анализу устойчивости линейных систем дифференциальных уравнений // Современные проблемы науки и образования, 2015, № 1; URL: http://www.science-education.ru/121-19424 (дата обращения: 29.05.2015).
12. Буланов С.Г. Разработка и исследование методов программного моделирования устойчивости систем линейных дифференциальных уравнений на основе матричных мультипликативных преобразований разностных схем. Автореферат диссертации на соискание ученой степени канд. техн. наук. Изд-во ТРТУ, Таганрог, 2006, 20 с.
13. Ромм Я.Е., Буланов С.Г. Численный эксперимент по компьютерному анализу устойчивости линеаризованных систем нелинейных дифференциальных уравнений / ДЕП в ВИНИТИ 14.07.2016, № 102, В-2016, 18 с.
14. Катрич С.А., Ромм Я.Е. Численный эксперимент по выбору шага и оценкам временной сложности для компьютерного анализа устойчивости решений обыкновенных дифференциальных уравнений на основе преобразований разностных схем / ДЕП в ВИНИТИ 14.07.2016, № 103, В-2016, 63 с.
15. Ромм Я.Е., Джанунц Г.А. Компьютерный метод варьируемой кусочно-полиномиальной аппроксимации функций и решений обыкновенных дифференциальных уравнений // Кибернетика и системный анализ, 2013, № 3, с. 95-112.
16. Ромм Я.Е., Джанунц Г.А. Кусочно-интерполяционное решение двухточечной задачи Коши для обыкновенных дифференциальных уравнений с итерационным уточнением // Фундаментальные исследования, 2016, № 6 (часть 2), с. 308-317.
17. Арнольд В.И. Математические методы классической механики. - М.: Наука, 1989, 472 с.
18. I.V. Astashova. Application of dinamical systems to the study of asymptotic properties of solutions to nonlinear higher-order differential equations // Journal of Mathematical Sciences, 2005, Vol. 126, N 5, p. 1361-1390.
19. Robert Morgan. Linearization and Stability Analysis of Nonlinear Problems // Rose-Hulman Undergraduate Mathematics Journal, 2015, Vol. 16, N 2, p. 68-91.
20. Okereke R. Lyapunov Stability Analysis of Certain Third Order Nonlinear Differential Equations. Applied Mathematics, 2016, 7, p. 1971-1977.
21. K. Wright. Adaptive methods for piecewise polynomial collocation for ordinary differential equations // BIT Numerical Mathematics, 2007, 47, p. 197212.
22. P. Giesl and S. Hafstein. Computation of Lya-punov functions for nonlinear discrete time systems by linear programming // J. Difference Equ. Appl, 2014, 20, p. 610-640.
23. T. Zhaolu and G. Chuanqing. A numerical algorithm for Lyapunov Equations // Applied Mathematics and Computation, 2008, 202, p. 44-53.
24. L. Xiao-Lin and J. Yao-Lin. Numerical algorithm for constructing Lyapunov functions of polynomial differential systems // J. Appl. Math. Comput, 2009, 29, p. 247-262.
25. Peter A. Giesl , Sigurdur F. Hafstein. Revised CPA method to compute Lyapunov functions for nonlinear systems // Journal of Mathematical Analysis and Applications, February 2014, Vol. 410, Issue 1, p. 292306.
АНАЛИЗ ТЕРМИЧЕСКИХ ХАРАКТЕРИСТИК ВЛАЖНЫХ МАТЕРИАЛОВ КАК ОБЪЕКТОВ СУШКИ
Сажин В.Б.
доктор технических наук, профессор, академик, директор Фонда Российский инвестиционно-инновационный Фонд «Научная Перспектива»
Сажин Б. С.
доктор технических наук, профессор, академик, советник Фонда Российский инвестиционно-инновационный Фонд «Научная Перспектива»
ANALYSIS OF THE THERMAL PROPERTIES OF THE MATERIALS TO BE DRIED
Sazhin V.
Doctor of Technical Sciences, Professor, Academician, Director Russian "Scientific Perspective" Investment and Innovation Fund, Moscow
Sazhin B.
Doctor of Technical Sciences, Professor, Academician, Adviser Russian "Scientific Perspective" Investment and Innovation Fund, Moscow
АННОТАЦИЯ
Приведены результаты исследования взаимного влияния тепловых характеристик (теплопроводности, температуропроводности и теплоёмкости), а также их зависимость от внешних параметров (температуры, влажности), с учётом структуры материала, характеристик остова и видов связи влаги с остовом для материалов органической и неорганической природы (строительные материалы, почвогрунты и др.). Предложены наилучшие подходы к анализу теплофизических характеристик, рекомендованы модели, расчётные зависимости и установлены границы их рационального применения.
ABSTRACT
The results of research of mutual influence of thermal characteristics (thermal conductivity, thermal diffusiv-ity and specific heat), as well as their dependence on external parameters (temperature, humidity), taking into account the structure of the material, the characteristics of the core and the types of moisture due to the backbone for the materials of organic and inorganic nature (building materials , soils, etc.). Offer the best approaches to the analysis of thermal characteristics are recommended model, the calculated dependence and set the boundaries of their rational use.
Ключевые слова: сушка, материалы, свойства, структура.
Keywords: drying, structure, materials, properties.
