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

Решение задачи Коши для системы ОДУ методом Рунге-Кутта
Последовательный алгоритм
Последовательная сложность	$4 m n$
Объём входных данных	$n + 3$
Объём выходных данных	$m(n + 1)$
Параллельный алгоритм
Высота ярусно-параллельной формы	?
Ширина ярусно-параллельной формы	?

Основные авторы описания: Филимонова Юлия

1 Свойства и структура алгоритма

1.1 Общее описание алгоритма

Методы Рунге-Кутты (распространено неправильное название Методы Рунге-Кутта) — важное семейство численных алгоритмов решения обыкновенных дифференциальных уравнений и их систем. Данные итеративные методы явного и неявного приближённого вычисления были разработаны около 1900 года немецкими математиками К. Рунге и М. В. Куттой.

Формально, методом Рунге-Кутты является модифицированный и исправленный метод Эйлера, они представляют собой схемы второго порядка точности. Существуют стандартные схемы третьего порядка, не получившие широкого распространения. Наиболее часто используется и реализована в различных математических пакетах (Maple, MathCAD, Maxima) стандартная схема четвёртого порядка. Классический метод Рунге-Кутты четвёртого порядка столь широко распространён, что его часто называют просто методом Рунге-Кутты, опуская порядок.

1.2 Математическое описание алгоритма

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

$\dot{x} = f(t, x),\ t_0 \leqslant t \leqslant t_1,\ x(t_0) = x_0.$

Здесь $x(t), x_0 \in \mathbb{R}^n,\ t, t_0, t_1 \in \mathbb{R}, f: \mathbb{R} \rightarrow \mathbb{R}^n$.

Введем равномерную сетку

$t_i = t_0 + ih,\ i = \overline{1, n},\ h = \frac{t_1 - t_0}{m},$

$x(t_i) = x_i.$

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

$x_{i+1} = x_i + \frac{h}{6} (K_1 + 2 K_2 + 2 K_3 + K_4).$

Вычисление нового значения происходит в четыре стадии:

$K_1 = f (t_i, x_i),$

$K_2 = f (t_i + \frac{h}{2}, x_i + \frac{h}{2} K_1),$

$K_3 = f (t_i + \frac{h}{2}, x_i + \frac{h}{2} K_2),$

$K_4 = f (t_i + h, x_i + h K_3).$

1.3 Вычислительное ядро алгоритма

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

1.4 Макроструктура алгоритма

Макроструктура алгоритма представлена одним шагом итерационного процесса, описанного в пункте 1.2. Основной макрооперацией алгоритма является операция обращения к функции $f$, и основное внимание будет уделено распараллеливанию этой операции.

1.5 Схема реализации последовательного алгоритма

Приведем здесь псевдокод

начало;

ввод начальных параметров ($x_0, t_0, t_1, m$);

цикл по числу узлов сетки: $i = \overline{1,m-1}$

цикл по числу стадий: $j = \overline{1,4}$

цикл по числу компонент вектора $x:\ k = \overline{1,n}$

вычисление коэффициентов $K_{j, k}^{i}$;

конец цикла по $k$;

конец цикла по $j$;

цикл по числу компонент вектора $x:\ k = \overline{1,n}$

вычисление $y_{j}^{i+1}$;

конец цикла по $k$;

цикл по числу компонент вектора $x:\ k = \overline{1,n}$

сохранение $y_{j}^{i} = y_{j}^{i+1}$;

конец цикла по $k$;

вывод $y^{i+1}$;

конец цикла по $n$;

конец;

1.6 Последовательная сложность алгоритма

В данном алгоритме производится четыре операции обращения к функции $f$ (порядка $n m$ арифметических операций) и четырнадцать операций сложения векторов и умножения вектора на число (порядка $n$ арифметических операций). Поскольку операция обращения к функции является более сложной, то сложность последовательного алгоритма можно считать равной $4 n m$.

1.9 Входные и выходные данные алгоритма

На вход алгоритма подаются следующие данные:

1. вектор начальных значений $x_0$ размерности $n$;
2. границы временного интервала $t_0, t_1$;
3. частота дискретизации $n$.

Общий размер входных данных $n + 3$.

На выходе получаются следующие данные:

1. вектор времени $t$ размерности $m$;
2. матрица значений $x$ ($m$ векторов длины $n$) размерности $m n$.

Общий размер выходных данных $m + m n$.