Устойчивость

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

Рассмотрим пример неустойчивой разностной схемы для задачи Коши дифференциального уравнения u¢ = a u. Выберем следующее однопараметрическое семейство разностных схем:

. (53)

Исследуем рост ошибки dyn начальных данных уравнения (53). Поскольку уравнение (53) линейно, постольку ошибка dyn удовлетворяет тому же уравнению (53). Изучим специальный вид ошибки dyn = ln. Подставим это представление в (53), тогда

. (54)

Решение квадратного уравнения (54) при h ® 0 дает следующие оценки корней

. (55)

Из оценок корней в (55) следует, что при s < ½ второй корень |l2| > 1, т.е. за один шаг ошибка возрастает в несколько раз. Проверим это.

На листинге_№5 приведен код программы, иллюстрирующей расчет по неустойчивой при s = 0,25 схеме (53) и по устойчивой схеме при s = 0,75. В начальных данных выбирались малые возмущения. Далее проводились серии расчетов с уменьшающимся значением шага сетки h. На рис.11 приведены итоговые графики зависимости значения возмущения в начальных данных на правом конце отрезка интегрирования в зависимости от шага сетки. Отчетливо видно сколь разительно отличаются друг от друга расчеты по неустойчивой и устойчивой схемам. Используя данную программу можно убедиться в пороговом значении параметра s = 0,5: при s < 0,5 схема неустойчива, при s ³ 0,5 — устойчива.

 

Листинг_№5

%Программа расчета по неустойчивой схеме при

%sigma=0,25 и по устойчивой схеме при sigma=0,75

%очищаем рабочее пространство

clear all

%определяем константу уравнения u'=alpha*u

alpha=1;

%определяем значения sigma=0,25; 0,75

sigm=0.25:0.5:0.75;

for s=1:length(sigm)

sigma=sigm(s);

%определяем начальное значение шага сетки

h=0.5;

for i=1:8

h=h/2;

x=0:h:1; N=length(x);

%определяем возмущения начальных данных

dy(1)=1e-6; dy(2)=1e-6;

%осуществляем расчет возмущения начальных

%данных на правом конце отрезка интегрирования

for n=2:(N-1)

dy(n+1)=(2+(alpha*h-1)/sigma)*dy(n)+...

(1/sigma-1)*dy(n-1);

end

%запоминаем возмущение на правом конце и

%шаг сетки

deltay(i)=dy(N);

step(i)=h;

end

%рисуем график зависимости возмущения на

%правой границе от шага сетки

plot(step,deltay);

hold on

end

 

Рис.11. Графики зависимости возмущения при расчете по
схеме (53) на правой границе от шага сетки h

 

Разностная схема (51¢), (52¢) устойчива[4], если решение системы разностных уравнений непрерывно зависит от входных данных j, c и эта зависимость равномерна относительно шага сетки. Уточним непрерывную зависимость. Это означает, что для любого e > 0 найдется такое d(e), не зависящее от h, что

, (56)

когда

. (57)

Если разностная схема (51¢), (52¢) линейна, то разностное решение линейно зависит от входных данных. В этом случае можно положить, что d(e) = e /(M + M1), где M, M1 — некоторые неотрицательные величины, независящие от h. В итоге условие устойчивости для линейных разностных схем можно записать в виде:

. (58)

Непрерывную зависимость разностного решения от j называют устойчивостью по правой части, а от cустойчивостью по граничным данным.

В дальнейшем будем рассматривать двуслойные разностные схемы, т.е. такие схемы которые содержат один известный и один новый, неизвестный слой.

Двуслойная разностная схема называется равномерно устойчивой по начальным данным, если при выборе начальных данных с любого слоя t* (t0 £ t* < T) разностная схема устойчива по ним, причем устойчивость равномерна по t*. Для линейных схем условие равномерной устойчивости можно записать в виде

, (59)

где константа K не зависит от t* и h, — решения разностной схемы Ahy = j с начальными данными и с одной и той же правой частью.

Достаточный признак равномерной устойчивости. Для равномерной устойчивости по начальным данным достаточно, чтобы при всех m выполнялось

. (60)

Доказательство. Условие (60) означает, что если на некотором слое возникла ошибка dy, то при переходе на следующий слой норма возмущения ||dy|| возрастает не более чем в (1 + Сt) £ eCt раз. Согласно (59), при переходе со слоя t* на слой t требуется m = (t - t*)/t шагов по времени, т.е. ошибка возрастает не более чем в . В итоге имеем

,

что, согласно определению в (59), означает равномерную устойчивость по начальным данным.

Теорема.Пусть двухслойная разностная схема Ahy = j равномерно устойчива по начальным данным и такова, что если два разностных решения Ahyk = jk равны на некотором слое, т.е. , то на следующем слое выполняется соотношение

, (61)

где a = const. Тогда разностная схема устойчива по правой части.

Доказательство. Помимо решения y рассмотрим решение , соответствующее возмущенной правой части . В дальнейшем будем считать, что . Это можно предположить, т.к. исследуется устойчивость по правой части.

Определим последовательность сеточных функций wm(t) при t ³ tm -1 согласно условиям:

(62)

Функции wm, m = 0,1,2,… определены так, что w0(t) º y(t) и при . Функции wm(t) и wm+1(t) на слое tm совпадают по определению в (62). С учетом (61), (62) имеем

.

При t ³ tm +1 функции wm(t) и wm+1(t) удовлетворяют разностной схеме с одной и той же правой частью j, но с разными начальными данными на слое tm +1. В силу равномерной устойчивости исходной разностной схемы по начальным данным можно сделать следующую оценку на последнем временном слое tM:

.

Далее воспользуемся неравенством треугольника

Последняя цепочка неравенств доказывает утверждение теоремы об устойчивости по правой части.

В теории разностных схем рассматривается несколько способов исследования устойчивости:

· принцип максимума;

· метод разделения переменных;

· метод операторных неравенств

и некоторые другие.

Начнем с принципа максимума. Запишем двуслойную схему в следующем виде:

, (63)

где суммирование на каждом слое производится в пределах шаблона около n-го узла. Считаем, что коэффициенты ak таковы, что . В этом случае:

a) схема равномерно устойчива по начальным данным, когда

, C = const > 0; (64)

b) схема устойчива по правой части, если верно (64) и

, k = const > 0. (65)

Доказательство. Докажем первую часть утверждения a). Фиксируем правую часть jn в (63) и возмущаем решение dy на исходном слое. В этом случае ошибка на следующем слое удовлетворяет уравнению

,

т.е.

.

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

,

или в другой форме

. (66)

Согласно (64), имеем

. (67)

Комбинируя (66), (67), получаем

,

что соответствует обеспечению достаточного признака устойчивости по начальным данным (60). Первая часть утверждения доказана.

Докажем вторую часть утверждения b). Возмутим в (63) правую часть, оставляя неизменным решение на нижнем слое, тогда

.

Из последнего равенства можно записать следующее неравенство

.

Аналогично предыдущей части доказательства, выберем узел , в котором возмущение на следующем шаге максимально и заменим соответствующие величины своими максимумами, тогда

.

Учитывая теперь (65), получаем

,

что, согласно (61), означает устойчивость по правой части. Вторая часть утверждения доказана.

Рассмотрим пример решения нестационарной краевой задачи для уравнения теплопроводности с постоянным коэффициентом теплопроводности:

(68)

Запишем неявную разностную схему (24) для задачи (68) на равномерной сетке:

(69)

Если переписать (69) в форме (63), то

(70)

Остальные коэффициенты в (70) равны нулю. Непосредственно можно проверить, что при любых соотношениях шагов t и h условие (65) выполнено во всех регулярных узлах, а условие (64) — во всех узлах сетки. Это означает, что схема (69) безусловно устойчива по начальным данным, правой части и краевым условиям.

Метод разделения переменных

Метод разделения переменных широко используется для обоснования устойчивости многих линейных разностных схем. При помощи этого метода устанавливается устойчивость в норме .

Рассмотрим применение метода разделения переменных на примере явной разностной схемы решения уравнения теплопроводности (26), заданной на равномерной сетке 0 = x0 < x1 < … < xN = a:

. (71)

Погрешность решения уравнения (71) удовлетворяет такому же уравнению, т.е.

. (71¢)

Решение уравнения (71¢) будем искать методом разделения переменных, т.е. в виде

. (72)

Подставляя (72) в (71¢), находим

. (73)

Можно сформулировать следующий признак устойчивости. Схема (71) с постоянными коэффициентами устойчива по начальным данным, если для всех q выполняется неравенство

, C = const. (74)

Доказательство. Система функций , q = 0,1,…,N - 1 полна и ортогональна на равномерной сетке xn, n = 0,1,…,N. Разложим произвольную ошибку начальных данных в ряд Фурье по приведенной выше системе функций, тогда

.

В силу линейности уравнения (71¢) и представления (72), можно записать

.

Используя ортогональность гармоник, находим

(75)

Учитывая (74), (75), получаем выражение

,

которое совпадает с признаком устойчивости по начальным данным (60), т.е. утверждение доказано.

Применим признак устойчивости по начальным данным (74) к оценке (73). При C = 0 можно найти, что при любых q = 0,1,…,N - 1, когда , т.е. схема (71) является условно устойчивой.