, ,

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

Положим (метод трапеций), тогда

, (6.3.2.1)

то есть значение производной «подправляется» значением в предварительно определенной точке.

В методе прямоугольников , тогда ,

В этом случае

(6.3.2.2)

6.4. Методы Рунге-Кутта высших порядков

В методе Рунге- Кутта третьего порядка точности:

Разлагая в ряд по h до h и сравнивая с рядом Тейлора (6.1.1 ) получим следующую систему из шести уравнений для восьми неизвестных:

Наиболее употребительна в этом случае симметричная разностная схема ( аналог метода парабол при численном интегрировании ): , тогда:

, , , , , .

,

,

,

.

В методе Рунге-Кутта точности порядка получается система из 11 уравнений для 13 неизвестных.

Наиболее употребительны две вычислительные схемы:

1. Аналог метода 3/8 в численном интегрировании.

, где

,

,

,

,

1. Аналог метода парабол.

, (6.4.1)

где ,

,

,

.

Проблема выбора той или иной вычислительной схемы при заданной точности зависит от вида (х,у), так как от этого зависит величина остаточного члена.

6.5. Задание к теме и пример решения ОДУ

Найти решение задачи Коши для ОДУ:

, на интервале . K и L параметры из табл. 4.3

Решить пятью методами:

1. Метод вариации постоянных (точное решение).

2. Разложение в ряд Тейлора до четвертого порядка.

1. Метод Эйлера (6.3.1.1).

2. Метод трапеций (Коши-Эйлера) (6.3.2.1).

3. Метод Рунге-Кутта (6.4.1).

Построить графики и сравнить точность различных методов, шаг .

1. В методе вариации постоянных решение ищется в виде , Однородное уравнение имеет очевидное решение . Подстановка в неоднородное уравнение дает уравнение для коэффициента: . После интегрирования и подстановки начального условия получим: .

2. Разложение в ряд Тейлора проводится в точке х=0. Все производные в этой точке известны ;

, ;

, ;

, .

3. Метод Эйлера. Расчет ведется по формуле (6.3.1.1):

4. Метод Коши-Эйлера (метод трапеций). Вначале рассчитывается значение , которое затем используется в окончательном выражении (6.3.2.1).

1. Метод Рунге-Кутта. Последовательно вычисляются значения производной в промежуточных точках и используются в окончательном выражении с заданными весами (6.4.1).

Пример. К=3, L=2. , .

Результаты расчетов представлены в таблице:

x	0	0,5	1	1,5	2
Точное решение	2	3,34488	5,87313	10,6768	19,5562
Ряд Тейлора	2	3,3438	5,83333	10,3438	18
Метод Эйлера	2	3	4,625	7,4375	12,2812
Метод Коши- Эйлера	2	3,3125	5,72656	10,2432	18,4889
Метод Рунге- Кутта	2	3,34440	5,87111	10,6710	19,5423

1. Численное решение начально-краевых задач для дифференциальных уравнений в частных производных

Рассматриваются простейшие уравнения математической физики гиперболического, параболического и эллиптического типов с начально-краевыми условиями для дифференциальных уравнений в частных производных (ДУЧП). Задачи такого типа возникают в физике, технике, других прикладных науках.

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

1. Конечные разности.

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

где шаги сетки по координатам и соответственно, целые числа.

Неизвестная функция , участвующая в краевой задаче, заменяется искомой сеточной функцией на узлах сетки. Частные производные по координатам заменяются соответствующими конечными разностями, которые могут быть различного порядка точности по шагу сетки вдоль координаты. Пусть шаг сетки вдоль рассматриваемой координаты, значение функции в рассматриваемой точке последующие и предыдущие значения сеточной функции по данной координате. Тогда первая производная по этой координате может быть заменена правой или левой конечной разностью порядка : (7.1.2)

или центральной конечной разностью порядка : (7.1.3)

Можно записать также для первой производной конечную разность порядка :

Центральные конечные разности для второй производной порядка и выглядят следующим образом: (7.1.4)

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

1. Гиперболические уравнения

В качестве примера гиперболического уравнения рассматри-вается волновое уравнение колебаний эластичной струны , (7.2.1) с граничными условиями (закрепление струны на концах) , (7.2.2) и начальными условиями первого рода для функции (отклонение от положения равновесия) и второго рода для её производной (скорости отклонения)

, для (7.2.3)

На сетке (7.1.1) сеточная функция удовлетворяет следующим соотношениям, следующим из (7.2.1) и (7.1.4)

Обозначим . После небольших преобразований получаем явное выражение для значения сеточной функции на слое по времени через её значения на и слоях:

(7.2.4) для

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

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

На рис. 1 показан полупериод колебаний струны рассчитанный по формулам (7.2.4), (7.2.5) при следующих исходных данных, граничных и начальных условиях: , , , , , . Шаги сетки по пространственной координате и по времени ,

Число шагов по и по соответственно

рис.1.

1. Параболические уравнения

Данный тип уравнений рассмотрим на примере одномерного нестационарного уравнения теплопроводности (7.3.1) с граничными (7.3.2) и начальными условиями (7.3.3), описывающего процесс установления температуры в изолированном стержне, имеющем на концах постоянную температуру и и заданное начальное распределение температуры вдоль стержня : , (7.3.1)

, , (7.3.2)

(7.3.3)

Для аппроксимации уравнения (7.3.1) используем конечные разности (7.1.2) и (7.1.4)

Обозначим . После преобразований получаем явную четырехточечную сеточную схему, в которой значение функции на слое по времени выражается через три соседних значения на нижнем, -ом слое:

(7.3.4)

Формула (7.3.4) позволяет последовательно найти все значения сеточной функции, начиная со слоя , на котором заданы начальные условия (7.3.3). Однако вычисления по этой формуле устойчивы только в том случае, если выполняется условие . Это накладывает жесткие ограничения на шаг сетки по времени, обязывая выбирать этот шаг намного меньшим, чем шаг по пространственной координате, что существенно увеличивает время расчета и ограничивает применимость явной схемы.

Для аппроксимации уравнения (7.3.1) может быть использована левая конечная разность (7.1.2)

, что приводит к неявной четырёхточечной разностной схеме , (7.3.5) которая устойчива при любых соотношениях шагов сетки.

Из (7.3.5) следует, что для каждого слоя по времени значения неизвестной сеточной функции , связаны СЛАУ с трехдиагональной матрицей. В этой матрице на главной диагонали находится значение , а на двух соседних диагоналях - . Значение на главной диагонали близко к , т.к. значение , как правило, . Вектор в правой части (7.3.5)(при постоянном значении ) известен из вычислений на предыдущем шаге по времени и входит в правую часть СЛАУ.