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

Конкретный метод считается устойчивым, если кумулятивный эффект всех ошибок округления, возникших в процессе применения данного алгоритма, является пренебрежимо малым. Говоря более определенно, рассмотрим ошибки

?? = Г;?-Т?, (4.15)

вносимые в узловых точках (/, я), где /= 2, 3, 1 и

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

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

нии данных о *Т1и У1 и yl+u так что

ТГ = S (Т i) + (1 - 2s) (Т?) + S (Ti). (4.16)

Вспоминая, что точное решение алгебраических уравнений Г/ удовлетворяет уравнению (3.41), и пользуясь формулой (4.15),

8 К- Флетчер, т. 1

С учетом (4.16) получим однородное алгебраическое уравнение Г = sll-i + (1 - 2s) i; + slhi. (4.17>

Если предположить, что заданы начальные и граничные значения, то все начальные ошибки / = 2, 3, 1, а также граничные ошибки и g, /г = О, 1, 2, ..., в соответствии с уравнением (4.17) будут равны нулю. Если при расчете значения Г? в каком-то внутреннем узле не будут внесены те или иные ошибки округления, то и результирующие ошибки решения останутся равными нулю.

Двумя наиболее общеупотребительными методами анализа устойчивости являются матричный метод и метод Неймана В основе обоих методов лежит возможность предсказания, будет ли иметь место рост ошибки как разности между истинным решением численного алгоритма и фактически вычисляемым решением, т. е. разности, учитывающей погрешности округления.

Альтернативный вариант анализа устойчивости связан с предположением о том, что начальные условия представляются с помощью рядов Фурье. Каждая гармоника, или мода, ряда Фурье будет возрастать или затухать в зависимости от вида дискретизированного уравнения, которое, как правило, сама порождает конкретную форму выражения для показателя роста (или затухания) соответствующей моды. Если данная мода может неограниченно возрастать, то дискретизированное уравнение имеет неустойчивое решение. Эта интерпретация понятия устойчивости [Richtmyer, Morton, 1967] непосредственна используется в методе анализа устойчивости, предложенном Нейманом (см. п. 4.3.4 и 4.3.5). Неограниченный рост некоторой конкретной моды остается возможным, даже если дискре-тизированные уравнения решаются точно, т. е. если ошибки округления отсутствуют. Если же ошибки округления будут внесены, то сама неустойчивая природа дискретизированных уравнений будет вызывать неприемлемый рост ошибок. Следовательно, методика анализа устойчивости дискретизированных уравнений будет оставаться одной и той же независима от проявления присущей им устойчивости или неустойчивости.

4.3.1. Матричный метод: схема ВВЦП

Данный метод состоит в том, что система уравнений, определяющих распространение ошибки представляется в матричной форме, после чего мы анализируем собственные значения соответствующей матрицы. Для начала этот метод будет про-

иллюстрирован на примере схемы ВВЦП, реализуемой согласно (3.41).

Подставляя в (4.17) значения / = 2, 3, /-1 и замечая, что на границах ошибки равны нулю, т. е. что при всех п -будет = = 0, получим

+ = (1-25)з + 5?з ,

ir = s? , + (l-2s)i; + s+b iy+/=s?-2 + (l-2s)?-,..

(4.18)

В матричной форме это можно записать как Г+ = АГ, п = 0, 1,

(4.19)

где А - квадратная матрица порядка / - 2, а % - вектор длины / - 2; то и другое определяются выражениями

(l-2s) S

(1 - 2s) s s . s

s (l-2s)

Можно показать, что при увеличении п ошибки остаются ограниченными, если все собственные значения А. tn матрицы А различны и имеют абсолютные значения, меньшие или равные единице; иначе говоря, если

I К 1 для любых X. (4.20)

Далее формула (9.48) позволяет найти собственные значения трехдиагональной матрицы А в виде

=l-4ssin2(2(7), m=l, 2, /-2. (4.21)

Условие устойчивости (4.20) позволяет использовать только такие значения 5, которые удовлетворяют двойному неравенству

(4.22)