Московский физико-технический институт
Опубликован: 25.10.2007 | Доступ: свободный | Студентов: 3914 / 1195 | Оценка: 4.50 / 4.33 | Длительность: 24:00:00
ISBN: 978-5-9556-0065-9
Специальности: Программист, Математик
Лекция 9:

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

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

8.1. Базовые понятия

Рассмотрим численные методы решения задачи Коши для обыкновенных дифференциальных уравнений (ОДУ) вида

\begin{gather*}
 \frac{d u(t)}{d t} = f(t, u), \quad t > 0, \\ 
u(0) = u_0, 
\end{gather*}
 $ ( 8.1)

а также систем ОДУ

\begin{gather*}
\frac{d {\mathbf{u}}(t)}{d t} = {\mathbf{f}}(t, {\mathbf{u}}), \quad \\ 
{\mathbf{u}}(0) = {\mathbf{u}}_0, 
\end{gather*}

где

\mathbf{u} = {(u_1, u_2, \ldots , u_m)}^{T}, \mathbf{f} = {(f_1, f_2, \ldots , 
f_m)}^{T}

- векторы столбцы искомых функций и правых частей соответственно.

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

\begin{gather*}
\frac{d^m u}{d t^m} = g\left({t, u, \frac{d u}{d t}, \ldots , \frac{d^{m - 1} u}{d t^{m - 1}}}\right), \quad t > 0, \quad \\ 
u(0) = a_0, \quad \\ 
\frac{d u}{d t}(0) = a_1, \ldots , \frac{{d^{m - 1} u}}{{d t^{m - 1}}}(0) = a_{m - 1}, 
\end{gather*}

если положить

\begin{gather*}
u_1 = u, u_2 = \frac{d u_1}{d t}, u_3 = \frac{d u_2 }{d t}, \ldots , u_m = \frac{d u_{m - 1}}{d t}, \\ 
\frac{{d u_m}}{{d t}} = g(t, u_1, u_2, \ldots , u_m), \\ 
u_i (0) = a_{i - 1}, i = 1, 2, \ldots , m.
\end{gather*}

Введем в расчетной области t \in \left[{0, T}\right] точки (узлы расчетной сетки) \{ t_{n} = n\tau , n = 0, 1, \dots , N\}, в которых вычисляется искомое решение. Совокупность узлов называется расчетной сеткой, (сеточной областью), \tauшагом интегрирования. Здесь для простоты введена равномерная сетка. В реальных расчетах применяются и неравномерные сетки.

Введем сеточную функцию u^{\tau}, определенную в узлах сетки и представляющую собой совокупность приближенных значений искомой функции, U^{\tau} — проекцию точного решения искомой задачи на сетку и f^{\tau} — значения правой части в узлах сетки.

Введем вслед за [8.1] также операторное обозначение дифференциальной задачи

{L(u) = F, } ( 8.2)

где

$ L(u) = \left\{ \begin{array}{cc}
 {\frac{d u}{d t} - f(t, u), } & {t > 0;} \\
 {u(0), } & {t = 0;} \\
\end{array} \right. 
F = \left\{ \begin{array}{cc}
 {0, } & {t > 0;} \\
 {u_0, } & {t = 0;} \\
\end{array} \right. $

и аппроксимирующей разностной задачи

{L_{\tau}(u^{\tau}) = F_{\tau}, } ( 8.3)

где L_{\tau } — обозначения разностного оператора, F_{\tau }проекция F на расчетную сетку. Заметим, что u и u^{\tau } являются элементами соответственно функционального и конечномерного пространств. Определим основные понятия теории разностных схем [8.1], [8.2].

Определение. Решение задачи (8.3) u^{\tau } сходится при {\tau}\to 0 к решению исходной задачи (8.2), если

\left\|{u^{\tau}- U^{\tau}}\right\| \to 0

при {\tau}\to 0.

При этом, если имеет место оценка

\left\|{u^{\tau}- U^{\tau}}\right\| \le C {\tau}^{p}(C \ne C({\tau})),

то имеет место сходимость порядка p.

Определение. Говорят, что задача (8.3) аппроксимирует задачу (8.2) на ее решении, если невязка

\left\|{r_{\tau}}\right\| \to 0

при {\tau}\to 0, где r_{\tau}\equiv L_{\tau}(U^{\tau}) - F_{\tau} ; при этом, если имеет место оценка

\left\|{r_{\tau}}\right\| \le C_1 {\tau}^{p}(C_1 \ne C_1 ({\tau})),

то говорят, что имеет место аппроксимация порядка p.

Определение. Задача (8.3) устойчива, если из соотношений

L_{\tau }(u^{\tau }) - F_{\tau } = \xi _{\tau } , 
\\
L_{\tau }(v^{\tau }) - F_{\tau } = \eta _{\tau }

следует

\left\|{u^{\tau}- v^{\tau}}\right\| \le C_2 \left({\left\|{\xi_{\tau}}\right\| + 
\left\|{\eta_{\tau}}\right\|}\right), C_2 \ne C_2 ({\tau}).

Теорема 1 (В.С.Рябенького - П. Лакса). Решение задачи (8.3) сходится к решению исходной задачи (8.2), если задача (8.3) устойчива и аппроксимирует задачу (8.2); если аппроксимация имеет порядок p, то сходимость также имеет порядок p.

Доказательство.

В силу аппроксимации имеем оценку: \left\|{r_{\tau}}\right\| \le C_1 {\tau}^{p} . Тогда из определения устойчивости, положив v^{\tau } = U^{\tau }, получим

\left\|{u_{\tau}- U_{\tau}}\right\| \le C_2 \left\|{r_{\tau}}\right\| \le C_2 C_1 {\tau}^{p} = 
C{\tau}^{p},

поскольку в данном случае

|| \eta ^{\tau } || = 0

и, кроме того,

|| r_{\tau } || = || \xi _{\tau } ||.

Приведем примеры простейших разностных уравнений, аппроксимирующих (8.1):

\begin{gather*}
\frac{u_{n + 1} - u_n}{{\tau}} = f(t_n, u_n), 0 \le n \le N - 1, \\ 
\frac{u_{n + 1} - u_n}{{\tau}} = f(t_n, u_{n + 1}), 0 \le n \le N - 1, \\ 
\frac{u_{n + 1} - u_{n - 1}}{2{\tau}} = f(t_n, u_n), 1 \le n \le N - 1.
 \end{gather*}

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

Для реализации третьей схемы необходимо задание функции un в двух точках: t0 и t1 . Один из возможных вариантов — решение на первом шаге нелинейного уравнения вида:

$ 
\frac{u_{n + 1} - u_{n - 1}}{2{\tau}} = \frac{1}{2} [f(t_{n - 1}, u_{n - 1}) + f(t_{n + 1}, u_{n + 1})] $

при n = 1.

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

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

Дифференцируя по t исходное уравнение (8.1), получим

u'' = f't(t, u) + f't(t, u) u', 
u''' = f''tt(t, u) + 2f''tu(t, u) u' + f''uu(t, u) (u')2 + f'u(t, u) u'', ...