6.3. Метод разложения в степенной ряд.

Рис. 6.4. Алгоритм метода Эйлера

К методу Эйлера можно прийти и иными путями, не исполь­зуя геометрических построений. Простейший из них таков. Ап­проксимируем производную у', входящую в уравнение (6.1), пра­восторонним конечно-разностным отношением в точке (см. формулу (5.21)), и получим . Здесь - точное значение искомого решения в точке . Обозначая через прибли­женное значение решения, получим: — форму­ла (6.7).

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

(6.9)

Число известно из начального условия (6.2), значение про­изводной также известно, поскольку из дифференциально­го уравнения (6.1) следует, что . Ограничиваясь в разложении (6.9) первыми двумя членами ряда и положив немедленно получаем формулу метода Эйлера (6.7).

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

(6.10)

(использовано правило дифференцирования сложной функции само уравнение (6.1)) Положив в (6.9) и используя (6.10), получим во втором порядке разложения:

(6.11)

Первые два слагаемых в первой части (6.11) – формула метода Эйлера, а следующее слагаемое – поправка к ней.

В связи с таким подходом становится понятно, почему метод Эйлера называют методом первого порядка аппроксимации; метод связанный с формулой (6.11), - второго порядка. Продолжая разложение в ряд, можно, в принципе, достичь более высокого порядка аппроксимации решения.

Кроме того, становится очевидным, что оценка погрешности применения метода Эйлера на одном шаге такова ( - точное значение решения в точке ):

(6.12)

где К – некоторая константа, вид которой следует из сравнения (6.11) и (6.12). Тем не менее, следует помнить, что оценка погрешности любого многошагового метода на одном шаге не дает возможности учесть накопление погрешностей при многих шагах, на каждом из которых сказывается погрешность предыдущего шага. Формула, аналогичная формуле (6.11), но для любого узла, имеет вид:

(6.13)

Пример 6.2. Получить рабочую формулу для интегрирования дифференциального уравнения, рассмотренного в примере 6.1, во втором порядке аппроксимации.

Имеем: Таким образом, с помощью (6.13) получаем:

Метод разложения в степенной ряд на практике почти не применяется, однако он создает теоретическую базу для построения более удобных, практически используемых методов. Рассмотрим ниже один из них.

6.4 Метод Рунге-Кутта.

Отмеченное выше неудобство решения дифференциальных уравнений путем разложения в степенной ряд связано, прежде всего, с тем, что в расчетные формулы типа (6.13) входят наряду с функцией ее частные производные. В формулу второго порядка (6.13) входят производные первого порядка, в формулу третьего порядка – производные второго порядка и т.д. Поскольку нахождение формул для производных – задача не для стандартно применяемых при решении математических задач систем программирования, то производные приходится искать вручную и конструировать с их помощью достаточно громоздкие формулы.

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

Проиллюстрируем это на примере одного из возможных методов второго порядка. Для метода Эйлера путем простейшей правосторонней аппроксимации производной имеем: .

В это выражение входят значения функции в двух точках. Если взять в правой части уравнения (6.1) значение в точке с номером , то придем к методу Эйлера. Однако можно рассудить и иначе: раз для аппроксимации производной взяты две точки, то лучшей аппроксимацией правой части уравнения будет полусумма . Тогда для нахождения получим равенство . Поскольку оно представляет собой уравнение для , то решать его можно методом итераций, причем в качестве первого приближения взять то значение , которое определяется методом Эйлера. В итоге получим:

(6.14)

где

(6.15)

- расчетные формулы нового метода, обеспечивающего пошаговое интегрирование дифференциального уравнения (6.1). Описанная процедура получения формул (6.14) и (6.15) опирается лишь на некие правдоподобные утверждения. Докажем теперь, что (6.14) согласуется с (6.13) с точностью до второго порядка по включительно. Для этого разложим величину , входящую в (6.14), в ряд вплоть до членов первого порядка по

(6.16)

(учитывая, что входит в оба аргумента). Подставим (6.16) в (6.14)

в точности совпадает с (6.13). Таким образом, несмотря на различный вид, формулы (6.14) и (6.13) согласуются друг с другом во втором порядке по включительно. Это, а также то обстоятельство, что в формулы (6.14) и (6.15) входят значения лишь самой функции , дает основания называть (6.14) формулой Рунге-Кутта 2-го порядка.

На самом деле, (6.14) – лишь одна из множества возможных формул Рунге-Кутта 2-го порядка. Попробуем отыскать такие формулы в семействе формул вида:

(6.17)

где некоторые коэффициенты, значения которых следует подобрать так, чтобы (6.17) согласовывалось с (6.13) с точностью до второго порядка по .

Разложим в ряд аналогично тому, как это было сделано при получении (6.16). Получим:

(6.18)

Сопоставляя (6.18) с (6.13), получаем ограничения на коэффициенты , при выполнении которых будет иметь место совпадение (6.17) с (6.13):

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

В частности формула (6.14) получается при .

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

Чем выше порядок формул Рунге-Кутта, тем более точные значения они дают. На практике соблюдается некоторый компромисс между высоким порядком формул и их громоздкостью, с одной стороны, и объемом вычислений по ним для достижения заданной точности, с другой. Ниже приведена одна из самых популярных формул четвертого порядка, часто называемая просто, без уточнений, формулой Рунге-Кутта:

(6.19)

где

(6.20)

причем сначала последовательно вычисляются а лишь затем по формуле (6.19).

Общий недостаток методов Рунге-Кутта - отсутствие простых способов оценки погрешности метода. Погрешность на одном шаге оценить сравнительно нетрудно, гораздо труднее оценить накопление погрешностей на протяжении многих шагов. Широко используемый на практике для этих методов полуэмпирический способ контроля точности – двойной счет. Допустим, что - точное решение уравнения при . Тогда для рассмотренного выше метода четвертого порядка где верхний индекс, означает шаг, с которым вычислено приближенное значение. Точно так же при решении с шагом : Отсюда получаем: Следовательно, ошибка при вычислении с шагом есть:

(6.21)

Из формулы (6.21), в частности, следует, что при достаточно малом и малых погрешностях вычислений решение уравнения (6.1), полученное методом Рунге-Кутта по формулам (6.19)-(6.20), будет близким к точному.

Пример 6.3. Решить дифференциальное уравнение на отрезке [0; 0,5] с начальным условием и шагом методом Эйлера и Рунге-Кутта 4-го порядка.

Начнем с того, что получим аналитическое решение уравнения . При заданных начальных условиях это решение выглядит так:

Пользуясь этой формулой, с помощью табулирования можно получить таблицу "точного" решения заданного уравнения.

Порядок вычислений по методу Эйлера был продемонстрирован в примере 6.1. Покажем в деталях один шаг вычислений по методу Рунге-Кутта 4-го порядка. Проведем последовательные вычисления по формулам (6.20):

Далее действуем по формуле (6.19):

Поскольку вычисления достаточно громоздки и трудоемки, то лучше выполнять их используя возможности компьютера. Для сравнения результаты численного решения уравнения указанными методами и значения "точного" решения, вычисленные в соответствующих точках (с шестью знаками после запятой), сведены в табл. 6.2. Из табл. 6.2. видно, что результаты, полученные по меиоду Рунге-Кутта, практически совпадают с "точным" решением – в отличие от соответствующих значений, полученных методом Эйлера.Подтвердить истинную точность результатов, полученных методом Рунге-Кутта, можно также путем повторения решения с мелким шагом.

Таблица 6.2.

x	y
	Метод Эйлера	Метод Рунге-Кутта	"точное" решение
0,00	1,000000	1,000000	1,000000
0,05	1,050000	1,049958	1,049958
0,10	1,099875	1,099659	1,099659
0,15	1,149369	1,148837	1,148837
0,20	1,198218	1,197217	1,197217
0,25	1,246146	1,244520	1,244520
0,30	1,292877	1,290462	1,290462
0,35	1,338127	1,334758	1,334758
0,40	1,381617	1,377128	1,377128
0,45	1,423065	1,417295	1,417295
0,50	1,462199	1,454991	1,454991