Численные методы решения задачи Коши для систем обыкновенных дифференциальных уравнений
8.1. Базовые понятия
Рассмотрим численные методы решения задачи Коши для обыкновенных дифференциальных уравнений (ОДУ) вида
![]() |
( 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*}](/sites/default/files/tex_cache/12e77c30288c02f7487f80ae83b08984.png)
где
![\mathbf{u} = {(u_1, u_2, \ldots , u_m)}^{T}, \mathbf{f} = {(f_1, f_2, \ldots ,
f_m)}^{T}](/sites/default/files/tex_cache/8773a0de7e540c8d6d9456f961f07901.png)
- векторы столбцы искомых функций и правых частей соответственно.
К аналогичной форме приводится задача Коши для обыкновенного дифференциального уравнения (системы уравнений) порядка выше первого, вида
![\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*}](/sites/default/files/tex_cache/56cee5145121f44a506fa7df50038fda.png)
если положить
![\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*}](/sites/default/files/tex_cache/fb0ba04ce8e443c1208f88e0a78819bf.png)
Введем в расчетной области точки (узлы расчетной сетки)
в которых вычисляется искомое решение. Совокупность узлов называется расчетной сеткой, (сеточной областью),
— шагом интегрирования. Здесь для простоты введена равномерная сетка.
В реальных расчетах применяются и неравномерные сетки.
Введем сеточную функцию определенную в узлах сетки и
представляющую собой совокупность приближенных значений искомой функции,
— проекцию точного решения искомой задачи на сетку и
— значения правой части в узлах сетки.
Введем вслед за [8.1] также операторное обозначение дифференциальной задачи
![]() |
( 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. $](/sites/default/files/tex_cache/c251a7fbbe920746571c0ddd8cbc6210.png)
и аппроксимирующей разностной задачи
![]() |
( 8.3) |
где — обозначения разностного оператора,
— проекция F на
расчетную сетку. Заметим, что u и
являются элементами соответственно функционального и конечномерного пространств. Определим основные понятия теории разностных схем [8.1], [8.2].
Определение. Решение задачи (8.3) сходится при
к решению исходной задачи (8.2), если
![\left\|{u^{\tau}- U^{\tau}}\right\| \to 0](/sites/default/files/tex_cache/6354a1b4f059fa7854eeaaf9178db2db.png)
при
При этом, если имеет место оценка
![\left\|{u^{\tau}- U^{\tau}}\right\| \le C {\tau}^{p}(C \ne C({\tau})),](/sites/default/files/tex_cache/39d7c3f5939311f63de76c05c2d7922a.png)
то имеет место сходимость порядка p.
Определение. Говорят, что задача (8.3) аппроксимирует задачу (8.2) на ее решении, если невязка
![\left\|{r_{\tau}}\right\| \to 0](/sites/default/files/tex_cache/922bef7c01a5b600ead9ed1ecc48103c.png)
при где
; при этом, если имеет место оценка
![\left\|{r_{\tau}}\right\| \le C_1 {\tau}^{p}(C_1 \ne C_1 ({\tau})),](/sites/default/files/tex_cache/f7b3a30a3a4e3902868f0e9a4bb0faf0.png)
то говорят, что имеет место аппроксимация порядка p.
Определение. Задача (8.3) устойчива, если из соотношений
![L_{\tau }(u^{\tau }) - F_{\tau } = \xi _{\tau } ,
\\
L_{\tau }(v^{\tau }) - F_{\tau } = \eta _{\tau }](/sites/default/files/tex_cache/d6dd688a107b959ef5e2b1c4ae55785b.png)
следует
![\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}).](/sites/default/files/tex_cache/5cf83124ea2dcd27314095c1c6d57e4f.png)
Теорема 1 (В.С.Рябенького - П. Лакса). Решение задачи (8.3) сходится к решению исходной задачи (8.2), если задача (8.3) устойчива и аппроксимирует задачу (8.2); если аппроксимация имеет порядок p, то сходимость также имеет порядок p.
Доказательство.
В силу аппроксимации имеем оценку: Тогда из определения устойчивости, положив
получим
![\left\|{u_{\tau}- U_{\tau}}\right\| \le C_2 \left\|{r_{\tau}}\right\| \le C_2 C_1 {\tau}^{p} =
C{\tau}^{p},](/sites/default/files/tex_cache/931469ca3ca788c6e6e45daff601f1dc.png)
поскольку в данном случае
![|| \eta ^{\tau } || = 0](/sites/default/files/tex_cache/a804f4d89225d997241012346bb5bece.png)
и, кроме того,
![|| r_{\tau } || = || \xi _{\tau } ||.](/sites/default/files/tex_cache/5b5d84715d2c018b72875f540af71a25.png)
Приведем примеры простейших разностных уравнений, аппроксимирующих (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*}](/sites/default/files/tex_cache/6d81c4b5c2d9606d5e7fbbbb801f9348.png)
Первая из схем называется явной (явная схема Эйлера), вторая — неявной (неявная схема Эйлера). Алгоритмическая реализация первой схемы — бегущий счет (рекуррентная формула), второй — решение нелинейного алгебраического уравнения на каждом временном шаге.
Для реализации третьей схемы необходимо задание функции 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})] $](/sites/default/files/tex_cache/1e6c85a1c977f5b85507cbd1e814a2d3.png)
при 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'', ...