Найти решение дифференциального уравнения dy/dx=F(x,y). Составить таблицу значений аргумента X и функции Y=Φ(Х) при условии, что X изменяется на интервале X0£X£Xmax с шагом DX. Построить график решения дифференциального уравнения.
Как известно, решением дифференциального уравнения dy/dx=F(x,y) при заданных начальных условиях X0 и Y0 (т.е. решением задачи Коши) является функция, проходящая через точку с координатами X0, Y0. Численные методы решения дифференциальных уравнений (например методы Эйлера, Эйлера-Коши, Рунге-Кутта) позволяют получить таблицу значений этой функции на заданном интервале изменения аргумента Х (Х0<=X<=Xmax) с неким шагом h.
Рассмотрим решение дифференциального уравнения методом Рунге-Кутта, дающим наиболее высокую точность (из трех названных). Алгоритм решения дифференциального уравнения представляет собой итерационный вычислительный процесс, т.е. значения, полученные на одном шаге, являются начальными для следующего.
Формулы численного метода Рунге-Кутта:
yk+1=yk+(k1+2k2+2k3+k4)/6
k1= h*F(xk,yk)
k2= h*F(xk+h/2,yk+k1/2)
k3= h*F(xk+ h/2,yk+k2/2)
k4= h*F(xk+h,yk+k3)
Применим метод Рунге-Кутта для решения дифференциального уравнения dy/dx=x+y/x при начальных значениях x0=1, y0=0.
Решение будем искать на интервале [1;3] c шагом 0,1. Значение шага h занесем в ячейку A5. Используя в формулах значение шага, будем делать абсолютную ссылку на эту ячейку. Это важно для получения верного решения.
Шапку таблицы заполним следующим образом:
| A | B | C | D | E | |
| h | x | x+h | x+h/2 | y | |
| 0,1 | 1,1 | 1,05 |
| F | G | H | I | |
| k1 | k2 | k3 | k4 | |
| h*f(xk,yk) | h*f(xk+h/2,yk+k1/2) | h*f(xk+h/2,yk+k2/2) | h*f(xk+h,yk+k3) |
Верхний ряд этих таблиц – названия столбцов, левый столбец – номера строк рабочего листа.
Часть результата приведена на рисунке ниже.
Таблица решения дифференциального уравнения
В ячейку С5 следует ввести формулу для вычисления x+h: =B5+$A$5. В ячейку D5 введите формулу расчета x+h/2: B5+$A$5/2. В ячейку E6 поместите формулу для вычисления следующего значения решения дифференциального уравнения: =E5+(F5+2*G5+2*H5+I5)/6 (см. формулы метода Рунге-Кутта). В ячейки F5 – I5 поместите формулы для расчета К1, К2, К3, К4. Для их записи пользуйтесь формулами метода, правой частью дифференциального уравнения (x+y/x) и не забывайте, что формула будет содержать не имена переменных, а ссылки на ячейки с нужной информацией. Затем с помощью меню Правка–Заполнить–Прогрессия…( или с помощью контекстного меню) заполните столбец значений х.
Решающий этап – заполнение таблицы. Предварительно выделите ячейки С5 и D5 и скопируйте эти формулы на соседнюю строчку с помощью маркера Автозаполнения. Затем выделите ячейки с F5 по I5 и проделайте то же самое.
Теперь можно быстро заполнить всю таблицу. Для этого выделите ячейки с С6 по I6, поместите указатель мыши на маркер Автозаполнения и протащите мышь, удерживая ее левую кнопку до строки, в которой кончается ряд значений x. При этом произойдет копирование всех формул и автоматический расчет по ним. Столбцы x и y дают решение дифференциального уравнения (таблицу значений искомой функции).
Для построения графика выделите несмежные диапазоны значений x и y=Φ(x), расположенные в столбцах B и E. (Тип диаграммы - Точечная диаграмма со значениями, соединенными сглаживающими линиями). Затем добавьте график точного решения уравнения по формуле в таблице.
Задания для самостоятельного решения
Найдите решение дифференциального уравнения первого порядка dy/dx=F(x,y) методом Рунге-Кутта на интервале [xн,xк] с шагом h=0.1 при начальном условии y0=y(xк).
Постройте численное и точное решение на одном графике. Сделайте выводы о совпадении (или не совпадении) решений уравнения. (Вариант – в журнале преподавателя).
| Исходные данные для расчета | Данные для проверки | ||||
| № | F(x,y) | Начальное условие | [xн;xк] | Точное решение | Значение y в 2-х точках |
| y(2)=1 | [2;4] | 
| y(3)=1.0606 y(4)=2.8284 | |
| y(1)=1 | [1;1,9] | 
| y(1.5)=1.9365 y(1.9)=5.8643 | |
 
| y(1)=1 | [1;3] | 
| y(2)=1.4577 y(3)=2.1343 | |
| y(1)=-1 | [1;3] | 
| y(2)= -0.5 y(3)= | |
| y(0)=1 | [0;2] | 
| y(1)=1.1451 y(2)=3.1822 | |
| y(0)=1 | [0;2] | 
| y(1)=1.7321 y(2)=2.2361 | |
| y(0)=0 | [0;2] | 
| y(1)=1.8508 y(2)= -4.806 | |
| y(0)=1 | [0;0.9] | 
| y(0.5)=1.366 y(0.9)=1.3359 | |
| y(1)=1 | [1;3] | 
| y(2)=0.25 y(3)=0.111 | |
| y(1)=2 | [1;3] | 
| y(2)=6.4261 y(3)=13.621 | |
| y(0)=-1 | [0;2] | y=sinx-1 | y(1)=-0.1585 y(2)=-0.0907 | |
| y(1)=1 | [1;3] | 
| y(2)=1.8402 y(3)=2.5091 | |
| y(0)=1 | [0;2] | 
| y(1)=0.8921 y(0)=0.1856 | |
| y(1)=1 | [1;1.7] | 
| y(1.4)=2.69 y(1.7)=30.91 | |
| y(0)=1 | [0;0.9] | 
| y(0.5)=3
y(0.9)=19
| |
| y(1)=1 | [1;3] | 
| y(2)=0.8 y(3)=0.7143 | |
| y(1)=-1 | [1;3] | y=-x | y(2)=-2 y(3)=-3 | |
| y(1)=1 | [1;3] | 
| y(2)=0.6581 y(3)=0.8406 | |
| y(2)=2 | [2;4] | 
| y(3)=2.4495 y(4)=2.8284 | |
| y(1)=0 | [1;3] | 
| y(2)=1.1287 y(3)=2.3239 | |
| y(  )=1
| [ ;  ]
| y=2sinx-1 | y  =-1
y( )=-3
| |
| y(0)=1 | [0;1] | 
| y(0.5)=1.1802 y(1)=3.5305 | |
| 4x-2y | y(0)=0 | [0;2] | y=2x-1-e-2x | y(1)=1.1353 y(2)=3.0183 | |
| y(0)=1 | [0;2] | 
| y(1)=0.5518 y(2)=0.0549 | |
| y(1)=1 | [1;2] | 
| y(1.5)=0.866 y(2)=0 | |
| y(0)=1 | [0;2] | y=(x+1)(x2+1) | y(1)=4 y(2)=15 | |
| cos x – y | y(0)=1 | [0;2] | 
| y(1)=0.8748 y(2)=0.3142 | |
| 0.5xy | y(0)=1 | [0;2] | 
| y(1)=1.284 y(2)=2.718 | |
| ex+y | y(0)=1 | [0;1.5] | 
| y(1)=5.4365 y(1.5)=11.2 | |
| e2x - y | y(0)=1 | [0;1.5] | 
| y(1)=2.7082 y(2)=6.844 | |
| y(1)=1 | [1;3] | y=2x-1 | y(2)=3 y(3)=5 | |
| y(0)=0 | [0;2] | y=ex-1 | y(1)=1.7182 y(2)=6.389 | |
| y(0)=0 | [0;2] | y=sinx +cosx | y(1)=1.3817 y(2)=0.4931 | |
| y(1)=1 | [1;3] | 
| y(2)=0.5 y(3)=1/3 | |
| ytgx-y2cosx | y(0)=1 | [0;1.5] | 
| y(1)=0.9254 y(1.5)=5.655 | |
| y(1)=1 | [1;3] | 
| y(2)=0.5906 y(3)=0.4765 | |
| y(1)=1 | [1;1.9] | 
| y(1.5)=3 y(1.9)=19 | |
| x-y | y(0)=0 | [0;2] | y=x-1+e-x | y(1)=0.368 y(2)=1.1353 | |
| y(1)=1 | [1;3] | y=x | y(2)=2 y(3)=3 | |
| y(1)=0 | [1;3] | 
| y(2)=-3.5 y(3)=-8.66 | |
| y(1)=1 | [1;3] | 
| y(2)=1.4832 y(3)=1.6124 | |
| y(1)=0 | [1;3] | 
| y(2)=-1.3863 y(3)=-3.295 | |
| y(0)=0 | [0;1] | 
| y(0.5)=1.459 y(1)=6.8731 | |
| 2y-x2 | y(0) 
| [0;2] | 
| y(1)=1.25 y(2)=3.25 | |
| 2x-y | y(0)=-1 | [0;2] | y=2x-2+e-x | y(1)=0.3678 | |
| y(0)=0 | [0;2] | 
| y(1)=0.3651 y(2)=0.3408 | |
| y(1)=2 | [1;1.9] | y= 
| y(1.5)=4 y(1.9)=20 | |
- 
| y(1)=0 | [1;3] | 
| y(2)=-0.75 y(3)=-1/3 | |
| y(1)=2 | [1;3] | 
| y(2)=2.33 y(3)=2.5 | |
| y(1)=1 | [1;3] | 
| y(2)=3.0895 y(3)=5.3642 |