1.1. Сеточная область
Для построения разностной схемы необходимо построить сетку
Gh-конечное множество точек, принадлежащих G, плотность
распределения которых характеризуется параметрами h — шагом сетки.
Пусть область изменения аргумента х есть отрезок. Разобьем этот
отрезок точками на
n равных частей длины каждая.
Множество точек называется равномерной сеткой на
отрезке и обозначим , а число h-расстояние между точками
(узлами) сетки называется шагом сетки. Разбиение отрезка
точками можно производить произвольным образом -. Тогда
получаем сетку с шагами , которое зависит от номера узла сетки.
Если хотя бы в одной точке, то сетка называется неравномерной
и такую сетку обозначают . Точки х0и хnназовем граничными узлами и
обозначим их. Остальные узлы назовем внутренними и обозначим их wh.
Узлы соседние с граничащими назовем приграничными. Тогда имеем
[4].
1.2. Сеточная функция. Пространство сеточных
функций. Нормы сеточных функций
Функция y=y(xi) дискретного аргумента хiназывается сеточной
функцией, определенной на сетке
. Сеточные функции можно
рассматривать как функции целочисленного аргумента, являющегося
номером узла сетки, т.е.. Далее мы будем писать .
Сеточная область whзависит от параметра h. При различных значениях
параметра hимеем различные сеточные области. Поэтому и сеточные
функции yh(x) зависят от параметра h.
Функции u(х) непрерывного аргумента являются элементами
функционального пространства H. Множество сеточных функций yh(x)
образует пространство Hh. Таким образом, в методе сеток
пространство Н, заменяется пространством Hhсеточных функций
yh(x).
Так как рассматривается множество сеток {wh}, то мы получаем
множество {Hh}пространств сеточных функций, определенных на
{wh}.
Пусть u(х) — решение исходной непрерывной задачи (1.1), uH;
yh-решение разностной задачи. yhHh. Для теории приближенных
вычислений представляет большой интерес оценка близости u(х) и
yh(x), но u(х) и yh(х) являются элементами из различных
пространств. Пространство Н отображается на пространство Hh.
Каждой функции u(х)Н ставится в соответствие сеточная функция
yh(x), х wh, так что yh=Phu Нh, где Ph— линейный
оператор из Н в Hh. Это соответствие можно осуществить различными
способами, т.е. зависит от выбора оператора Ph. Теперь, имея
сеточную функцию uh, образуем разность yh-uh, которая является
вектором пространства HhБлизость yhи uh, характеризуется числом ,
где — норма на Hh.
Соответствие функций u(х) и uhможно установить различными
способами, например,
uh=u(x), х wh.
В дальнейшем мы будем пользоваться этим способом соответствия.
В линейном пространстве Hhвведем норму , которая является аналогом
нормы || • ||н в исходном пространстве Н. Обычно принято выбирать
норму в пространстве Hhтак, чтобы при стремлении к нулю hона
переходила в ту или иную норму функций, заданных на всем отрезке,
т.е. чтобы выполнялось условие
,
(1.2)
где — норма в пространстве функций, определенных на отрезке,
которому принадлежит решение.
Условие (1.2) называют условием согласования в пространствах
Hh, и Н.
Рассмотрим простейшие типы норм в Hhдля случая сеток wh={xi=ih} на
отрезке .
1.Норма
удовлетворяет условию (1.2), если в качестве Н рассматривать
пространство непрерывных функций с нормой
а сеточную функцию определять в виде (1.2), т.е.
2. Норма
удовлетворяет условию (1.2), если за Н принять пространство
непрерывных функций с нормой
а сеточную функцию определять в виде
[4, 6].
продолжение
1.3. Аппроксимация дифференциальных операторов
Пусть имеем дифференциальный оператор . Этот оператор можно
аппроксимировать несколькими способами. Например,
— правая разностная
производная;
(1.3)
— левая разностная
производная;
(1.4)
— центральная разностная производная.
(1.5)
Можно взять их линейную комбинацию
,
(1.6)
где — вещественный параметр.
При =1 из (1.6) получаем аппроксимацию (1.3); при =0 —
аппроксимацию (1.4), а при =0.5 — аппроксимацию (1.5).
Чтобы показать погрешность аппроксимации, разложим по формуле
Тейлора
,
(1.7)
предполагая, что функция v(x) достаточно гладкая в некоторой
окрестности
(х- ho,x+ h) точки х, h
Подставляя это разложение в (1.3), (1.4), (1.5), получим:
Отсюда видно, что
Пусть L— дифференциальный оператор, Lh— разностный оператор,
заданный на сетке wh. Говорят, что разностный оператор Lh:
1) аппроксимирует
дифференциальный оператор Lв узле, если, где v(x) — достаточно
гладкая функция, стремится к нулю при ;
2) аппроксимирует Lс порядком
n>0 в узле если , т.е.
, M= const>0.
В качестве следующего примера рассмотрим оператор Lv= v"(x).
Для аппроксимации этого
оператора используем трехточечный шаблон
(x-h, х, x+h).
Замечая, имеем
Отсюда
Пользуясь разложением (1.7), покажем, что порядок аппроксимации
равен двум, т.е.
так как [1].
продолжение
1.4. Разностная схема
Как правило, дифференциальное уравнение решается с некоторыми
дополнительными условиями — начальными (задача Коши), краевыми
(краевая задача) либо и с начальными, и с краевыми
условиями (смешанные задачи). Эти дополнительные условия при
переходе к разностным уравнениям надо также аппроксимировать.
Пусть имеем некоторую дифференциальную задачу, записанную в
виде
Lu= f(x), xG
(1.8)
с дополнительным условием
lu= (х), хГ.
(1.9)
Введем в области сетку
и поставим в соответствие задаче (1.8), (1.9) разностную задачу
Lhyh=fh, xwh,
(1.10)
, xyh.
(1.11)
Функции Уh(х), fh(x), зависят от шага сетки. Меняя h, получаем
множества функций , зависящих от параметра h. Таким образом, мы
рассматриваем не одну разностную задачу, а семейство задач,
зависящее от параметра h. Это семейство задач называется разностной
схемой.
Рассмотрим примеры разностных схем, аппроксимирующих
дифференциальные задачи.
Пример 1. Имеем задачу Коши
Используем аппроксимации:
После этого имеем разностную схему:
Расчетный алгоритм имеет вид
Пример2. Рассмотрим задачу Коши
Воспользуемся следующими аппроксимациями:
После этого имеем разностную схему
[4].
1.5. Корректность разностной схемы
Пусть имеем дифференциальную задачу
Lu= f(x),
xG,
(1.12)
lu= (х), хГ
(1.13)
и на сетке аппроксимируем ее разностной схемой
Lhyh=fh,
xwh,
(1.14)
,
xyh.
(1.15)
Задача (1.12), (1.13) поставлена корректно, если выполнены
условия:
1) задача однозначно разрешима при любых правых частях fH,
Н;
2) решение задачи непрерывно зависит от правых частей fH, Н,
т.е.
Аналогично определяется понятие корректности разностной схемы
(1.14), (1.15). Говорят, что разностная схема (1.14), (1.15)
корректна, если при всех достаточно малых |h|
1) решение уhразностной схемы существует и единственно для всех
входных данных fhHh, Hh;
2) существуют постоянные М1>0, М2>0 не зависящие от hи такие,
что
при любых fhHh,Hhсправедлива оценка
(1.16)
Свойство 2), означающее непрерывную зависимость, равномерную
относительно h, решения разностной схемы от правых частей,
называется устойчивостью разностной схемы. Рассмотрим примеры.
Пример 1. Пусть имеем задачу:
(1.17)
Точным решением задачи (1.17) является функция
Если ввести новую функцию , то получим задачу
(1.18)
Решением задачи (1.18) является функция
Задачу (1.18) аппроксимируем на равномерной сетке схемой:
(1.19)
Перепишем схему (1.19) в виде
Отсюда имеем
Рассмотрим фиксированную точку и выберем последовательность
сеток таких, чтобы ,т.е. является узлом сетки
при
Вычислим значение у в этой точке .Taк как при и любых
h, то при любом h. Из этого неравенства видно, что решение
разностной схемы (1.19) непрерывно зависит от входных данных. В
таких случаях говорят, что разностная схема устойчива по входным
данным (по начальным условиям и по правой части).
Пример2. Имеем уравнение
.
(1.20)
Точным решением задачи (1.20) является функция
.
Отсюда следует неравенство
|u(x)||uo|
(1.21)
при .
Для устойчивости вычислительных алгоритмов
решения задачи (1.20) должно быть выполнено условие вида
(1.21), т.е.
(1.22)
Задачу (1.20) аппроксимируем явной схемой Эйлера
(1.23)
y(0)=u(0).
Выражая решение схемы (1.23) через начальное условие, имеем
Неравенство (1.22) будет выполнено, если
т.е. .
Таким образом, явная схема Эйлера условно устойчива.
Пример 3. Для численного
решения задачи (1.20) используем неявную
схему Эйлера
Отсюда
т.е.
при .
Схема (1.24) абсолютно устойчива, ибо выполнено условие (1.22) при
любом h.
Пример 4. Задачу (1.20) аппроксимируем схемой с весом
(1.25)
Отсюда имеем
Условие (1.22) будет выполнено, если
т.е.
Отсюда получаем
Схема абсолютно устойчива при
и
т.е. схема (7.25) условно устойчива при [4].
продолжение
1.6. Аппроксимация и сходимость
Для того, чтобы выяснить, с какой точностью приблизили функцию
u=u(x) с помощью функции у(х), мы должны их сравнить. Пусть
uhзначение функции u(х) на сеточной области , т.е. uhHh.
Рассмотрим погрешность решения разностной схемы (1.14), (1.15),
которая аппроксимирует на сетке дифференциальную задачу
(1.12), (1.13).
Введем функцию погрешности решения
zh=yh-uh,
где yh— решение схемы (1.14), (1.15), uh— решение задачи (1.12),
(1.13) на сетке . Подставив yh= zh+ uhв линейную задачу (1.14),
(1.15), получим для zhзадачу того же вида, что и (1.14),
(1.15):
(1.26)
(1.27)
где
(1.28)
Функции (1.28) называются погрешностью аппроксимации задачи (1.12),
(1.13), схемой (1.14), (1.15) на решении задачи (1.12), (1.13).
Будем говорить, что решение разностной схемы (1.14), (1.15)
сходится к решению задачи (1.12), (1.13), если
при
Разностная схема сходится со скоростью O(hh) или имеет n-ый порядок
точности, если при достаточно малом hhвыполняется неравенство
где М>0, не зависит от h, n>0.
Говорят, что разностная схема имеет n-ый порядок аппроксимации,
Если
т.е.
Теорема Лакса. Пусть дифференциальная задача (1.12), (1.13)
поставлена корректно, разностная схема (1.14), (1.15) является
корректной и аппроксимирует исходную задачу (1.12), (1.13). Тогда
решение разностной схемы (1.14), (1.15) сходится к решению
исходной задачи (1.12), (1.13), причем порядок точности совпадает с
порядком аппроксимации.
Доказательство. Если схема (1.14), (1.15) корректна, то нетрудно
получить оценку погрешности решения через погрешность
аппроксимации (1.28).
Задача (1.26), (1.27) аналогична задаче (1.14), (1.15), поэтому для
нее пользуясь априорной оценкой вида (1.16), получим оценку
(1.29)
Таким образом, если схема (1.14), (1.15) корректна и
аппроксимирует задачу (1.12), (1.13), то она сходится при h0.
Норма погрешности | при
h0, если и при h0.
Из оценки (1.28) видно, что порядок точности схемы (1.14),
(1.15) определяется порядком аппроксимации, и чтобы схема сходилась
со скоростью O(hn), n>0 достаточно, чтобы она имела
аппроксимацию того же порядка, т.е.
Рассмотрим примеры.
Пример 1. Рассмотрим явную схему Эйлера
которая аппроксимирует дифференциальную задачу (1.20). Покажем
порядок погрешности аппроксимации и сходимость.
Рассмотрим функцию погрешности решения
Для получаем схему:
(1.30)
Разложим по формуле Тейлора в точке , имеем
(1.31)
Подставляя (1.31) в , получим т.е. имеем первый порядок
аппроксимации. Из (1.30) имеем
При имеем Выражая через , получим:
Отсюда видно, что при .
Для точности схемы имеем
т.е. схема имеет первый порядок точности.
Пример 2. Рассмотрим неявную схему Эйлера
которая аппроксимирует дифференциальную задачу (1.20). Для
погрешности решения получаем разностную схему:
Подставляя разложение (1.31) в , получим
Отсюда имеем
т.е. первый порядок аппроксимации. Для сходимости рассмотрим
решение задачи для
Множитель при . Выражая через , имеем
Отсюда , т.е. схема имеет первый порядок точности. Таким же образом
можно показать, что схема с весом
имеет первый порядок аппроксимации и при выполнении условий
устойчивости имеет место сходимость и притом порядок точности
совпадает с порядком погрешности аппроксимации[4].
продолжение
§2. Основные понятия и история вопроса экономичных разностных
схем
1. Одним из основных критериев оптимальности в теории численных
методов является требование минимума арифметических операций.
Для одномерных задач математической физики особых затруднений в
этом плане нет. Разностные схемы в этом случае реализуются
экономичным алгоритмом прогонки, который на слое сетки допускает
количество арифметических операций пропорциональное количеству
узлов.
Особую остроту приобретает вопрос об экономичности вычислительных
алгоритмов в численном решении многомерных задач математической
физики. Многомерные краевые задачи моделировать многомерными
разностными схемами и решать непосредственно эти схемы
нецелесообразно, так как алгоритм становится сложным,
неэкономичным и нереализуемым на ЭВМ.
Пусть в цилиндре ищется решение уравнения
(2.1)
удовлетворяющее условиям
(2.2)
(2.3)
где -мерный единичный куб, — боковая поверхность .
Задачу (2.1)-(2.3) моделируем разностной схемой с
весом: (2.4)
где — вещественный параметр.
Разностное уравнение запишем в
операторном
где Е — единичный_ оператор.
Нахождение при требует обращения многомерного
оператора, что связано с весьма трудоемкой вычислительной работой,
поскольку соответствующая этому оператору матрица порядка не
имеет специального вида. Обращение такой матрицы производится по
методу Гаусса, что требуетарифметических операций или методом
итераций, который также требует много машинного времени. По этой
причине непосредственное решение схемы (2.4)-(2.6) даже при
нецелесообразно, так, как возникают серьезные проблемы
памяти ЭВМ и резко возрастает количество арифметических операций на
слое сеткой. Схема может оказаться нереализуемой на современных
ЭВМ.
Если в уравнении (2.4) , то получаем
явную схему
которую перепишем в виде
В этом случае обращается диагональная матрица, требующая лишь
0(N) арифметических операций.
Условием устойчивости схемы (2.4)-(2.6) будет
Отсюда видно, что условием устойчивости явной схемы будет
Из условия (2.7) видно, что схемы о , в том числе чисто неявная
схема при и симметричная схема при , абсолютно устойчивы.
Если имеем уравнение с переменными коэффициентами, т.е.
(2.8)
то
Условие устойчивости схемы (2.4)-(2.6) с эллиптическим оператором
(2.8) имеет вид
Отсюда получаем условие устойчивости явной схемы в виде
(2.9)
Условие устойчивости (2.9) налагает весьма жесткие ограничения на
шаг по времени и выход на заданный момент времени по
явной схеме требует неоправданно большого числа временных
шагов.
Из условия (2.9) видно, что выбор шага существенно зависит
от количества пространственных переменных р
и от величины М. С ростом числа р и для
быстро меняющегося коэффициента теплопроводности
K=K(x,t) шаг становится еще мельче.
В этом случае схема (2.4)-(2.6) при становится
неэффективной, так как выход на заданный момент времени
t=T требует слишком большого числа временных
шагов.
Чисто неявная схема () (2.4)-(2.6) абсолютно устойчива, т.е. на
параметры сетки h и
не налагаются ограничения. Поэтому счет можно
вести более крупными шагами hи и тем самым удается
значительно понизить порядок системы, уменьшить количество
временных шагов, необходимых для достижения момента времени t = Т
.
В неявных схемах придется решать системы алгебраических уравнений,
как правило, высокого порядка и с разреженными матрицами. Итак,
реализация систем алгебраических уравнений для многомерных задач в
общем случае представляется невозможной из-за громоздкости порядка
системы (проблема памяти ЭВМ), большого объема арифметических
операций.
Таким образом, лучшими качествами явной и неявной схем являются
количество арифметических операций, равное 0(N) в явной
схеме, и абсолютная устойчивость неявной схемы, недостатками —
условная устойчивость явной схемы и большое количество
арифметических операций в неявной схеме. Отсюда ясно, что если
построим разностные схемы, сочетающие в себе лучшие качества
обычных явных и неявных разностных схем, то можно эффективно решать
многомерные задачи математической физики.
Итак, безусловная (абсолютная) устойчивость и независимость
количества арифметических операций, требуемые для вычисления
приближенного решения задачи в отдельной точке сетки от общего
количества узлов сеточной области, определяют класс экономичных
схем.
Уравнение (2.1) можно аппроксимировать по-другому, стремясь
упростить вычислительный процесс. Уравнение (2.1) аппроксимируем
разностным уравнением
В этой схеме аппроксимация по первому направлению неявная, а по
всем остальным — явная.
Реализация этой схемы идет по направлению , как в одномерном
случае методом линейной прогонки, затратой арифметических операций
порядка . Однако схемы (2.10), (2.5), (2.5) условно
устойчивы. Таким образом, из требований экономичности разностной
схемы первое условие не выполнено, а второе — выполнено. Поэтому
схемы типа (2.10), (2.5), (2.6) не относятся к классу
экономичных схем, они конструируются специальным образом.
2. Достижением вычислительной математики является разработка
экономичных методов решения многомерных краевых задач
математической физики. Первые экономичные схемы были схемами в
дробных шагах по времени t. Они предложены и обоснованы в 1955 г.
одновременно американскими учеными D.W. Peacemdn,
H.H.Rachford и J.Douglas.
Характерной особенностью экономичных схем этого периода является
то, что все они основывались на идее ведения дробных моментов
времени и поэтапном решении р задач в промежутках
. Решением исходной задачи будет решение последней
р -й задачи в момент времени . Решения остальных
(р-1) задач являются вспомогательными в моменты времени
. Таким образом, наряду с основной сеткой по рассматривается
вспомогательная сетка. Все экономичные схемы этого периода
относятся к схемам последовательного перехода по времени t. Следуя
Н. Н. Яненко, назовем их методом дробных шагов.
В 1965г. А. А. Самарский предложил и обосновал экономичный метод
без привлечения вспомогательной сетки , т.е. экономичная схема
конструируется на исходной сетке . Таким образом, характерной
областью экономичных схем этого периода является то, что сетка
не вводится, вспомогательные функции рассматриваются на
верхнем слое. Составная схема конструируется на исходной сетке .
Как будет показано ниже, отказ от вспомогательной сетки дает
широкие возможности для построения различных экономичных схем и тем
самым удается значительно расширить их класс. Экономичные методы,
построенные на такой основе, назовем методом целых шагов. На основе
метода целых шагов строятся экономичные схемы последовательного и
параллельного (одновременного) перехода с нижнего слоя на верхний
слой по t.
Любая разностная схема, моделирующая исходную дифференциальную
краевую задачу, должна удовлетворять данным условиям устойчивости,
аппроксимации на решение исходной задачи и простоты. Если эти
требования для одномерной разностной схемы выполняются сравнительно
легко, то при переходе к двумерной (многомерной) задаче
возникают значительные трудности.
Перед вычислительной математикой встала сложная задача построения
экономичных методов решения многомерных задач математической
физики. Стало невозможным разрешить эту проблему на базе однородных
разностных схем, где при переходе от одного временного слоя к
другому одновременно удовлетворяются условия устойчивости и
аппроксимации. При этом, конечно, формула получается более
простой, но схема становится менее гибкой и имеет в своем
распоряжении небольшое количество произвольных параметров, что
создает основные трудности выполнения вышеописанных требований.
Экономичные схемы, расчленяя переход от нижнего слоя к верхнему на
ряд промежуточных этапов и не требуя на каждом этапе обязательного
выполнения свойств аппроксимации исходного уравнения и
устойчивости, имеют в своем распоряжении набор параметров, что
дает возможность выбора наиболее эффективного вычислительного
алгоритма.
Пусть имеем уравнение теплопроводности
Разрабатывая в 1955г. первые экономичные схемы переменных
направления, Писмен, Рэкфорд и Дуглас имели в виду упростить
решение алгебраической системы уравнений высокого порядка,
сохранить абсолютную устойчивость и приемлемую точность, тем самым
удовлетворили до некоторой степени вышеописанные требования. Идея
метода заключается в следующем. Уравнение (2.2) аппроксимируется
разностным уравнением
(2.12)
Нетрудно заметить, что схема (2.12) является условно устойчивой и
симметризуется так, чтобы и
менялись ролями от шага к шагу
(2.13)
Если в схему (2.23) счет повторяется при переходе
с j–го слоя на (j+2)-й, то (j+1)-й слоя выступает
в роли вспомогательного слоя.
Поэтому, введя дробные моменты времени , схему (2.13) можно
переписать в виде
(2.14)
Путем исключения вспомогательных функций в дробные моменты времени
получаем эквивалентную ей однородную схему
(2.15)
Схема (2.15) и эквивалентная ей схема (2.14) аппроксимируют
уравнение (2.11), как и схема
.
Безусловная устойчивость схемы (2.15) или (2.14) устанавливается
методом Неймана. Схема (2.14) решается двумя одномерными
прогонками, что не только упрощает алгоритм счета, но и уменьшает
объем вычислений. Из (2.13) или (2.14) видно, что методы Писмена,
Рэкфорда и Дугласа предполагают обязательную аппроксимацию
дифференциального уравнения по каждому направлению. Они
представляют неявный метод переменных направлений.
Заметим, что этот метод не годится для пространства нечетного
числа измерений. Например, при р=3
интегрирование в каждом направлении происходит один раз
неявно, а явно. Тогда возрастание ошибки в явной схеме не
компенсируется убыванием ее в неявной схеме. Это подтверждается
точным анализом устойчивости при р=3. Заметим также, что метод
переменных направлений не годится для уравнений со смешанными
производными даже при р=2.
Для решения трехмерного уравнения (2.2) в работе J.Douglas,
H.Rachfozd была предложена следующая
схема:
§3 Схема расщепления с последовательным переходом3.1.Постановка
задачи.
Требуется найти функцию u(x1,x2,t), удовлетворяющую следующим
условиям:
, 0
,
x = (x1, x2 )
u(x,0) = u0(x) , ut(x,0) = u1(x) ,
(2)
u(0, x2,t) = μ1- , u(l1, x2,t) =
μ1+
(3)
u( x1,0,t) = μ2- , u(x1, l2,t) = μ2+
l2
X2
II
I
l1
0
X1
3.2. Последовательный переход
(4)
(5)
Приближенным решением исходной задачи является решение задачи (5),
т.е.
продолжение
3.3. Разностная схема.
Задача (4) аппроксимируется схемой
(6)
(7)
Суммарная погрешность аппроксимации
Задача (6) решается методом правой прогонки по направлению .
Решением задачи на каждом слое будет решение второй схемы (7)
Задача (6) – вспомогательная.
3.4. Стандартный вид разностных схем.
Схемы (6) и (7) сведем к стандартному виду
,
(*)
где
по х1: ,1 , 0 ,
, ,
,
0 , 1 ,
по х2: , , ,
, ,
,
, ,
Построим сетку по времени и по пространству с шагами и h:
Метод правой прогонки
, ,
, ,
(**)
,
X2
x2
Схема(6)
Схема(7)
продолжение
0
x1
0
x1
§4 Схема расщепления с параллельным переходом. 4.1.Постановка
задачи.
Требуется найти функцию u(x1,x2,t), удовлетворяющую следующим
условиям:
, 0
,
x = (x1, x2 )
u(x,0) = u0(x) , ut(x,0) = u1(x) ,
(2)
u(0, x2,t) = μ1- , u(l1, x2,t) =
μ1+
(3)
u( x1,0,t) = μ2- , u(x1, l2,t) = μ2+
4.2. Параллельный переход на дифференциальном уровне.
Задачи (1) — (3) расщепляем на две автономные задачи по
направлениям осей координат.
(8)
(9)
Решение задачи определяется в виде
4.3. Параллельный переход на разностном уровне
(или
0)
(10)
(или
)
(11)
Схема (10) решается методом правой прогонки по направлению . Схема
(11) решается методом правой прогонки по направлению . Схемы (10) и
(11) свести к стандартному виду (*) и решить методом правой
прогонки (**).
, , , Т=1.
x2
l1
x2
Схема(10)
Схема(11)
продолжение