book1989

на границе единичного круга комплексной плоскости, причем на границе единичного круга нет кратных корней. Справедлива

Т е о р е м а 1. Условие	корней необходимо и достаточно для
устойчивости уравнения (8)	по начальным данным.

Д о к а з а т е л ь с т в о . Докажем сначала необходимость. Пусть уравнение (8) имеет корень q, для которого | ^ | > 1 . Задавая в качестве начальных данных функции Vj—q1, j = 0, 1, .... m—1, по­ лучим решение vn = qn, n^.m , неограниченно возрастающее при п—>-оо. Для такого решения невозможна оценка вида (18) с кон­

стантой Ми не зависящей от п. Следовательно, условие |(?й|^ 1 , k= 1, 2, .... m, необходимо для устойчивости.

Пусть уравнение (9) имеет корень q кратности т> 1, для кото­ рого \q\ = 1. Тогда разностное уравнение (8) имеет решение nr~lqn, растущее при п—>-оо как пг~'у и, следовательно, в этом случае оцен­ ка (18) также невозможна.

Прежде чем переходить к доказательству достаточности усло­ вий теоремы 1, необходимо провести некоторые вспомогательные построения.

Запишем (8) в виде эквивалентной системы уравнений

— V n-m + 1> ■ • • > Vn-1 — Vn- i ,

V n

II 1

atn Oo

C2 a 1

a m -1	a i
Vn—m-ri	. . --------- -v,
a0

и представим эту систему в векторной форме
l/n=SVn_1,		п=пг,	т + 1 , . . . ,	(19)
где Vn {vn—m+i, V n —	• • • » V n )	,
	0	1	0	0
	0	0	1	0
S =	0	0	0	1
	am	ani~1	am-2 _____
_	a0	a0	a0

Начальный вектор Vm~ —(^o> V it . . ■, Dm-,y задан.

Нетрудно проверить, что множество собственных чисел матри­ цы S совпадает с множеством корней характеристического уравне­ ния (9).

Л е м м а 1. Если выполнено условие корней,	то существует
норма ||'- Н, вектора такая, что для подчиненной нормы матрицы S
справедливо неравенство ||S ||,^1 .
Д о к а з а т е л ь с т в о . Заметим, что данная	лемма	уточняет
лемму 1 из § 3 гл. 2, поэтому их доказательства	похожи	(см. [2,
с. 338]). С помощью преобразования подобия
S==QSQ~l		(21)
		241

п р и в едем S к м о д и ф и ц и р о в а н н о й ж о р д а н о в о й ф о р м е

			~Si			О
	5 =
			О			J
-где S„ —либо число, либо жорданова клетка
				_		О
	s A— 1
			О
•g„ —собственное число		матрицы		S,		е > 0 —любое число. Оценим
норму			т
						max (| su | + | h u 1 1).
	IS ||с = шах		2 1 S‘7 I =
	i<:i<+i .				u£i<m
		l=i
Здесь	совпадает	с одним		из корней характеристического
уравнения, а внедиагональнын элемент
	_[ 0 в случае			простого корня sl£,
	s£,i+i = |			кратного корня
	| е в случае

Если s(i+1 = 0 для некоторого i, то по условию леммы |s,i|^l. Если же sii+1 = e, то si£ —кратный корень, и согласно условию кор­

ней выполняется строгое неравенство |s„|<E Но тогда при доста­ точно малом е имеем

| Sи | + | sii,-i | <С 1 ■

Таким образом, выбирая е достаточно малым, получим |[S|lc^ -s^l. Введем норму вектора

ll</ll.= IIQyllic,	(22)
где Q определено согласно (21). Тогда получим	||S||,= ||S||C^ 1.
■Лемма 1 доказана.

Завершим теперь доказательство теоремы 1. Покажем, что ус­ ловие корней достаточно для устойчивости уравнения (8) по на­ чальным данным.

Из уравнения (19), учитывая лемму 1, получим неравенство

I * п- 1 II * -

и, следовательно,
IIKJI-OIEm-il!., п=т, т+ 1,...	(23)

По определению (22) нормы || -1|. имеем

l!V|j. = ||Q E ||c ^ llQ I |c |!E!ic.

242

С другой стороны, для любой невырожденной матрицы Q спра­ ведливо тождество

V=Q~lQV,

из которого следует оценка
Таким образом, если норма ||-||. определена равенством		(22),
то выполняются оценки
(IIQ-'y-'l^clKIIFII.^IIQIicilFllc.		(24)
Из (23) и (24) получаем
ilVJc^ATIIV^IIc,		(25)
где Afi= llQ-'llcIlQIlc. Заметим, что константа Aft не зависит		от п.
Из (25) следует неравенство
| Уп | ^ 44j шах	|и/|,	(26)
означающее устойчивость уравнения	(8) по начальным данным.
Теорема 1 полностью доказана.

4. Оценка решения неоднородного уравнения. Рассмотрим за­

дачу Коши для неоднородного разностного уравнения
ЯаУп+ GiUti-i + • • ■+	(27)

где п = т, т+ 1........ величины у0, уи . . . , ym- t заданы, правая часть gh, &= 0, 1,... — заданная функция. Если ааФ0, то для каждой пра­ вой части решение задачи Коши существует и единственно. Оно мо­

жет быть найдено по рекуррентной формуле
У п ~ — а0 Уп-т	Яо Уп—т+1	а1	I	Т£п-т
		я0	Уп~1 “1	До	у
	п = т, т+ 1,				(28)

исходя из заданных начальных условий у0, уи . . ., y m - i и известной правой части qk.

В предыдущем пункте получена оценка решения однородного уравнения через начальные данные, означающая устойчивость по начальным данным. Получим теперь аналогичную оценку решения неоднородного уравнения (27) через начальные данные и правую часть.

Основным результатом здесь будет доказательство того, что

если однородное уравнение (8) устойчиво по начальным данным,

то для неоднородного уравнения (27)	справедлива оценка
	п—т
\Уп\< -иа шах | г/, ] +	Мг V т | gk |,	(29)
	k—Q

где Mi и Мг не зависят от п.

Выполнение оценки (29) означает по определению устойчивость уравнения (27) по правой части. Таким образом, будет показано,

243

что из устойчивости по начальным данным следует устойчивость по правой части.

Представим уравнение (27) в векторном виде
Е«= 5У„_1+ т:(7п-1,	п = т, т+ 1, . .. ,	(30)
ГД6 Уп~ (i/n—m+i, Ул—m+lj *• *1 Уп) j
Gn~i — (О, 0, ..	о, 8п~т ) ,	(31)
	йо /

матрица S определена согласно (20).

Пусть уравнение (8) устойчиво по начальным данным. Тогда согласно теореме 1 выполнено условие корней. При этом условии в соответствии с леммой 1 для некоторой нормы матрицы 5 спра­ ведливо неравенство ||5 ||,^ 1 . Таким образом, из (30) получаем не­ равенства

I I I I •+ т||G^li,, k = m, т+ 1, ...

Суммируя эти неравенства по k от т до п, получим

\Yn II.< 1 ^ - 4 + 2 *11 б*II..		(32)
k=rn-i
Используем далее неравенство (24), устанавливающее эквива­
лентность норм М |. и IIНе. Тогда из (32) получим
(«<Г ЦсГ1I У„Цс ^ IIQlie ( IJV * ||с +	2 * h)1 •	(33)

\k—m-i J

Учитывая специальный вид (31) вектора Gh, имеем
	GkНс - I &khl—П
	I	«0 I
и, следовательно,		n—m
/г—l	t|G.lc= -f		T|S.|.
2		S
k—m-1	0	k=Q
Отсюда и из (33) получаем требуемое неравенство
1Уа|с М Мх max I у/1 +			n—m
			М, 2 * Iff*|,	(34)

где Ali=[|Q_ti!cllQIIc, М2= /И^о1. Напомним, что Q — матрица, пре­ образующая S согласно (21) к модифицированной жордановой

форм е.

5.Оценки погрешности разностного метода. Вернемся к урав­ нению для погрешности (3), полученному в п. 1. Нам нужно оце­

нить погрешность метода zn через погрешность аппроксимации

/г = 0, 1, ... , п—т. Если бы функция (рд равнялась нулю, k = 0, 1, ...

244

п—m, то достаточно было бы воспользоваться оценкой (34). Однако наличие функции фА, зависящей от решения z и характери­ зующей нелинейность задачи, усложняет получение требуемых оце­ нок. Следуя [2, с. 484], докажем, что справедлива

Те о р е м а 2. Пусть все корни характеристического уравнения

(9)лежат внутри или на границе единичного круга, причем на гра­

нице нет кратных корней. Пусть				| fu(t, и) |			для /е(О , Т]. Тогда
при			=						(35)
					2 I Ь0 |	L
для решения уравнения (3)			справедлива оценка
| гпК	max	|z,-| +	М2	V х | фй|,			п = m, m +		1, ... , (36)
где		M ^W Q U Q - Ъ
=	Cl =	MlVL, 0 =		^ - 2		(1‘	\bk\ + M		\b0
				a:о k=i
Д о к а з а т е л ь с т в о . Согласно формуле конечных приращений
Лагранжа имеем
f (^n—hi Уп—h)			f {.in—ft) Ип—ь)				C—	й,
где
ln—h /и {tn—к) Un—ft)i Un—h //n—ft4~ 02-т,—hi 0-C;0-C~ 1,								&= 0, 1,2,
Представим функцию ф„-m, определенную согласно									(5), в виде
	фп—т — Ь0102п			^ b fjn ^ ft^ n — k'					(37)
					k=i
Подставляя (37) в уравнение (3), получим
		т
(а0—Ьа10%)г„= 2			(—ak + xbkln.k)				+ тф„_т .		(38)
		k — \
Из условия (35) следует, что
		к - л м ; ^ >0’							(39)
поэтому уравнение (38) можно разрешить относительно zn:
„	v\	ak'5C		k	,	| ^			(40)
=	2 ------------— ------Zn- k +					T -	QQ— T&o^O
	A=l	a0	tb0l0
Добавляя к правой части уравнения (40)							и вычитая из нее выра­
жение			т
				4	7
			2	а0	Zfl-ky
									245

перепишем (40) в виде

2п — — 2		а		1		-л	(41)
		— Zfl~k “Ь т 2			Vnk2n-k + т — ———
где			.	а ,		ао - tv»
				aob^n-kr~ аМ>
		Vak					(42)
				а„ (ао— т&0/0)
Представим уравнение (41) в векторной форме. Для этого вве-
дем векторы	Zn		(^п-т+1, 2n—m-i-2t • • • , ^n)			,
	Yn- ! =		fo,	.... 0,	%~m
			4		flo — TOQ'O
и матрицы					1	0	0 -1
~o	.. . 0 "			г 0
0	.. . 0	,	s =	0	0	0	1
				am	V i	am-2	____£r_
	' *
У run				L an	«0	aQ	a0 J

Тогда получим, что уравнение (41) эквивалентно векторному урав­ нению

Z n = SZ,_i+ TVn_ iZ n_1+ T Y „ _ t.

(43)

Для оценки решения уравнения (43) используем лемму 1 из п. 3. Согласно этой лемме, ||S||,=^1 в некоторой норме || -1|.. Поэтому из

(43) получим
l!ZJ,<||Z„_1||,+ T||Vn_1!!.l|Zn_1H.+ Tll1F„_1||..	(44)

Покажем, что ||У„_,||, ограничена числом, не зависящим от п и

т. Согласно (39), (42) имеем
I ао11 bk1+	| ak1| b0			1 , 2 , . . . ,	m,
\vnk | < 2 -
поэтому		m
I Vn-1 lie =			1v«b I ^	vL>
		2
где		k=\
о m	(K l	I	1+ M	|ь0j).	(45)
V= — 2

Q0 fe=l

Для любого вектора x имеем

||V _1*||. = ||QVn_1x||c= || (QVn-.Q-1) (Q*) ||c<

<IIQ Vn_1Q-,||cllQJc||c< ^ .IIV n-1||c|i*||.t

246

где Mi = ||Q||c||Q_1|lc. Следовательно, lll'»-ll|.^All||l'n_l||c<A l10L.

Подставляя эту оценку в (44), приходим к неравенству

H 2 J .< (l + TC1)||Zt_1||. + T||4f4- i||.l

(46)

где &= m, т+ 1, . . . , c ^ M ^ L . Из (46)

получим

j|2n|l.^ (l +T Cir

m||Zm.1I. + т 2

k —m

<

ес‘‘™

I!

I. +

т

2

II

II.,

(47)

k—m

где tn- m=x{n—m ) ^ T .

Далее, учитывая (24), (39), имеем

II Znl >

(«ОТ1Нс)"11Z„ ||c >

(1QT1 lie)"11zn I,

IZm_!1. ^

|i Q||c I Zm_!||c=

IQ lie

max

| г,-1,

1*A-i II.^

I! Q1U

2 II

Q|U

I ^

I•

IIQlc 1 Ъ-г ||c = -— — 7 7 7 1tu n I

I

a0

I

I aa — Jb0l0 \

Подставляя

эти неравенства

в

(47)

и обозначая A^HIQIIcX

2

XlIQ-'lic, Мг=М1-— - , получим оценку

I «о I

п—т

lanK /H ie^7

шах

| z, | +

М2 V т | ф* |,

совпадающую с (36). Теорема 2 доказана.

условие

корней,

С л е д с т в и е .

Если

О ^ и т ^ Г ,

выполнено

| Zj|—>-0 при т-»-0, /= 0, 1,

. . . , т—1, и разностное

уравнение

(2)

аппроксимирует исходное уравнение (1), то решение разностной за­ дачи (2) сходится при т-vO к решению исходной задачи (1).

Доказательство следует немедленно из оценки (36), если учесть,

что

п—т

УI т | ф * | < 7 шах | ф*|-

иo^h^n-m

§5. Численное интегрирование жестких систем обыкновенных дифференциальных уравнений

1.Условно устойчивые и абсолютно устойчивые разностные ме­ тоды. Для задачи Коши

= f(t,u), / > 0, и (0) = и0

(1)

at

247

будем рассматривать разностные методы вида

т	ak	т	*/«-*). п = т < т + 1. ■• •	(2)
2	•— Уп-к =	2 bkf
*=О		А=О

В § 4 показано, что устойчивость и сходимость метода определя­ ются расположением корней характеристического уравнения

т

2 akq'n-* = 0 .

(3)

k=0

А именно, требуется, чтобы все корни удовлетворяли условию | q \ ^ =^1, причем корни q, для которых |^| = 1, не должны быть крат­ ными.

Эти условия устойчивости являются очень общими и не могут учесть многие характерные свойства решений исходной дифферен­ циальной задачи (1) и аппроксимирующего ее разностного метода

(2). Они означают лишь, что все решения однородного разностно­ го уравнения, соответствующего (2), остаются ограниченными при п—»-оо.

В частности, при таком подходе коэффициенты bh,k = 0, 1, . .. , m, входящие в правую часть уравнения (2), никак не влияют на устой­ чивость.

Предположим, однако, что заранее известна та или иная харак­ терная особенность в поведении решения исходной дифференциаль­ ной задачи. Тогда естественно требовать, чтобы эта особенность со­ хранялась и у решения разностного уравнения. Такое требование приведет к сужению класса допустимых разностных методов. В на­ стоящем параграфе будут рассмотрены методы, предназначенные для расчета асимптотически устойчивых решений уравнения (1).

Рассмотрим сначала характерный пример. Уравнение

~ — Xu, t > 0, и (0) = н0,

(4)

at

где Я<0, имеет решение

u(t) =иаеи,

монотонно убывающее при t-*-oo. При любых т> 0 для решения это­ го уравнения справедливо неравенство

\u ( t+ x ) \^ \u ( t) \,

(5)

означающее устойчивость решения u(t).

Естественно требовать, чтобы и для решения разностной задачи,

аппроксимирующей	(4), выполнялось бы неравенство, аналогичное
(5). Рассмотрим с этой точки зрения метод Эйлера
^	— = Ц п ,	п = 0 , 1 , . . .	(6)

т

Из уравнения (6) получаем

yn+i=qyn, q= 1+тЯ.

248

Оценка вида (5), т. е. неравенство
I Уп+i | | уп| » м = 1, 2, . ..	(7 )
для метода (Ь) будет выполнено тогда и только тогда, когда	^

^ 1 . В случае Я<0 это условие эквивалентно следующему ограни­ чению на шаг т:

0 (8)

Таким образом, разностный метод (6) устойчив в смысле выпол­ нения оценки (7), если шаг т удовлетворяет неравенству (8).

Разностный метод (2) называется абсолютно устойчивым, если он устойчив при любых т>0, и условно устойчивым, если он устой­ чив при некоторых ограничениях на шаг т. Мы видели, что метод Эйлера (6) условно устойчив при условии (8). Примером абсолютно устойчивого метода для уравнения (4) с Я<0 является неявный ме­ тод Эйлера

— fi-Упы,

1

для которого \ q\ = | (1—тХ)~‘| <1 при любых т>0.

Приведенные здесь простые примеры характерны тем, что и для более общих асимптотически устойчивых систем дифференциаль­ ных уравнений явные разностные методы являются условно устой­ чивыми, а среди неявных методов существуют абсолютно устойчи­ вые методы.

Условная устойчивость является недостатком явного метода, так как вынуждает брать слишком мелкий шаг т. Например, если %= = —200, то условие (8) выполнено при т^0,01, и для того чтобы вычислить решение «(/) при / = 1, надо сделать сто шагов по методу Эйлера. Неявный метод лишен этого недостатка, однако его приме­ нение к задаче (1) приводит к необходимости решения на каждом шаге системы алгебраических уравнений, в общем случае нелиней­ ной.

2. Понятие жесткой системы дифференциальных уравнений. Многие из рассмотренных в § 1—4 численных методов интегрирова­ ния обыкновенных дифференциальных уравнений переносятся без изменений на системы дифференциальных уравнений. Однако в случае численного решения систем уравнений могут появиться до­ полнительные трудности, связанные с разномасштабностью процес­ сов, описываемых данной системой.

Поясним характер возникающих трудностей на примере систе­ мы, состоящей из двух независимых уравнений

- - L	+ aiU i = о,	t>- 0 ,	(9)
dt	' х 1
где at и аг— положительные постоянные.
Система (9) имеет решение
	их (/)= «!(0)е~а' \	«2(/)= м2 (0)e~a’J ,
			249

монотонно убывающее с ростом t.			Предположим, что а2 гораздо
больше, чем аи Тогда компонента			u 2 ( t )	затухает гораздо быстрее,
чем щ(0,	и, начиная	с некоторого	t ,	поведение решения	u ( t ) =
= { u i ( t ) ,	u 2 ( t ) } почти	полностью определяется компонентой			U i ( t ) .

Однако оказывается, что при решении системы (9) разностным ме­ тодом шаг интегрирования т определяется, как правило, компонен­ той u z { t ) , не существенной с точки зрения поведения решения си­ стемы. Например, метод Эйлера

+ djU* = О, + аги\ = О, ( 10)

где и" = u,(tn) , г = 1, 2, будет устойчив, если шаг т удовлетворяет одновременно двум неравенствам Tfli^2, %a2^L2. Поскольку а{^> ^>а,>0, условие устойчивости приводит к ограничению т==:2/й2.

Приведенный пример может показаться искусственным, так как ясно, что каждое из уравнений системы (9) следует решать незави­ симо одно от другого со своим шагом интегрирования т , /= 1, 2, т1^ 2 /а 1, т2=£^2/аг. Однако аналогичные трудности возникают и при решении любой системы обыкновенных дифференциальных урав­ нений

— = Аи,	(И)
dt

если матрица этой системы имеет большой разброс собственных чи­ сел.

Предположим, например, что матрицу А системы (11) можно привести пре­ образованием подобия Q - ‘AQ к диагональному виду. Тогда замена u= Qv пре­ образует систему (И ) в систему независимых уравнений

dv

Q-UQv,

dt

матрица которой имеет те же собственные числа, что и матрица А.

Сформулируем теперь определение жесткой системы уравнений. Рассмотрим сначала систему (11) с постоянной, т. е. не зависящей от t матрицей А. Система дифференциальных уравнений (11) с по­ стоянной матрицей А(т Хт ) называется жесткой, если

1)	ReXbCO, /г=1,2, ...,п г

(т. е. система асимптотически устойчива по Ляпунову), 2) отношение

max	| Re Xk \
s __	_______
1m i n	[ R e X 4 |

велико.

Число s называется числом жесткости системы (11). Второе тре­ бование не указывает границу для s, начиная с которой система ста­ новится жесткой.

250