Устойчивость
Под устойчивостью (неустойчивостью) разностной схемы понимается то, что малые ошибки, возникающие в процессе счета (или внесенные с входными данными), при последующих расчетах уменьшаются (возрастают).
Рассмотрим пример неустойчивой разностной схемы для задачи Коши дифференциального уравнения 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) является условно устойчивой.