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

(2.35)

Считаем, что задача поставлена корректно [1]. На  задана равномерная сетка , как в примере 5.

Затем проинтегрируем уравнение на промежутке . Получим

.   (2.36)

Теперь достаточно указать эффективный способ вычисления интеграла в (2.36), чтобы получить вычислительное правило решения задачи (2.35). Так, если заменить интеграл в (2.35) по формуле левых прямоугольников, то получается классический метод Эйлера:

 

– задано.

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

Существует общий способ, предложенный Рунге и Кутта, построения явных одношаговых вычислительных правил [2, с. 151-156], [3, с. 21-32], [4, с. 246-250], [6, с. 450-459], [7, с. 218-230]. Рассмотрим его. Пусть

.

Тогда (2.36) примет вид

.                               (2.37)

 

Кадр 27

Чтобы перейти к промежутку интегрирования , введём в (2.37) новую переменную .

Тогда вместо (2.37) получим

.               (2.38)

Задаем три группы параметров:

;

 

.

С помощью и  определим величины :

;

;

--------------------------------------------------------

.

Все они могут быть последовательно вычислены, если параметры  известны.

Величины

 

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

(2.39)

Численный метод общего вида (2.39) для решения задачи (2.35) фактически построен. Осталось обсудить способ выбора параметров ,,.

Для этих целей введём в рассмотрение погрешность приближённого равенства (2.39):

.                                              (2.40)

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

,             (2.41)

где . Из формулы (2.41) вытекает следующее: если параметры ,, подобраны так, что справедливо

,                                  (2.42)

то погрешность приближённого равенства (2.39) будет величиной порядка не ниже , т.е.

.

Число  из равенства (2.42) называют порядком точности соответствующего численного метода.

Из проведённых выше рассуждений следуют два эквивалентных приёма получения параметров ,,.

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

2^ой приём (эквивалентен первому). Параметры ,, подбирают так, чтобы разложения

и

по степеням  совпадали до членов с наиболее высокими степенями .

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

I. Метод первого порядка точности (). При из (2.39) имеем:

(2.43)

где .

Согласно формуле (2.40) погрешность в этом случае имеет вид:

.

Для определения значения  воспользуемся первым приёмом и вычислим ,.... Получим

,

.

Т.к.   не зависит от , то не имеет смысла вычислять,.... Ясно, что можно лишь требовать. Тогда и  определяется из соотношения

.

Следовательно, . Из (2.43) получаем численный метод первого порядка точности для решения задачи (2.35):

Согласно формуле (2.41) локальная погрешность построенного правила имеет вид:

.

Итак, в случае способ Рунге-Кутты построения одношаговых вычислительных правил решения задачи Коши приводит к явному методу Эйлера.

Кадр 28

II. Метод второго порядка точности (). При  имеем приближённое равенство

.      (2.44)

Для определения параметров воспользуемся вторым  приёмом и разложим  и  в ряд по степеням . Будем иметь:

 

или

(2.45)

Аналогично

(2.46)

Сравниваем в (2.45) и (2.46) коэффициенты при и . Получим следующие соотношения:

(2.47)

При  за счёт выбора не удаётся добиться совпадения всех членов с   в разложениях (2.45), (2.46). Следовательно, локальная погрешность любого построенного правила при  есть , а сам метод будет второго порядка точности.

Из системы (2.47) находим:

(2.48)

Получается, что  можно задавать числом, не равным нулю. Возможны следующие случаи.

Случай 1. Пусть . Тогда  и

,                                           (2.49)

причём , а локальная погрешность имеет вид:

(2.50)

Вычислительное правило (2.49) называют методом Хьюна – это аналог квадратурной формулы трапеций.

Случай 2. Положим . Тогда  и

, где .

Получили вычислительное правило, аналогичное квадратурной формуле средних прямоугольников.

Случай 3. Пусть . Тогда  и

 

где . В этом случае локальная погрешность имеет представление:

.

Общая схема построения вычислительных правил третьего и более высокого порядка точности остается прежней, однако выкладки становятся громоздкими. Так, например, при выводе методов третьего порядка точности () для определения параметров метода получается 6 уравнений на 8 неизвестных, а в случае методов четвёртого порядка () – 11 уравнений на 13 неизвестных. Другими словами, при q > 2 для определения неизвестных параметров получаются системы, в которых число уравнений меньше числа неизвестных. Такие системы имеют бесконечное множество решений. Из них выбирают решения, которые приводят к наиболее простым (экономичным) алгоритмам методов Рунге-Кутты. Приведём без вывода примеры подобных методов [7, с. 224-230], [8, с. 294-304].

 

 

 

III. Методы третьего порядка точности ().

Случай 1.

,

где

 

Случай 2.

,

где

 

Случай 3.

 

где

 

IV. Методы четвёртого порядка точности ().

Случай 1.

,                   (2.51)

где

 

Это одна из самых распространённых формул Рунге-Кутты, используемых в практических расчетах.

Случай 2.

 

где

 

Представленные формулы можно считать аналогом квадратурной формулы Симпсона.

Случай 3.

,

где

 

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

Кадр 29

Увеличения точности методов можно достичь при . С примером правил Рунге-Кутты пятого порядка точности можно ознакомиться, например, в [3, с. 31], [9, с. 78-82]. В [9, с. 56-64] наглядно представлена геометрическая интерпретация методов Рунге-Кутты первого, второго и четвёртого порядков точности. Повышение порядка точности численных одношаговых методов Рунге-Кутты приводит к быстрому возрастанию трудоёмкости вычислений, т.к. на одном шаге многократно приходится вычислять значения функции при разных значениях аргументов. Заметим, однако, что методы более высокого порядка точности позволяют использовать больший шаг . Этот факт может уменьшить общие вычислительные затраты. Кроме того, методы Рунге-Кутты дают возможность использования переменного шага сетки. Так как эффективная оценка погрешности методов Рунге-Кутты затруднительна, на практике обычно на каждом шаге применяют двойной пересчёт. Суть его заключается в следующем. Исходя из точного значения , вычисляют двумя способами: с шагом  и с двойным шагом . Если полученные значения различаются в пределах допустимой точности, то шаг  для данного этапа расчёта выбран правильно. Полученное с его помощью значение решения можно принять за . В противном случае шаг уменьшается в два раза.

Проиллюстрируем сказанное на примере.

Пример 7. Используя метод Рунге-Кутты (2.51), выбрать шаг для приближённого решения задачи Коши:

 

Решение. Покажем начало процесса. Пусть . Вычислим  с шагом   и . Последовательно получаем при :

 

Отсюда

и

.

Тогда

.

Получили приближённое значение решения в точке .

Получим теперь значение .

Имеем

 

Тогда

 

Здесь .

Получим это значение теперь с шагом . Имеем

 

Следовательно,  и . Сравнивая это значение с ранее полученным , находим, что они совпадают с точностью .

Заметим, что точное решение поставленной задачи Коши есть . Его значение при   равно . Таким образом, дальнейшие расчёты на промежутке [0;0,5] можно проводить с шагом . Как доказано в [7, с. 221-224], если метод Рунге-Кутты аппроксимирует исходное уравнение, то он сходится при , причём порядок точности совпадает с порядком погрешности аппроксимации. Здесь под погрешностью аппроксимации понимается выражение

 

где .

Все методы Рунге-Кутты обобщаются на системы обыкновенных дифференциальных уравнений.

 

Кадр 30

 

Пусть дана система дифференциальных уравнений

(2.52)

и начальные условия

.                                           (2.53)

Здесь

,   , .

Выберем  и построим равномерную сетку

.

Поставим задачу определения значения приближённого решения

, по формулам

,

где  вычисляется, например, по следующей формуле (обобщение формулы (2.51) на случай систем дифференциальных уравнений):

,

(2.54)

 

Итак, зная , по формулам (2.54) вычисляем , а затем находим  . Принимая  за исходные данные и повторяя тот же процесс, находим  и т.д.

Подобным образом любая вычислительная схема метода Рунге-Кутты для одного уравнения переносится на систему уравнений вида (2.52), (2.53).

Пример 8. Методом Рунге-Кутты получить численное решение уравнения колебаний маятника в сопротивляющейся среде (см. [2, с. 122]):

,                       (2.55)

при начальных условиях

(2.56)

Здесь  - угол отклонения маятника,  - время.

Решение. Положим . Тогда уравнение (2.55) может быть представлено в виде системы:

 

Начальные условия (2.56) перепишутся так:

 

Пусть . Тогда имеем

(2.57)

Положим   и поставим задачу определения приближенных значений  в точках  с помощью формул (2.54).

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

1)   по исходным данным   сначала вычислялись  с помощью соответствующих формул из (2.54),

2)   по формулам  , , определялись приращения для  и ;

3)   за значения принимались соответственно значения ;

4)   полученные в точке  значения принимались за исходные данные и по ним вычислялись , затем новые значения  и т.д. до тех пор, пока не получились значения решения в точке .

Таблица2.3

 

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

С другими примерами численного решения задачи Коши в случае одного ОДУ, систем ОДУ первого порядка, дифференциального уравнения и систем дифференциальных уравнений высших порядков с применением методов Рунге-Кутты можно ознакомиться в [10, с.133-149].