Опубликован: 25.10.2007 | Уровень: специалист | Доступ: платный | ВУЗ: Московский физико-технический институт

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

< Лекция 10 || Лекция 11: 123456789

10.6.2. Метод квазилинеаризации (метод Ньютона)

Рассмотрим простейшую разностную аппроксимацию краевой задачи для ОДУ второго порядка

$ \frac{d^2u}{dt^2} = f(u), u(0) = U_1, u(L) = U_2  $

следующего вида:

$ \frac{u_{n - 1} - 2u_n + u_{n + 1}}{\tau^2} = f(u_n), n = 1, \ldots , N - 1, u_0 = U_1, u_n = U_2, \tau = L/N . $

Зададим некоторые начальные приближения u_n^0 = \varphi (t_n) к искомой функции и построим итерационный процесс, воспользовавшись линейным приближением функции правой части

\begin{gather*}
f(u_n^{i + 1})  \approx  f(u_n^{i}) + f^{\prime}_u (u_n^{i})(u_n^{i + 1} - u_n^{i})), \\ 
\frac{{u_{n - 1}^{i + 1} - 2u_n^{i + 1} + u_{n + 1}^{i + 1}}}{{\tau ^2 }} = f(u_n^{i}) + 
f^{\prime} (u_n^{i})(u_n^{i + 1} - u_n^{i}), u_n^0 = \varphi (t_n), i = 0, 1, \ldots
\end{gather*}

где i является итерационным индексом, а начальное приближение \varphi (t_n) удовлетворяет граничным условиям \varphi(0) = U_1, \varphi(L) = U_2.

На первой итерации, т.е. при i = 0, получаем первое приближение к искомому решению, решая методом трехточечной прогонки разностное уравнение

$ \frac{{u_{n - 1}^1 - 2u_n^1 + u_{n + 1}^1 }}{{\tau ^2 }} = f(u_n^0 ) + f^{\prime}_u (u_n^0 )(u_n^1 - u_n^0 ), u_n^0 = \varphi_n .  $

Полученное первое приближение вновь подставляем в линеаризованную разностную задачу и методом прогонки получаем второе приближение. Аналогично поступаем и далее, до достижения заданной точности в соответствии с условием \left\|u^{i + 1} - u^{i}\right\| \le \varepsilon.

10.6.3. Аппроксимация граничных условий

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

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

\begin{gather*}
\frac{d^2 u}{dt^2} + h\frac{du}{dt} + su = f, \\ 
u^{\prime}(0) = a, u^{\prime}(1) = b, 
\end{gather*}

и соответствующую ей разностную задачу

$ \frac{u_{n + 1} - 2u_n + u_{n - 1}}{\tau ^2 } + h\frac{u_{n + 1} - u_{n - 1}}{2\tau } + su}_n = f_n , $

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

u_{1} - u_{0} = \tau  a,  u_{n} - u_{ N - 1} = \tau  b.

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

Наиболее распространены три способа.

  1. Использование формул одностороннего дифференцирования более высокого порядка точности. Эти формулы разбирались в "Предмет вычислительной математики. Обусловленность задачи, устойчивость алгоритма, погрешности вычислений. Задача численного дифференцирования" . Подход очевиден, однако имеет недостаток — матрица системы уравнений для определения решения будет уже не трехдиагональной, следовательно, алгоритм прогонки неприменим. Можно избавиться от этого недостатка, проведя перед численным решением системы соответствующие алгебраические преобразования. Они в случае применения конкретной формулы численного дифференцирования для аппроксимации граничных условий свои, но достаточно очевидны.
  2. Использование фиктивной ячейки (фиктивного узла). Рассмотрим расширение сеточной области. Введем в рассмотрение узлы с индексами 0 и N + 1, находящиеся формально за пределами рассматриваемой области. Пусть в них тоже определено значение сеточной функции. Тогда узел с индексом 0 выступает как внутренний узел, и в нем может быть записано соотношение

    $ \frac{{u_1 - 2u_0 + u_{- 1}}}{{\tau ^2 }} + h\frac{{u_1 - 
u_{- 1}}}{{2\tau }} + su}_0 = f_0   $,

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

    $ \frac{{u_1 - u_{- 1}}}{\tau } = a  $.
    Из этих двух уравнений осталось только исключить значение в фиктивном узле (фиктивной ячейке).

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

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

  3. Использование ряда Тейлора. Для формулы на левой границы области интегрирования u_{1} - u_{0} = \tau  a в явном виде выпишем главный член невязки:

    $ \frac{{u_1 - u_0 }}{\tau } = u^{\prime} (0) + 
\frac{\tau }{2}u^{\prime\prime}(0) + o(\tau)  $.

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

    $ \frac{{d^2 u(0)}}{{dt^2 }} = f(0) - h\frac{{du(0)}}{dt} - 
su}(0), $

    откуда следует

    $ \frac{{u_1 - u_0 }}{\tau } = a + \frac{\tau }{2}(f(0) - 
h\frac{{du(0)}}{dt} - su}(0)) + o(\tau ).  $

    Заменив в последнем соотношении производную на конечную разность, получим

    $ \frac{{u_1 - u_0 }}{\tau } = a + \frac{\tau }{2}(f_0 - h\frac{{u_1 - 
u_0 }}{\tau } - su}_0 ) . $

    Осталось только привести это соотношение к виду, удобному для использования в методе прогонки.

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

< Лекция 10 || Лекция 11: 123456789
Эдуард Макаров
Эдуард Макаров
Россия, Челябинск, Челябинский политехнический институт, 1966
Иван Кузнецов
Иван Кузнецов
Россия, г. Новосибирск