§ 8. Методы решения сеточных эллиптических уравнений

В этом параграфе будут рассмотрены методы решения систем сеточных эллиптических уравнений. Рассмотрим простейший пример (см. § 6). Пусть

система сеточных уравнений относительно неизвестных , получающаяся в результате аппроксимации краевой задачи Дирихле уравнения Пуассона в квадрате . Будем для простоты считать, что . Для решения системы уравнений (1) можно предложить большое число методов. Поэтому для того, чтобы их сравнивать, необходимо выбрать один или несколько критериев, по которым будет проводиться сравнение. Условимся критерием качества метода считать количество арифметических операций, которые необходимы либо для получения точного решения, либо для получения решения с некоторой заданной точностью. Так как часто не удается вычислить точное число арифметических операций либо этот подсчет достаточно громоздок, оценивают лишь порядок числа арифметических операций по отношению к числу узлов сетки.

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

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

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

Выясним, какое потребуется количество арифметических операций, если вычисления проводить только над ненулевыми элементами. Напомним, что матрица А системы уравнений (1) при «естественной» нумерации компонент вектора неизвестных имеет вид (с точностью до множителя )

где

а — диагональные матрицы с элементами на диагонали, равными — 1. Матрицы имеют размерность .

Таким образом, матрица А оказывается блочно-трехдиагональной и в каждой строке матрицы не более пяти элементов отличны от нуля; кроме этого, матрица А является ленточной с шириной ленты равной : все элементы при . Задача решения системы (1) является частным случаем решения системы m уравнений с неизвестными

с -диагоиальной матрицей.

Для решения такой системы могут быть применены, например, методы Гаусса, квадратного корня, отражений и вращений с исключением несодержательных операций во всех случаях. Под исключением несодержательных операций имеется в виду следующее. Поскольку при , то не надо проводить операций по обнулению этих элементов. Поэтому производятся операции лишь по обнулению элементов. Кроме того, при реализации каждого шага обнуления все элементы при остаются равными нулю и поэтому каждый шаг требует арифметических операций. Таким образом, общее число операций при решении системы уравнений (2) оказывается величиной порядка . В данном случае и поэтому общее число арифметических операций . Матрица А симметричная и положительно определенная. Как отмечалось ранее, в случае при реализации этих методов не возникает операции деления на нуль.

Еще один возможный метод решения этой системы — это метод прогонки в блочной форме. Исходную систему уравнений записываем в виде

где — вектор с компонентами .

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

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

В конкретном случае системы (2) при можно предложить метод, требующий арифметических операций. Для этого следует применять метод Гаусса в векторной форме при следующем порядке исключения неизвестных: сначала исключаются векторы , т.е. все векторы , затем векторы , т.е. все векторы , затем векторы , т.е. все векторы и т.д. Оказывается, что при таком способе исключения многократно умножаются и обращаютя одни и те же матрицы. Если вместо этого запомнить произведения соответствующих матриц, то общее число операций окажется .

Описанные методы требуют в общем случае для запоминания ленты матрицы слов памяти ЭВМ.

Задача 1. Показать, что при нумерации неизвестных

количество арифметических операций в методе Гаусса по порядку также равно .

Задача 2. Для случая разностной аппроксимации уравнения Лапласа в произвольной двумерной области указать порядок исключения неизвестных, при котором решение системы получается за арифметических операций; здесь М — общее число узлов.

Рассмотрим прямой метод решения системы уравнений (1), основанный на использовании дискретного преобразования Фурье в случае прямоугольника. Здесь мы не будем предполагать, что . Преобразуем вначале систему уравнений (1) к операторной форме. Пусть — сеточная функция, определенная при , которая соответствует решению (1), т.е. . Аналогично, сеточная функция, соответствующая правой части (1). Через обозначим операторы

Тогда система уравнений (1) может быть представлена в операторной форме

Пусть для определенности . Заметим, что этого всегда можно достичь, умножая обе части (1) на соответствущий коэффициент. Функция и как функция переменного представима в виде дискретной суммы Фурье по переменной (см. § 4.3):

Аналогично представим в виде

Подставляя полученные разложения в (3), получим

Напомним, что функции являются собственными для оператора :

Поэтому выражение (5) может быть заменено эквивалентным

Функции являются ортогональными в скалярном произведении

Поэтому (6) представляет собой систему независимых уравнений

Система уравнений (7) может быть получена из (6) также умножением обеих частей (6) скалярно на . Выражение (7) при фиксированном j представляет собой систему линейных алгебраических уравнений с трехдиагональной матрицей относительно неизвестных , которая может быть решена, например, методом прогонки.

Таким образом, алгоритм решения задачи (1) заключается в следующем:

а) находим из при каждом , коэффициенты ;

б) при решаем методом прогонки систему уравнений (7). В результате получаем функции ;

в) из формулы (4) при вычисляем значения функции .

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

Рассмотренный метод даст решение намного быстрее, чем метод Гаусса. Однако этот метод применим лишь в случае, когда исходная область является прямоугольником, в то время как метод Гаусса применим и в случае областей общей формы.

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

Рассмотрим другие приближенные методы решения системы уравнений (1), допускающие обобщение на случай более общих, чем прямоугольник, областей. В основу этих методов положено то свойство системы уравнений (1), что результат применения матрицы системы к вектору вычисляется по простым формулам и не требует запоминания матрицы. Количество арифметических операций, затрачиваемое на вычисление результата применения матрицы к вектору, по порядку равно , т. е. пропорционально длине вектора. Ранее мы упоминали, что матрица А системы уравнений (1) симметрична и положительно определена и в рассматриваемом случае ее собственные значения лежат в промежутке

Поэтому для решения системы уравнений

можно применить, например, метод простой итерации

где — вектор начального приближения. Согласно § 6.6 параметр целесообразно выбрать из соотношения

При этом метод (8) сходится со скоростью геометрической прогрессии и показатель скорости сходимости метода равен

При имеем , поэтому

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

Задача 3. Показать, что при расчетные формулы (8) приобретают вид:

как в случае прямоугольника, так и в случае произвольной области.

Поэтому для выполнения неравенства достаточно выбрать так, чтобы выполнялось неравенство . Отсюда . Так как , то . При малых h имеем требование

Так как нас интересует случай малых , т. е. когда точность решения системы достаточно высока, то второй член в правой части (9) может быть отброшен. Поэтому в данном случае количество итераций по порядку должно быть равно .

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

Поэтому матрицу А запоминать не надо.

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

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

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

Задача 4. Показать, что для итерационного процесса с чебышевским набором параметров требуемое число операций .

Задача 5. Получить такую же оценку числа действий для оптимального линейного итерационного процесса.

Задача 6. Получить такую же оценку числа действий для трехслойного итерационного процесса (см. § 6.6, задача 3) с фиксированным итерационным параметром .

Задача 7. Показать, что средняя трудоемкость трехслойного итерационного процесса (см. § 6.6) с фиксированным параметром со может быть снижена вдвое за счет следующего обстоятельства. При нахождении при четном вычисляются и запоминаются только значения в точках с четной суммой , а при нечетном — в точках с нечетной суммой .

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

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

а трехслойные — как аппроксимацию уравнения

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

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

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

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

Метод Федоренко (называемый также многосеточным методом). К числу наиболее эффективных и употребляемых методов решения сеточных эллиптических задач (включая краевые задачи для систем уравнений Навье—Стокса) относится многосеточный метод, предложенный в шестидесятые годы. Сначала практическое использование этого метода носило эпизодический характер из-за неприспособленности существовавшего тогда программного обеспечения к использованию методов такого типа.

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

Таким образом, решение исходной задачи сводится к решению задачи с относительно более гладким решением. Решение задачи с гладким решением на сетке с шагом h близко к решению задачи на сетке с более крупным шагом, например с шагом . Предлагается решить уравнение для погрешности на более грубой сетке и затем, проинтериолировав на исходную сетку, получить существенно лучшее приближение к решению. Для приближенного решения задачи на сетке с шагом применяется аналогичная процедура с переходом к решению задачи на сетке с шагом и так далее.

Часто одна итерация на сетке с шагом h состоит в дву- или трехкратном применении описанной процедуры перехода к решению задачи на сетке с шагом .

Рассмотрим этот метод на примере краевой задачи

Соответствующая сеточная краевая задача имеет вид:

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

Опишем подробнее процедуру сведения решения на сетке с шагом h к решению задачи на сетке с шагом . Рассмотрим итерационный процесс

Вычитая это соотношение из равенства

получим уравнение относительно погрешности :

В качестве начального приближения возьмем ; тогда

Далее в рассуждениях об уменьшении величины погрешности мы имеем в виду уменьшение погрешности решения относительно погрешности приближения .

Функции образуют полную ортонормированную систему в пространстве функций относительно скалярного произведения

и являются собственными функциями для оператора :

Поэтому

Отсюда видно, что итерационный процесс (10) сходится при . Если , где не зависит от h, то составляющие погрешности, соответствующие функциям со значениями q порядка N, то есть сильно колеблющиеся, умножаются на множители существенно меньшие единицы. Вследствие этого и происходит определенное сглаживание погрешности.

Положим далее , тогда проводимые выкладки имеют наиболее простой вид. В частности,

и соотношение (11) приобретает вид

Обозначим . Таким образом,

Определим оператор перехода с сетки с шагом h, на сетку с шагом :

и оператор перехода с сетки с шагом на сетку с шагом :

Лемма. Справедливы неравенства

(12)

Доказательство. Запишем неравенство Коши—Буняковского для скалярного произведения векторов :

Поэтому

Если , то после суммирования по четным и домножения на , получаем первое из утверждений леммы.

Аналогично имеем

Просуммируем по нечетным и добавим сумму . После умножения на h получим второе утверждение леммы. Если бы удалось решить точно уравнение

(13)

то, прибавив к , мы получили бы точное решение задачи .

Предположим, что известен алгоритм приближенного решения уравнения

такой, что приближенное решение может быть записано в виде

где и — некоторые линейные операторы и

Перенесем правую часть (13) на сетку с шагом и применим этот алгоритм к уравнению

Точное решение этого уравнения записывается в виде

Проинтерполировав получившееся приближенное решение на сетку с шагом h, получим приближение

положим

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

где

Функции являются собственными для операторов :

Справедливо равенство

Поэтому

Совпадение функции и в узлах сетки с шагом носит случайный характер и не имеет места в других случаях.

Справедливы равенства

Поэтому

Согласно равенству Парсеваля

Из неравенства Коши—Буняковского следует

Вследствие имеем

Поэтому

После суммирования по q получаем

Отсюда и из (12), (14) получаем оценку

Исходя из равенств

и равенств (16), подберем так, чтобы при всех выполнялось равенство

Получим уравнения

и следовательно, . Поэтому

Подставляя это разложение в (15), получим равенство

Согласно неравенству Коши—Буняковского

Первый множитель записывается в виде

Если ввести обозначение , то при всех q будет справедливо неравенство

и поэтому

Отсюда получаем оценку , где . Из этой оценки и оценки для следует, что

Здесь — ошибка начального приближения: .

Оценим величину . Справедливо равенство , где . Согласно формуле дифференцирования сложной функции

Отсюда получаем, что производная функции обращается в нуль в точках ; поэтому, исследуя график функции , получаем .

Задача 10. Показать, что

при .

Подведем итог проведенных построений. Начав с приближения , мы получили новое приближение , где ,

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

Далее фиксируем и возьмем . Из оценки (17) следует, что при использовании описанной процедуры норма погрешности решения на сетке с шагом h, умножится на множитель не больший 0,5.

После повторного применения этой процедуры норма погрешности решения на сетке с шагом h умножится на множитель не больший 0,25.

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

Обозначим число арифметических операций, достаточное для уменьшения погрешности в 4 раза на сетке с шагом через .

Полученный выше результат можно записать в виде неравенства

Функция удовлетворяет уравнению

При рассматриваемую сеточную задачу можно решить, совершив не более чем 3 арифметические операции. Поскольку , то .

Индукцией по можно получить оценку при .

Если требуется добиться уменьшения погрешности в М раз, то будет достаточно итераций, и, таким образом, общее число действий, требуемое для решения задачи, будет .

Эта оценка хуже, чем оценка числа действий методов прогонки или стрельбы.

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

В результате этого число операций, требуемое для уменьшения нормы погрешности на сетке с шагом h в 4 раза снижается до .

Можно предложить другой вариант уменьшения нормы погрешности, например, в 4 раза на сетке с шагом с затратой арифметических действий. Положим при . Применим описанный выше алгоритм последовательного сведения решения задачи на сетке с шагом к решению задачи на сетке с шагом при этом для каждого производим итераций по формуле (10) при . Задачу на сетке с шагом решаем точно.

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

При решении задачи существенно используется то, что при всех j выполняются наравенства .

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

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

В рассматриваемом случае при некотором видоизменении итерационного процесса (10) в результате одного полного цикла итерации — спуска от шага до шага и возвращения к шагу с общей затратой арифметических операций получается точное решение сеточной задачи. Однако этот факт носит случайный характер и не относится к более сложным задачам.

Рассмотрим случай размерности задачи . Пусть используется рассмотренный алгоритм сведения решения задачи на сетке с шагом h к решению задачи на сетке с шагом и число таких сведений при каждом шаге равно .

Для весьма широкого класса задач можно доказать существование , удовлетворяющих следующим соотношениям. После -кратного применения описанной выше процедуры, состоящей из итераций (10), перехода к задаче на сетке с шагом , решения задачи на этой сетке с помощью алгоритма и перехода на сетку с шагом , некоторая норма погрешности решения задачи на сетке с шагом h уменьшается в раз.

Естественно, что операторы перехода с одной сетки на другую , будут иными, чем выше.

Обозначим число арифметических операций, достаточное для уменьшения погрешности в раз на сетке с шагом , через . Справедливо неравенство .

Отсюда при можно получить неравенство .

Таким образом, общее число действий, достаточное для уменьшения погрешности на сетке с шагом h в раз, будет порядка , то есть порядка , где N - общее число узлов сеточной задачи.

Высказанные выше утверждения относятся и к случаю аппроксимаций метода конечных элементов.