50192 (Решение дифференциальных уравнений. Обзор), страница 2

набор точек в которых нужно найти решение;

само дифференциальное уравнение, записанное в некотором специальном виде, который будет описан ниже.

Один из наиболее эффективных алгоритмов интегрирования ОДУ основан на численном методе Рунге-Кутты четвертого порядка. Функция, реализующая этот метод, имеет вид rkfixed (y,x₁,x₂, npoints,D)

Здесь:

y-вектор начальных условий размерности n, где n- порядок дифференциального уравнения или число уравнений в системе (если решается система уравнений);

x₁, x₂ – граничные точки интервала, на котором ищется решение дифференциального уравнения. Начальные условия ,заданные в векторе y,- это значение решения в точке x₁;

npoints- число точек (не считая начальной точки), в которых ищется приближенное решение. При помощи этого аргумента определяется число строк (1+npoints) в матрице, возвращаемой функцией rkfixed;

D(x,y) – функция,возвращающая значение в виде вектора n элементов, содержащих первые производные неизвестных функций.

2.1 Метод Эйлера

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

Численное решение на отрезке [a, b] задачи Коши

y' = f(x, y), y(a) =y₀

состоит в построении таблицы приближенных значений

y₀, y₁, ..., y_i, ... y_N

решения y(x) в узлах сетки

a=x₀ < x₁ < ... < x_i < ...< x_N=b, y(x_i)@ y_i.

Если x_i = a+ i h, h=(b-a)/ N, то сетка называется равномерной.

Численный метод решения задачи Коши называется одношаговым, если для вычисления решения в точке x₀ + h используется информация о решении только в точкеx₀.

Простейший одношаговый метод численного решения задачи Коши - метод Эйлера. В методе Эйлера величины y_i вычисляются по формуле

y_i+1 = y_i + h f(x_i , y_i), i = 0, 1

Найдем методом Эйлера на [0, 1] с шагом h=0.2 приближенное решение задачи Коши

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

Определим правую часть уравнения

Расчетные формулы метода Эйлера для решения этой задачи имеют вид

x0=0, y0= 1, xi+1 = xi + 0.2, yi+1 = yi + 0.2(sinxi - cosyi), i =0, 1, ..., 4.

Изобразим приближенное решение графически.

y' = sin x – cos y, y(0)=1.

Определим диапазон изменения номера точки i=0,1, ..., 4

Знак присваивания можно ввести щелчком по соответствующей позиции в панели Evaluation.

Определим начальное условие - решение в начальной точке

Для того чтобы ввести символ диапазона изменения индекса , щелкните по соответствующей позиции в панели Matrix или введите с клавиатуры символ ("точка с запятой")

Определим шаг формулы Эйлера - шаг интегрирования

Для того чтобы ввести нижний индекс переменной, щелкните по соответствующей позиции в панели Matrix или в панели Calculator

Определим по формулам Эйлера значения приближенного решения в узлах сетки

Выведем в рабочий документ вычисленные значения решения

Построим график найденного решения y(x)

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

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

2.2 Метод Эйлера с шагом h/2.

Метод Эйлера допускает простую геометрическую интерпретацию. Пусть известна точка (x_i,y_i) интегральной кривой уравнения y'=f(x, y).

Касательная к интегральной кривой уравнения, проходящая через эту точку, определяется уравнением

y = y_i + f(x_i , y_i)(x-x_i).

Следовательно, вычисленная методом Эйлера точка (x_i+1 ,y_i+1 ),

Где x_i+1=x_i+h, y_i+1=y_i + h f(x_i , y_i), лежит на этой касательной.

Найдем методом Эйлера на [0, 1] с шагом h=0.2 и с шагом h=0.1 приближенное решение задачи Коши

y' = sin x – cos y, y(0)=1.

Изобразим приближенные решения графически.

Расчетные формулы метода Эйлера для решения этой задачи имеют вид

x0=0, y0= 1, xi+1 = xi + 0.2, yi+1 = yi + 0.2(sinxi - cosyi), i =0, 1, ..., 4

xi+1 = xi + 0.2, yi+1 = yi + 0.2(sinxi - cosyi), i =0, 1, ..., 9

Определим правую часть уравнения

Знак присваивания можно ввести щелчком по соответствующей позиции в панели Evaluation.

Определим диапазон изменения номера точки i=0,1, ..., 4 для вычислений с шагом h=0.2

Для того чтобы ввести символ диапазона изменения индекса , щелкните по соответствующей позиции в панели Matrix или введите с клавиатуры символ ("точка с запятой")

При решении задачи с шагом h=0.2 назовем шаг h1, аргумент - x1, а решение - y1.

Определим начальное условие

Для того чтобы ввести нижний индекс переменной, щелкните по соответствующей позиции в панели Matrix или в панели Calculator

Определим шаг формулы Эйлера - шаг интегрирования

Определим по формулам Эйлера значения приближенного решения в узлах сетки

Выведем в рабочий документ вычисленные значения решения

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

Построим график найденного решения y1(x1)

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

Определим диапазон изменения номера точки i=0,1, ..., 9 для вычислений с шагом h=0.1

Для того чтобы ввести символ диапазона изменения индекса , щелкните по соответствующей позиции в панели Matrix или введите с клавиатуры символ ("точка с запятой")

При решении задачи с шагом h=0.1 назовем шаг h2, аргумент - x2, а решение - y2.

Определим начальное условие

Для того чтобы ввести нижний индекс переменной, щелкните по соответствующей позиции в панели Matrix или в панели Calculator

Определим шаг формулы Эйлера - шаг интегрирования

Определим по формулам Эйлера значения приближенного решения в узлах сетки

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

Построим график решения y2(x2)

Построим на одном графике оба приближенные решения

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

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

Когда функции определены, щелкните по рабочему документу вне поля графиков.

2.3 Метод Рунге – Кутты

Методом Рунге-Кутты четвертого порядка точности называют одношаговый метод, относящийся к широкому классу методов Рунге-Кутты. В этом методе величины y_i+1 вычисляются по следующим формулам:

y_i+1 = y_i + h (k₁ + 2k₂ + 2k₃ + k₄)/6 , i = 0, 1, ...

k₁ = f(x_i , y_i),

k₂ = f(x_i+h/2, y_i+hk₁/2),

k₃ = f(x_i+h/2, y_i+hk₂/2),

k₄ = f(x_i+h, y_i+hk₃).

Найдем на [0, 1] приближенное решение задачи Коши y' = sin x – cos y, y(0)=1методом Рунге-Кутты 4-го порядка с шагом h=0.2 и методом Эйлера с тем же шагом.Изобразим оба приближенные решения графически

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

Определим начальное условие - решение в начальной точке

Для того чтобы ввести нижний индекс переменной, щелкните по соответствующей позиции в панели Matrix или в панели Calculator

Определим правую часть уравнения

Знак присваивания можно ввести щелчком по соответствующей позиции в панели Evaluation.

Вычислим приближенное решение на отрезке [0,1], выполнив n=1/h=5 одинаковых шагов, методом Рунге-Кутты 4-го порядка; обозначим приближенное решение Y