Введение в вычислительную физику - Федоренко Р.П.
ISBN 5-7417-0002-0
Скачать (прямая ссылка):
Устойчивость разностных схем. Рассмотрим схемы следующей формы:
*-~*n- = F(Xa), п = §, \, ..., N — \, N=T/x.
Как мы увидим в дальнейшем, такая форма охватывает важные классы схем Рунге—Кутты. Функция F(x), конечно, связана с правой частью уравнения /(х, t) (эта связь будет конкретизирована ниже). Все дальнейшее не претерпит никаких изменений, если вместо F(x) будут использоваться функции F(x, t, х), но ради простоты записи мы выбросим несущественные аргументы. Во все выкладки их при желании можно вписать механически, ничего не меняя.
Проверка устойчивости связана со сравнением решений двух «почти совпадающих» систем уравнений:
= F(Xn) -He;, X0=X*,
ПУп) +С Уо=У'о¦
72
ОСНОВЫ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ
[Ч. I
Относительно E1tl и є'' предположим, что они ограничены общей константой: HeI1II < є, ||е"И $ ?. Еще раз подчеркнем, что на самом деле ниже речь пойдет не о двух системах, а о семействах систем с параметром т. Чем меньше т, тем больше уравнений в системах, а нас будут интересовать оценки, равномерные по х. Ради простоты мы не пишем Xjl, F(xxn, г) и т.д., но «скрытый» аргумент х следует всегда иметь в виду.
Теорема. Пусть функция F(x) удовлетворяет условию Липшица с постоянной С: Ц/Х*) “-^(-У)И ^ С||х — у||. (Заметим, что это равномерная по г оценка, С от г не должна зависеть.) Пусть шаг х мал: Cr -«1. Тогда система разностных уравнений устойчива и имеет место оценка
Wxn-уп\\<ест \\х'а-Уі\\ + 2гест/с, Vn^T/x. (I)
(Таким образом мы докажем устойчивость разностных уравнений по начальным данным и правым частям.)
Доказательство. Перепишем уравнения в виде
*„+1 = *„ + xF(xJ + Т<’
(2)
Уп+i = Уп + x^yn) +
Вычитая второе уравнение из первого, получаем
IK-H - У„+іІІ < IK - yjl + *l№„) - ^(yjll + 2хє. Используя условие Липшица, преобразуем оценку
ІК+1 - yn+iH < (1 + Ct)IK - yjl + 2т?-
Применим эту основную оценку последовательно:
11*1 - УіІІ sS О + Сх)||х0 - у0|| + 2хє,
11*2 - JbJI < (1 + Ст)11*1 - УіІІ + 2хе <
« (1 + Сх)2||х0 - у0|| + 2хе[1 + (1 + Cl)],
11*з - УзІІ * О + Сг)11*2 - Уг\\ + 2х? ^
< (1 + Сх)3||х0 - yj + 2хе[1 + (I + Cx) + (I + Cx)2]. Легко угадывается и доказывается по индукции общая формула: ll*e-yJ«(l'+Cx)“||x0-y0|| +
+ 2хе[1 + (I + Cx) + ... + (I + Cx)"-1],
§ 7] ИССЛЕДОВАНИЕ СХОДИМОСТИ МЕТОДОВ рунге-кутты 73
или (после суммирования прогрессии)
уп\\ « (I + CxY + Ixz (1+Сх)” ,1 «
I +Cx — 1
2 г
(1 + Сх)п
I^O — З'о!! + Г
(3)
Заметим, что п H Т/х; следовательно, (I + Cx)'1 ^ (1 +• Сх)тн « ест при Сх« 1. Используя эту оценку в (3), получаем (1). Кстати, условие Cx-«1 не следует считать очень жестким, так как, напримрр, 1.5 « еол, а 1.1 « е0 095.
Теорема об устойчивости доказана. Применим ее к некоторым широко используемым на практике схемам.
Метод Эйлера с пересчетом. Перевод от известного X11 к новому хп+j делается в два этапа:
а) находится значение хп+ш = хп + 0.5т /(xn, tn);
б) вычисляется xn+1 = xn + X f(xn+m, tn + 0.5х).
Эту схему можно записать в общей форме:
*„+! = ** + т>-
ще F(x, t, х) = /(х + 0.5х f(x), t + 0.5х). Легко проверить, что если / удовлетворяет условию Липшица (по х) с константой С, то и F удовлетворяет этому условию с несущественно большей константой C(1 + Cx). Таким образом, схема метода Эйлера с пересчетом устойчива.
Проверим порядок аппроксимации, т.е. оценим выражение
T (*» + !--W,. *)•
Используем ряд Тейлора:
+ х) = агп + хІТа + ^Мга + O(Xi).
В силу уравнения 3? а = f(3?n, tn) имеем
= /Ж + A = /,(«'п. *„) /(*„. О + Л(^„> о-
Итак,
т (*n+i - *п) = / +1 Uj+/,I + 0(t2).
С другой стороны,
П*« ta, X) = / \хп +1 /(*„, tn), tn + Ij =
= A*V O n*»> 0 + X2ft(^n> O + O(^).
74 ОСНОВЫ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ [Ч. I
Объединяя оба результата, получаем
'„+I = T (*»+i - *п) - *», X) = 0(х2).
Таким образом, схема Эйлера с пересчетом имеет второй порядок аппроксимации и в силу теоремы § 6 — второй порядок точности (в предположении, что /(х, t) имеет две ограниченных производных и, следовательно, Sf (і) — три).
Методы Рунге-Кутты. Метод Эйлера с пересчетом является простейшим вариантом одной из наиболее распространенных в современной вычислительной практике схем численного интегрирования обыкновенных дифференциальных уравнений, объединяющей семейство методов с общим названием «методы Рунге—Кутты». Основу этих методов составляет ряд Тейлора. Связь 3?п+1 =
= af(tn + x) с Tn имеет форму
= + ^« + T *« + ••• + H 3fP + °Vk+i)- (4)
Приближенное решение находится из того же выражения, но с выброшенным остаточным членом:
*„+1 = +T+ ¦•¦ + Ti
Чтобы иметь вычислительную схему, нужно дать выражение для производных. В принципе здесь нет серьезных проблем:
~ f(Xn’ O’ = fх(хп> tJ f(xn’ *п) + ft(Xn' tn)
и т.д. Однако аналитические выражения начинают катастрофически усложняться в результате дифференцирования, и этот путь оказывается крайне неудобным. В методах Рунге—Кутты строится способ вычисления отрезка ряда Тейлора, требующий лишь вычисления /(х, t) в разных точках. Вот одна из распространенных схем.
