Конечно-разностные методы решения краевых задач

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

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

Процесс построения численного решения дифференциального уравнения можно представить в следующем виде:

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

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

(7.7)

При малых значениях , отбросив малые члены более высоких порядков относительно (пренебрегая ими) в разложении (7.7), можно получить выражения для аппроксимации производных. Погрешность аппроксимации определяется старшим членом в отброшенном остатке, который содержит , где определяет порядок точности аппроксимации. Для обозначения остаточного члена разложения используется специальная символика – «о малое от », которая определяется как «бесконечно малое относительно », т.е. , . Для обозначения остаточного члена разложения также используется символика – «О большое от », т.е. .

Для аппроксимации производных первого порядка в разложении (7.7) можно оставить два первых члена, т.е. линейную часть относительно :

.

Тогда конечно-разностные аналоги для аппроксимации производных первого порядка можно записать в виде:

, (7.8)

. (7.9)

Выражение (7.8) используется для аппроксимации правой производной, взятой с шагом вперед, выражение (7.9) – для левой производной, взятой с шагом назад. Определим порядок точности конечно-разностной аппроксимации производных первого порядка по формулам (7.8) и (7.9). Погрешность такой аппроксимации определяется старшим членом в отброшенном остатке. Заменяя в (7.8) значение разложением в ряд Тейлора (7.7), получим:

Таким образом, формула (7.8) имеет аппроксимацию первого порядка. Определим погрешность аппроксимации для левой производной по формуле (7.9):

Формула (7.9) для аппроксимации левой производной имеет первый порядок точности.

Для аппроксимации производных первого порядка также можно построить конечно-разностные аналоги более высокого порядка. Для этого разложение (7.7) запишем в виде:

(7.10)

(7.11)

Далее, если из разложения (7.10) почленно вычесть разложение (7.11), то из разложения разности исчезнут свободный член и члены четных степеней относительно :

Если в последнем разложении разности оставить линейную часть относительно , то можно получить более точную аппроксимацию первой производной, которая называется центральной разностью:

. (7.12)

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

Чтобы построить конечно-разностные аналоги для аппроксимации производных второго порядка, следует сложить разложение (7.10) почленно с разложением (7.11). Тогда из разложения суммы исчезнут члены нечетных степеней относительно :

Если в последнем разложении суммы оставить члены до второго порядка включительно относительно , то получим формулу для аппроксимации производной второго порядка ( ):

.(7.13)

Отброшенный остаток будет содержать член с во второй степени (после деления на ), т.е. формула (7.13) для аппроксимации второй производной имеет второй порядок точности.

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

Например:

, (7.14)

. (7.15)

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

,

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

получим следующее разложение:

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

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

(I)

(II)

(III)

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

Решение краевой задачи заключается в определении значения функции , удовлетворяющего уравнению (I) и краевым условиям (II)-(III). Основная идея численного решения краевой задачи (I)-(III) методом конечных разностей (МКР) заключается в следующем: в заданной области выбирается система узловых точек и задача считается решенной, если найдены значения искомой функции в узловых точках. При этом решение проводится в два этапа: 1) конечно-разностная аппроксимация производных дифференциального уравнения и краевых условий, т.е. составление конечно-разностного аналога задачи (дискретизация) – получается система уравнений; 2) решение полученной системы уравнений относительно значений искомой функции в узлах сетки.

На заданном отрезке выбирается сетка узлов: .

Для удобства расчетов можно выбрать и равномерную сетку: . Далее, дифференциальное уравнение (I)

рассматривается во внутренних узлах , где :

.

Заменяя первые и вторые производные во внутренних узлах центральными конечно-разностными аналогами (7.12) и (7.13), получим систему из -го уравнения с -им неизвестным относительно :

,

где .

Недостающие два уравнения получим из аппроксимации краевых условий (II) и (III). Производные первого порядка в краевых условиях заменяются конечно-разностными аналогами (7.8) и (7.9):

.

В итоге получим недостающие два уравнения:

и .

Таким образом, для определения -го неизвестного получили систему из -го линейного алгебраического уравнения:

(IV)

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

1. Использование формул одностороннего дифференцирования более высокого порядка точности.

2. Расширение сеточной области добавлением фиктивного узла.

3. Использование разложения искомой функции в ряд Тейлора в граничных точках.

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

Далее дифференциальное уравнение (I) и краевые условия (II) и (III) аппроксимируем в граничных точках и со вторым порядком. Замена первых и вторых производных дифференциального уравнения (I) в граничных точках и центральными конечно-разностными аналогами второго порядка (7.12) и (7.13) приводит к линейным уравнениям:

, (A)

. (B)

Аппроксимация краевых условий (II) и (III) центральной разностью второго порядка (7.12) в граничных точках и дает уравнения:

(С)

. (D)

Следующим шагом из уравнения (С) исключается значение сеточной функции в фиктивной точке , а из уравнения (D) – значение сеточной функции в фиктивной точке . Для этого из уравнения (A) выражается значение и подставляется это значение в уравнение (С). Аналогично, из уравнения (B) выражается значение и подставляется в уравнение (D). Теперь остается только первое и последнее уравнения системы (IV) заменить новыми уравнениями (С) и (D). В итоге построена разностная схема, которая имеет второй порядок аппроксимации во всех узловых точках сеточной области, включая и граничные.

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

Аналогичным образом, как и для функции одной переменной, используя разложение в ряд Тейлора функции многих переменных, можно получить конечно-разностные аппроксимации частных производных. Например, для функции двух переменных выражения для конечных разностей в предположении, что первый индекс относится к координате , а второй – , т.е. и принимая , , будут представлены в виде:

– центральная разность по оси ,

– центральная разность по оси ,

– вторая разность по оси , (7.16)

– вторая разность по оси ,

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

Конфигурацию узлов, используемую для составления разностной схемы (разностных алгебраических уравнений), называют шаблоном. Для решения одной и той же задачи можно использовать разные шаблоны. На рис.7.1 приводятся некоторые варианты вычислительных шаблонов для конечно-разностной аппроксимации частных производных функции двух переменных:

Рис. 7.1

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

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

2) Погрешность усечения (дискретизации) связана с тем, что для аппроксимации функции вместо бесконечных рядов используется лишь несколько первичных членов (не зависит от ЭВМ).

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

Указанные три источника погрешностей являются причиной наблюдаемых ошибок двух типов:

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

2) Глобальная ошибка – разность между вычисленным и точным решением, определяющая суммарную погрешность численного алгоритма, накопившуюся с момента начала вычислений.

Фундаментальным понятием при оценке точности численного метода является его порядок, который определяется в терминах локальной ошибки, получаемой при применении метода к задачам с гладким решением. Говорят, что численный метод имеет порядок , если существует число , такое, что для локальной ошибки выполняется условие , где – абсолютная величина локальной ошибки между теоретическим решением и приближенным решением, полученным с помощью численного метода на -м шаге, – не зависит от номера шага и величины шага , т.е. имеет место выражение . Расчетная формула численного метода согласуется (совпадает) с разложением решения в ряд Тейлора до членов порядка , тогда – порядок метода. Таким образом, порядок точности численного метода определяется порядком зависимости ошибки от величины шага. Например, метод Эйлера для решения обыкновенных дифференциальных уравнений (ОДУ) имеет первый порядок точности, так как для него зависимость ошибки от величины шага линейна, т.е. при уменьшении шага в раз ошибка также уменьшается в раз. Расчетная формула метода Рунге-Кутты для решения ОДУ содержит первые члены ряда Тейлора до четвертого порядка включительно относительно шага. Следовательно, он имеет четвертый порядок точности.

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