Ломаная эйлера для дифференциального уравнения

Численные методы решения обыкновенных дифференциальных уравнений

,

связывающее неизвестную функцию y(x), независимую переменную x и производные неизвестной функции, называется обыкновенным дифференциальным уравнением. Порядок n старшей производной называется порядком дифференциального уравнения.

В задаче Коши для дифференциального уравнения n-го порядка искомая функция y(x) кроме самого дифференциального уравнения должна удовлетворять начальным условиям

Методы решения задач для дифференциальных уравнений можно разбить на три типа: точные, приближенные и численные [9].

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

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

Численный метод решения дифференциального уравнения — это алгоритм вычисления значений искомого решения y(x) на некотором дискретном множестве значений аргумента x. При этом вычисляемые значения искомого решения y(x) являются приближенными, но могут быть и точными.

Задача Коши

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

(7.1)

. (7.2)

Требуется найти функцию , которая удовлетворяет уравнению (7.1) на интервале (x₀, X) и начальному условию (7.2) в точке x₀.

Приведем без доказательства теорему существования и единственности решения задачи Коши.

Теорема 7.1. Пусть в области функция f(x, y) непрерывна. Тогда на некотором отрезке существует решение уравнения (7.1), удовлетворяющее условию (7.2).

Если в области R функция f(x, y) удовлетворяет условию Липшица

,

то указанное решение единственно.

Произведем разбиение отрезка [x₀, X] на n частей

(7.3)

Найдем приближенные значения решения y(x) в точках x_i, i = 1, 2, …, n.

Метод Эйлера (метод ломаных)

Рассмотрим уравнение (7.1) в точках x_i, i = 0, 1, …, n – 1 и заменим производную разностной формулой

. (7.4)

Тогда получим рекуррентную формулу метода Эйлера для вычисления приближенных значений y(x_{i + 1}):

(7.5)

На рис. 7.1 дана геометрическая иллюстрация метода Эйлера (7.5). Уравнение касательной к графику решения y(x) в точке (x₀, y₀) имеет вид

(7.6)

так как . Интегральная кривая y(x) на отрезке [x₀, x₁] заменяется отрезком касательной (7.6), соединяющей точку (x₀, y₀) с точкой (x₁, y₁), где (см. рис. 7.1). Точка (x₁, y₁) уже не лежит на интегральной кривой y = y(x), удовлетворяющей начальному условию (7.2).

При i = 1 формула (7.5) дает точку (x₂, y₂), которая определяется с помощью касательной проведенной в точке (x₁, y₁) к интегральной кривой y(x), удовлетворяющей уравнению (7.1) и начальному условию .

Таким образом, с каждым шагом i метод Эйлера (7.5) дает точки (x_i, y_i), которые, вообще говоря, удаляются от интегральной кривой, соответствующей точному решению задачи Коши (7.1), (7.2). Вместо интегральной кривой метод Эйлера дает ломаную, изображенную на рис. 7.1, поэтому его часто называют методом ломаных.

Формулу (7.5) можно получить и другим способом. Рассмотрим разложение искомого решения y(x) в ряд Тейлора в окрестности точки x₀:

Ограничившись двумя слагаемыми и учитывая, что , при x = x₁ получим (7.5).

Мы также получили следующий результат — погрешность вычисления значения y₁ есть величина порядка O(h 2 ), а погрешность значения y_n — величина порядка O(h). Из-за большой погрешности метод Эйлера применяется редко.

Пример 7.1. Решить методом Эйлера на отрезке [1, 2] задачу Коши

Решение в программе Excel. Точным решением данной задачи является функция .

Выберем шаг h = 0,1.

В диапазоне A1:B12 заполняем ячейки, как показано в таблице 7.1.

В ячейке D2 записываем начальное значение y(0) = 1, а в ячейке C2 — формулу =(1+C2^2)*0,5/A2, которая соответствует выражению f(x₀, y₀). Затем выделим ячейку C2 и маркером заполнения протянем вниз до C11.

В ячейку D3 введем формулу =D2+0,1*C2 и протянем ячейку D3 маркером заполнения вниз до D12.

В столбце E2:E12 вычислим значения точного решения. Для этого в ячейку E2 введем формулу =TAN(LN(КОРЕНЬ(A2))) и протянем ячейку E2 маркером заполнения вниз до E12.

В столбце F2:F12 вычислим абсолютные ошибки с помощью формулы =ABS(C2-D2), введенной в ячейку F2 и скопированной в остальные ячейки.

Как видим, значения, полученные методом Эйлера, приближенно равны точным значениям. Абсолютная погрешность каждого значения равна в среднем ε ≈ 0,01. При этом видно, что погрешность y_i+1 увеличивается вместе с номером i.

y(x)Ошибки200,0476910,002309420,0914140,004154530,131940,005676640,1698420,006964750,2055560,00808860,2394260,009065970,271720,009951080,3026580,0107581190,3324180,0115031210

0,361150,0122

Создадим в программе Excel пользовательские функции для решения задачи Коши (7.1), (7.2) методом Эйлера, используя встроенный язык Visual Basic. Правую часть уравнения (7.1) для нашего примера определяет функция fxy(x, y):

Function fxy(ByVal x, ByVal y): fxy = (1 + y ^ 2) / (2 * x): End Function

Метод Эйлера реализует функция Eiler(x0, y0, x, n):

Function Eiler(ByVal x0, ByVal y0, ByVal x, ByVal n)

y1 = y0 + h * fxy(x0, y0): x0 = x0 + h: y0 = y1

Next i: Eiler = y1: End Function

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

Здесь исходными данными являются начальные значения x0, y0, значение переменной x и числовой параметр n, который означает число частей, на которые разбивается отрезок [x0, x], т.е. метод Эйлера применяется с шагом h = (x — x0) / n. Результатом является приближенное значение решения в точке x.

Чтобы ввести эти функции в Excel, выполним команду «Сервис – макрос – Редактор Visual Basic» и в открывшемся окне выбираем пункт меню «Insert – Module» и вводим тексты программ-функций (см. рис. 7.2). Теперь эти функции можно использовать в формулах, как и обычные функции программы Excel.

Переходим в Excel, и построим таблицу решения примера 7.1 с помощью созданных функций. В таблице 7.2 приведены результаты вычислений, а в таблице 7.3 показаны формулы. Уменьшение погрешности по сравнению с таблицей 7.1 объясняется тем, что метод Эйлера в таблице 7.2 применяется с шагом h = 0,01.

Отметим, что правая часть уравнения вызывается в описании функции Eiler. Если понадобится решить задачу Коши с другой правой частью, то достаточно переопределить функцию fxy(x, y).

Рис.7.2. Окно редактора Visual Basic

Создадим в программе Mathcad функции для решения задачи Коши методом Эйлера и вычислим решение в точке x = 2 с шагом h = 0,1 (n = 10) и шагом h = 0,01 (n = 100):

Как видим, результаты расчетов в программе Mathcad совпадают со значениями, полученными в Excel.

Метод Эйлера с уточнением

Для повышения точности метода Эйлера применяют следующий прием. Сначала находят приближенное значение решения по методу Эйлера

,

а затем уточняют его по формуле

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

(7.7)

Метод Эйлера-Коши имеет погрешность порядка O(h 2 ). Геометрическая иллюстрация метода Эйлера-Коши дается на рис. 7.2. Очередное приближение метода Эйлера-Коши соответствует точке пересечения диагоналей параллелограмма, построенного на двух звеньях ломаной метода Эйлера.

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

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

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

Выведем формулы метода Рунге-Кутта второго порядка точности. Для этого решение представим куском ряда Тейлора, отбрасывая члены с порядком выше второго. Тогда приближенное значение искомой функции в точке x₁ можно записать в виде:

. (7.8)

Вторую производную можно выразить через производную функции f(x, y), однако в методе Рунге-Кутта вместо производной используют разность

,

соответственно подбирая значения параметров . Тогда (7.8) можно переписать в виде

, (7.9)

где α, β, γ и δ — некоторые параметры. Рассматривая правую часть (7.9) как функцию аргумента h, разложим её по степеням h:

,

и выберем параметры α, β, γ и δ так, чтобы это разложение было близко к (7.8). Отсюда следует, что

С помощью этих трех уравнений выразим β, γ и δ через параметр α, тогда получим

(7.10)

Теперь, если вместо (x₀, y₀) в (7.10) подставить (x₁, y₁), то получим формулу для вычисления y₂ — приближенного значения искомой функции в точке x₂.

В общем случае метод Рунге-Кутта применяется на произвольном разбиении отрезка [x₀, X] на n частей, т.е. с переменным шагом

(7.11)

Параметр α выбирают равным 1 или 0,5. Запишем окончательно расчетные формулы метода Рунге-Кутта второго порядка с переменным шагом для α = 1:

(7.12)

(7.13)

Дадим геометрическую интерпретацию методов Рунге-Кутта (7.12), (7.13).

На рисунке 7.3 а) иллюстрируется формула (7.12). Сначала по методу Эйлера вычисляется приближенное значение в точке x_i + h_i/2. В этой точке определяется касательная к интегральной кривой, параллельно которой через точку (x_i, y_i) проводится прямая до точки пересечения с прямой x = x_{i + 1}. Ордината точки пересечения y_{i + 1} принимается за приближенное значение искомого решения в точке x_{i + 1}.

Рисунок 7.3 б) интерпретирует формулу (7.13). Методом Эйлера вычисляется приближенное значение в точке x_i+1. В этой точке тангенс угла наклона касательной к интегральной кривой равен выражению . Через точку (x_i, y_i) проводится прямая, тангенс угла наклона которой определяется как полусумма угловых коэффициентов касательных, проведенных через точки (x_i, y_i) и . За приближенное значение искомого решения в точке x_i+1 принимается ордината точки пересечения этой прямой с прямой x = x_i+1.

Отметим, что метод (7.13) совпадает с методом Эйлера-Коши (7.7).

Теорема 7.1. Если правая часть f(x, y) уравнения (7.1) непрерывна и ограничена вместе со своими производными до второго порядка включительно, то решение, полученное по формулам (7.12), (7.13), равномерно сходится к решению задачи (7.1), (7.2) с погрешностью .

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

(7.14)

Пример 7.2. Решить методом Рунге-Кутта (7.12) задачу Коши из примера 7.1:

Решение в программе Excel. В таблице 7.4 приведены результаты применения метода Рунге-Кутта второго порядка (7.12) с шагом h = 0,1.

Сравнивая таблицы 7.2 и 7.4 можно увидеть, что метод Рунге-Кутта второго порядка (7.12) с шагом h = 0,1 дает лучшие результаты, чем метод Эйлера с шагом 0,01.

Ломаная Эйлера и e-приближенное решение

2.5. Ломаная Эйлера и e-приближенное решение

Рассмотрим систему уравнений

(2)

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

Совокупность n функций z1(t), . zn(t) называется e-приближенным решением системы (2) на отрезке А, если каждая из этих функций непрерывна, имеет кусочно-непрерывную производную и

во всех точках tÎK, кроме точек разрыва непрерывности этой производной.

Пусть задана начальная точка (t0, x10, …, хn0) и пусть функции fi(t, xi. хn) непрерывны по t в области G и удовлетворяют в этой области условию Липшица по переменным t, x1, х2, . хn. Можно показать, что в этом случае функции fi(t, x1. хn) будут непрерывны по совокупности переменных t, x1. хn в области G. Из непрерывности функций fi (t, x1. хn) в замкнутой области G сле­дует их равномерная непрерывность. Таким образом, для любого e>0 найдется такое d>0, зависящее только от e, что при

будет справедливо неравенство

Построим e-приближенное решение системы (2). Для этого разобьем область G на кубы со сторонами, меньшими d (для случая n=1 построение проведено на рис. 2, в этом случае область разбивается на квадраты). Из точки (t, xlo, . хn0) проведем прямую

Эту прямую продолжим до пересечения с одной из сторон соответствующего куба. Обозначим точку пересечения (t1, x11. xn1). Из этой точки проведем прямую

которую продолжим до пересечения с одной из сторон куба; обозначим точку пересечения (t2, x12. xn2), через эту точку проводим новую прямую

В результате указанных действий получим ломаную xi=xi(t) (i=l, 2, . n), называемую ломаной Эйлера. Эта ломаная представляет собой непрерывную кусочно-линейную функцию. Ломаную Эйлера мы можем продолжить до границы области G.

Пусть xi(t) (i=l, 2, . n) — точное решение системы (2), удовлетво­ряющее начальным условиям. Обозначим через si(t) (i=1, 2, . n) e-приближенное решение системы (1) для тех же начальных условий. Тогда

Отсюда следует, что если |t–t0| 0 существует такое d(e, h)>0, что другое решение x=s(t, t0, z0), удовлетворяющее начальным условиям

Ломаная эйлера для дифференциального уравнения

Variant 19 (Sukach Maxim, BS17-03)

Найдем

В итоге, наше решение принимает вид:

Метод Эйлера дает возможность приближенно выразить функцию теоретически с любой наперед заданной точностью. Суть метода Эйлера в пошаговом вычислении значений решения y=y(x) дифференциального уравнения вида y’=f(x,y) с начальным условием (x0;y0). Метод Эйлера является методом 1-го порядка точности и называется методом ломаных.

Для вычисления используются следующие формулы:

Метод Эйлера и точное решение при x0 = 0, xf = 9, y0 = 1, h = 0.1

Метод Эйлера и точное решение при x0 = 0, xf = 3, y0 = 1, h = 0.1

Метод Эйлера и точное решение при x0 = 0, xf = 1, y0 = 1, h = 0.1

Усовершенствованный метод Эйлера

Суть усовершенствованного метода Эйлера в пошаговом вычислении значений решения y=y(x) дифференциального уравнения вида y’=f(x,y) с начальным условием (x0;y0). Усовершенствованный метод Эйлера является методом 2-го порядка точности и называется модифицированным методом Эйлера.

Разница между данным методом и методом Эйлера минимальна и заключается в использовании следующих формул:

Усовершенствованный Метод Эйлера и точное решение при
x0 = 0, xf = 9, y0 = 1, h = 0.1

Усовершенствованный Метод Эйлера и точное решение при
x0 = 0, xf = 3, y0 = 1, h = 0.1

Усовершенствованный Метод Эйлера и точное решение при
x0 = 0, xf = 1, y0 = 1, h = 0.1

Классический метод Рунге-Кутты

Суть метода Рунге-Кутты в пошаговом вычислении значений решения y=y(x) дифференциального уравнения вида y’=f(x,y) с начальным условием (x0;y0). Классический метод Рунге-Кутты является методом 4-го порядка точности и называется методом Рунге-Кутты 4-го порядка точности.

Ну и как обычно, формулы:

Классический метод Рунге-Кутты и точное решение при x0 = 0, xf = 9, y0 = 1, h = 0.1

Классический метод Рунге-Кутты и точное решение при x0 = 0, xf = 3, y0 = 1, h = 0.1

Классический метод Рунге-Кутты и точное решение при x0 = 0, xf = 1, y0 = 1, h = 0.1

Сравнение методов для заданной задачи

Размер ошибки всех методов на промежутке [0, 9] с шагом 0.1

Размер ошибки всех методов на промежутке [0, 3] с шагом 0.1

Размер ошибки всех методов на промежутке [0, 1] с шагом 0.1

Очевидно что, классический метод Рунге-Кутты справляется с задачей аппроксимации в случае данного уравнения намного лучше чем Метод Эйлера и Усовершенствованный метод Эйлера.

График глобальной средней ошибки

Глобальная ошибка в зависимости от размера шага H на промежутке от 0.01 до 0.91 для x0 = 1, xf = 9

🔍 Видео