Метод рунге кутта решения диф уравнений и их систем

Содержание

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

Существуют и другие явные одношаговые методы. Так, рассмотренные метод Эйлера (1.15) и его модифицированные варианты (1.22), (1.23) и (1.25), (1.26) являются частными случаями методов первого и второго порядков, относящихся к классу методов Рунге – Кутта. Эти методы используют для вычисления значения значение yi, а также значения функции f(x, у) при некоторых специальным образом выбираемых значениях и у. На их основе могут быть построены разностные схемы разного порядка точности.

Широко распространен метод Рунге – Кутта четвертого порядка.

Запишем алгоритм этого метода в виде

Таким образом, данный метод Рунге – Кутта требует на каждом шаге четырехкратного вычисления правой части f(x, Y) уравнения (1.9). Суммарная погрешность этого метода есть величина О(h4).

Метод Рунге – Кутта (1.27) требует большего объема вычислений по сравнению с методом Эйлера и его модификациями, однако это окупается повышенной точностью, что дает возможность проводить счет с большим шагом. Другими словами, для получения результатов с одинаковой точностью в методе Эйлера потребуется значительно меньший шаг, чем в методе Рунге – Кутта (1.27).

Проведем сравнительную оценку рассмотренных методов Рунге – Кутта на простом примере, позволяющем получить также и точное решение.

Пример. Решить задачу Коши

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

Теперь решим данную задачу численно с помощью рассмотренных выше методов. Результаты вычислений приведены в табл. 1.1. Анализ решения позволяет проследить рост погрешности с возрастанием xi. Как видно из табл. 1.1, самым точным является решение, полученное методом Рунге – Кутта четвертого порядка. При хi = 1погрешность составляет менее 0,003%. Для модифицированных методов Эйлера погрешность при хi = 1 составляет около 1%, а для самого метода Эйлера — почти 18%. Следовательно, при больших х метод Эйлера может привести к еще более существенным погрешностям, и в таких случаях предпочтительнее пользоваться численными методами высших порядков точности.

Результаты вычислений xi разными методами

С уменьшением шага hлокальная погрешность метода Эйлера снизится, однако при этом возрастет число узлов, что неблагоприятно повлияет на точность результатов. Поэтому метод Эйлера применяется сравнительно редко при небольшом числе расчетных точек. Наиболее употребительным одношаговым методом является метод Рунге – Кутта.

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

Начальные условия зададим в виде

По аналогии с (1.27) запишем формулы Рунге – Кутта для системы двух уравнений:

К решению систем уравнений сводятся также задачи Коши для уравнения высших порядков. Например, рассмотрим задачу Коши для уравнения второго порядка

Введем вторую неизвестную функцию . Тогда сформулированная задача Коши заменяется следующей:

В заключение еще раз отметим особенность одношаговых методов, состоящую в том, что для получения решения в каждом новом расчетном узле достаточно иметь значение сеточной функции лишь в предыдущем узле. Это позволяет непосредственно начать счет при i = 0 по известным начальным значениям. Кроме того, указанная особенность допускает изменение шага в любой точке в процессе счета, что позволяет строить численные алгоритмы с автоматическим выбором шага.

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

Обыкновенное дифференциальное уравнение первого порядка. Задача Коши, суть метода Рунге-Кутта. Выбор среды разработки. Программная реализация метода Рунге-Кутта 4-го порядка. Определение порядка точности метода. Применение языка программирования C++.

Рубрика	Программирование, компьютеры и кибернетика
Вид	курсовая работа
Язык	русский
Дата добавления	16.05.2020
Размер файла	163,4 K

Отправить свою хорошую работу в базу знаний просто. Используйте форму, расположенную ниже

Студенты, аспиранты, молодые ученые, использующие базу знаний в своей учебе и работе, будут вам очень благодарны.

Размещено на http://www.allbest.ru/

ФЕДЕРАЛЬНОЕ Государственное АВТОНОМНОЕ образовательное учреждение Высшего профессионального образования

«БЕЛГОРОДСКИЙ ГОСУДАРСТВЕННЫЙ НАЦИОНАЛЬНЫЙ

Институт инженерных технологий и естественных наук

Кафедра математического и программного обеспечения информационных систем

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

по дисциплине «Методы вычислений»

студента очного формы обучения

направления подготовки 010500.62

«Математическое обеспечение и администрирование

3 курса группы 07001302

Данькова Николая Алексеевича

1. Теоретическая часть

1.1 Обыкновенное дифференциальное уравнение первого порядка. Задача Коши

1.2 Суть метода Рунге-Кутта

1.3 Выбор среды разработки

2. Практическая часть

2.1 Программная реализация метода Рунге-Кутта 4-го порядка

Список использованных источников

Введение

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

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

Цель курсовой работы: изучить метод Рунге — Кутта 4-го порядка для решения обыкновенных дифференциальных уравнений.

Постановка задачи: необходимо составить программу, позволяющую решать обыкновенные дифференциальные уравнения методом Рунге — Кутта 4-го порядка.

Курсовая работа состоит из 3 разделов, содержит 6 рисунков, 3 листинга, 1 приложение и 18 страниц.

1. Теоретическая часть

= f(x, y) (1)

Согласно теореме существования и единственности для любой точки (x,y) ?G найдется решение у = у(х), определенное на некотором интервале (х -д, х +д), удовлетворяющее условию y(x) = y_0, такое, что точки (x,y(x)) ?G и y`_x ? f(x, y(x)), причем это решение будет единственным. Задача для уравнения (1) с начальным условием у(х) = y (задача Коши) состоит в нахождении функции у(х), обращающей и уравнение (1), и начальное условие в тождество. Допустим, что значения, которые принимает независимое переменное х, принадлежат интервалу (Х, X_N ) и запишем задачу Коши:

Разобьём отрезок [Х, X_N ] на N частей так, что x_n₊₁ — х_n = h_n _,

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

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

y(x0 ) = y (0) =(y (0)1 , . y (0)n ).

Методе Рунге-Кутта четвертого порядка вычисляет значение y_i k по формуле:

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

Описание алгоритма метода Рунге-Кутты решения систем дифференциальных уравнений первого порядка:

Задаем начальное x и конечное x_n значения отрезка интегрирования, массив начальных значений y, т.е. значения системы y0(i)… y0(n) в точке x, а также количество шагов интегрирования m.

Вычисляем шаг интегрирования h.

На основе известных x и y0(1)… y0(n)вычисляем правые части системы для всех n уравнений.

Вычисляем коэффициенты k₁ и для всех n уравнений.

Вычисляются коэффициенты k₂, k_3, k₄ для всех n уравнений системы в соответствии с формулами метода.

Вычисляем y₁… y_т для первого шага интегрирования в соответствии с главной формулой метода.

Итерации повторяются для вычисления состояний системы на каждом шаге интегрирования.

Конец цикла по х.

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

Определить Массив y(m)

Цикл по i от 1 до m с шагом 1

Определить Массив y0(m)

Определить Локальный Массив k1(m),k2(m),k3(m),k4(m),y1(m),y0_(m),f(m)

Локальные переменные i,x0,h

h=(xe-xn)/(n-1) * шаг интегрирования

Цикл по x0 от xn до xe с шагом h

Цикл по i от 1 до m

Определить Массив y(n),f(n)

Определить Массив y(m)

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

set decimals to 10

LOCAL ARRAY k1(m),k2(m),k3(m),k4(m),y1(m),y0_(m),f(m)

FOR x0 = xn TO xe STEP h

s= s + вывод_решения(x0+h,@y1,m) + CHR(10)

s = STR(x,12,2) + CHR(9) ;

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

Const PiNumber = 3.14159265358979

Dim k1(), k2(), k3(), k4(), y1(), y_tmp(), f()

x1 = 0: xn = 2 * PiNumber

ReDim k1(m): ReDim k2(m): ReDim k3(m): ReDim k4(m)

ReDim y1(m): ReDim y_tmp(m): ReDim f(m)

Debug.Print «X Y F(x) ABS(Y-F(x))»

For i = 1 To m Step 2

h = (xn — x1) / (n — 1)

For x0 = x1 To xn Step h

Call правые_части(x0, y(), f, m)

y_tmp(i) = y(i) + k1(i) / 2

Call правые_части(x0 + h / 2, y_tmp, f, m)

y_tmp(i) = y(i) + k2(i) / 2

Call правые_части(x0 + h / 2, y_tmp, f, m)

y_tmp(i) = y(i) + k3(i)

Call правые_части(x0 + h, y_tmp, f, m)

dk = 1 / 6 * (k1(i) + 2 * k2(i) + 2 * k3(i) + k4(i))

Call вывод_решения(x0 + h, y1, m, dec)

Sub правые_части(x, y, f, n)

For i = 1 To 200 Step 2

Sub вывод_решения(x, y, m, dec)

For i = 1 To dec

Debug.Print Format(x, dStr), Format(y(1), dStr), Format(f(x), dStr), Format(Abs(y(1) — f(x)), dStr)

Рис. 1 График решения системы дифференциальных уравнений при , на интервале x = 0 до . Сплошная линия – зависимость , пунктирная — .

Рис.2 Фазовый портрет решения системы дифференциальных уравнений при , на интервале x = 0 до . По горизонтальной оси по вертикальной оси

Контрольные вопросы

1. Какое уравнение называется дифференциальным?

2. Что называется обыкновенным дифференциальным уравнением?

3. Что называется уравнением в частных производных?

4. Что называется порядком дифференциального уравнения?

5. Что называется решением дифференциального уравнения?

6. Что называется интегрированием дифференциального уравнения?

7. Что называется задачей Коши?

8. Что называется краевой задачей?

9. Что называется общим решением дифференциального уравнения?

10. Что такое интегральные кривые?

11. Что такое однородное дифференциальное уравнение первого порядка?

12. Что такое линейное дифференциальное уравнение первого порядка?

13. Запишите вид динамической системы дифференциальных уравнений.

14. Что называется фазовым пространством?

15. Что называется фазовой траекторией?

16. Что называется линейной системой дифференциальных уравнений?

17. Запишите вид линейной системы дифференциальных уравнений с постоянными коэффициентами.

18. Запишите вид характеристического уравнения линейной однородной системы дифференциальных уравнений.

Задания

Решить методами Эйлера, Кранка-Николсона и Рунге-Кутты одно из следующих уравнений вида y’ = F(x,y) на интервале [x,x_n] с начальным условием y(x)=y, принимая h = 0,01.

y’ = F(x,y)	[x,x_n]	y(x)
№1.	y’ = x 2 +y;	[0; 0,2];	y0 = 1.
№2.	y’ = 2x 2 +y;	[0; 0,2];	y0 = 1.
№3.	y’ = 2x+y;	[0; 0,2];	y0 = 1.
№4.	y’ = x+2y;	[0; 0,2];	y0 = 1.
№5.	y’ = x 2 –y;	[1; 1,2];	y0 = 0.
№6.	y’ = x–2y;	[1; 1,2];	y0 = 0.
№7.	y’ = 2(x+y);	[1; 1,2];	y0 = 0.
№8.	y’ = 2x–3y;	[1; 1,2];	y0 = 0.
№9.	y’ = 2x+3y;	[0; 0,2];	y0 = 1.
№10.	y’ = x+3.5 y;	[0; 0,2];	y0 = 1.
№11.	y’ = 4x+y;	[0; 0,2];	y0 = 1.
№12.	y’ = 3x–y;	[0; 0,2];	y0 = 1.
№13.	y’ = 4x–y;	[0; 0,2];	y0 = 1.
№14.	y’ = 1+xy;	[1; 1,5];	y0 = 1.
№15.	y’ = x+y;	[0; 0,5];	y0 = 1.
№16.	y’ = 2x+y;	[0; 0,5];	y0 = 1.
№17.	y’ = 3x+y;	[0; 0,5];	y0 = 1.
№18.	y’ = 4x+y;	[0; 0,5];	y0 = 1.
№19.	y’ = y–2x;	[0; 0,5];	y0 = 1.
№20.	y’ = y–3x;	[0; 0,5];	y0 = 1.
№21.	y’ = x+y 2 ;	[0; 0,5];	y0 = 1.
№22.	y’ = x–y 2 ;	[1; 1,5];	y0 = 0.
№23.	y’ = x–2y 2 ;	[1; 1,5];	y0 = 0.
№24.	y’ = 2x–y 2 ;	[1; 1,5];	y0 = 0.

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

1. Для диссипативной колебательной системы (с силой сопротивления пропорциональной скорости)

при различных значениях коэффициентов k и g.

2. Для колебательной системы с зависящим от времени положением равновесия

при различных значениях коэффициентов k, g и s. Начальные условия и диапазон решения определите самостоятельно.

3. Для диссипативной колебательной системы с внешней силой

при различных значениях коэффициентов k, g, r и s. Начальные условия и диапазон решения определите самостоятельно.

Организация стока поверхностных вод: Наибольшее количество влаги на земном шаре испаряется с поверхности морей и океанов (88‰).

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

Опора деревянной одностоечной и способы укрепление угловых опор: Опоры ВЛ — конструкции, предназначенные для поддерживания проводов на необходимой высоте над землей, водой.

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

Рассмотрим задачу Коши для системы обыкновенных дифференциальных уравнений $$ \begin \tag <1>\frac &= f_i (t, u_1, u_2, \ldots, u_n), \quad t > 0\\ \tag <2>u_i(0) &= u_i^0, \quad i = 1, 2, \ldots, m. \end $$

Используя векторные обозначения, задачу (1), (2) можно записать как задачу Коши $$ \begin \tag <3>\frac> &= \pmb(t, \pmb), \quad t > 0, \\ \tag <4>\pmb(0) &= \pmb_0 \end $$ В задаче Коши необходимо по известному решению в точке $ t = 0 $ необходимо найти из уравнения (3) решение при других $ t $.

Численные методы решения задачи Коши

Существует большое количество методов численного решения задачи (3), (4). Вначале рассмотрим простейший явный метод Эйлера и его программную реализацию. Затем будут представлены методы Рунге—Кутта и многошаговые методы.

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

Идея численных методов решения задачи (3), (4) состоит из четырех частей:

1. Вводится расчетная сетка по переменной $ t $ (время) из $ N_t + 1 $ точки $ t_0 $, $ t_1 $, $ \ldots $, $ t_ $. Нужно найти значения неизвестной функции $ \pmb $ в узлах сетки $ t_n $. Обозначим через $ \pmb^n $ приближенное значение $ \pmb(t_n) $.

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

3. Аппроксимируем производные конечными разностями.

4. Формулируем алгоритм, который вычисляет новые значения $ \pmb^ $ на основе предыдущих вычисленных значений $ \pmb^k $, $ k 0 $ при $ \tau\to 0 $.

Явный метод Эйлера

Проиллюстрируем указанные шаги. Для начала введем расчетную сетку. Очень часто сетка является равномерной, т.е. имеет одинаковое расстояние между узлами $ t_n $ и $ t_ $: $$ \omega_\tau = \< t_n = n \tau, n = 0, 1, \ldots, N_t \>. $$

Затем, предполагаем, что уравнение выполнено в узлах сетки, т.е.: $$ \pmb^\prime (t_n) = \pmb(t_n, u(t_n)), \quad t_n \in \omega_\tau. $$

Заменяем производные конечными разностями. С этой целью, нам нужно знать конкретные формулы, как производные могут быть аппроксимированы конечными разностями. Простейший подход заключается в использовании определения производной: $$ \pmb^\prime(t) = \lim_ <\tau \to 0>\frac<\pmb(t+\tau) — \pmb(t)><\tau>. $$

В произвольном узле сетки $ t_n $ это определение можно переписать в виде: $$ \begin \pmb^\prime(t_n) = \lim_ <\tau \to 0>\frac<\pmb(t_n+\tau) — \pmb(t_n)><\tau>. \end $$ Вместо того, чтобы устремлять шаг сетки к нулю, мы можем использовать малый шаг $ \tau $, который даст численное приближение $ u^\prime(t_n) $: $$ \begin \pmb^\prime(t_n) \approx \frac<\pmb^ — \pmb^><\tau>. \end $$ Такая аппроксимация известна как разностная производная вперед и имеет первый порядок по $ \tau $, т.е. $ O(\tau) $. Теперь можно использовать аппроксимацию производной. Таким образом получим явный метод Эйлера: $$ \begin \tag <5>\frac<\pmb^ — \pmb^n> <\tau>= \pmb(t_n, \pmb^). \end $$

Четвертый шаг заключается в получении численного алгоритма. Из (5) следует, что мы должны знать значение $ y^n $ для того, чтобы решить уравнение (5) относительно $ y^ $ и получить формулу для нахождения приближенного значения искомой функции на следующем временном слое $ t_ $: $$ \begin \tag <6>\pmb^ = \pmb^n + \tau \pmb(t_n, \pmb^) \end $$

При условии, что у нас известно начальное значение $ \pmb^0 = \pmb_0 $, мы можем использовать (6) для нахождения решений на последующих временных слоях.

Программная реализация явного метода Эйлера

Выражение (6) может быть как скалярным так и векторным уравнением. И в скалярном и в векторном случае на языке Python его можно реализовать следующим образом

При решении системы (векторный случай), u[n] — одномерный массив numpy длины $ m+1 $ ($ m $ — размерность задачи), а функция F должна возвращать numpy -массив размерности $ m+1 $, t[n] — значение в момент времени $ t_n $.

Таким образом численное решение на отрезке $ [0, T] $ должно быть представлено двумерным массивом, инициализируемым нулями u = np.zeros((N_t+1, m+1)) . Первый индекс соответствует временному слою, а второй компоненте вектора решения на соответствующем временном слое. Использование только одного индекса, u[n] или, что то же самое, u[n, :] , соответствует всем компонентам вектора решения.

Функция euler решения системы уравнений реализована в файле euler.py:

Строка F_ = lambda . требует пояснений. Для пользователя, решающего систему ОДУ, удобно задавать функцию правой части в виде списка компонент. Можно, конечно, требовать чтобы пользователь возвращал из функции массив numpy , но очень легко осуществлять преобразование в самой функции решателе. Чтобы быть уверенным, что результат F будет нужным массивом, который можно использовать в векторных вычислениях, мы вводим новую функцию F_ , которая вызывает пользовательскую функцию F «прогоняет» результат через функцию assaray модуля numpy .

Неявный метод Эйлера

При построении неявного метода Эйлера значение функции $ F $ берется на новом временном слое, т.е. для решении задачи (5) используется следующий метод: $$ \begin \tag <7>\frac<\pmb^ — \pmb^n> <\tau>= \pmb(t_, \pmb^). \end $$

Таким образом для нахождения приближенного значения искомой функции на новом временном слое $ t_ $ нужно решить нелинейное уравнение относительно $ \pmb^ $: $$ \begin \tag <8>\pmb^ — \tau \pmb(t_, \pmb^) — y^n = 0. \end $$

Для решения уравнения (8) можно использовать, например, метод Ньютона.

Программная реализация неявного метода Эйлера

Функция backward_euler решения системы уравнений реализована в файле euler.py:

Отметим, что для нахождения значения u[n+1] используется функция fsolve модуля optimize библиотеки scipy . В качестве начального приближения для решения нелинейного уравнения используется значение искомой функции с предыдущего слоя u[n] .

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

Одношаговый метод Рунге—Кутта в общем виде записывается следующим образом: $$ \begin \tag <9>\frac<\pmb^ — \pmb^n> <\tau>= \sum_^s b_i \pmb_i, \end $$ где $$ \begin \tag <10>\pmb_i = \pmb\left( t_n + c_i\tau, \pmb^n + \tau \sum_^s a_\pmb_j \right), \quad i = 1, 2, \ldots, s. \end $$ Формула (9) основана на $ s $ вычислениях функции $ \pmb $ и называется $ s $-стадийной. Если $ a_ = 0 $ при $ j \geq i $ имеем явный метод Рунге—Кутта. Если $ a_ = 0 $ при $ j > i $ и $ a_ \ne 0 $, то $ \pmb_i $ определяется неявно из уравнения $$ \begin \tag <11>\pmb_i = \pmb\left( t_n + c_i\tau, \pmb^n + \tau \sum_^ a_\pmb_j + \tau a_ \pmb_i \right), \quad i = 1, 2, \ldots, s. \end $$ О таком методе Рунге—Кутта говорят как о диагонально-неявном.

Одним из наиболее распространенных является явный метод Рунге-Кутта четвертого порядка: $$ \begin \tag <12>\pmb_1 & = \pmb(t_n, \pmb^n), &\quad \pmb_2 &= \pmb\left( t_n + \frac<\tau><2>, \pmb^n + \tau \frac<\pmb_1> <2>\right),\\ \pmb_3 &= \pmb\left( t_n + \frac<\tau><2>, \pmb^n + \tau \frac<\pmb_2> <2>\right), &\quad \pmb_4 &= \pmb\left( t_n + \tau, \pmb^n + \tau \pmb_3 \right),\\ \frac<\pmb^ -\pmb^n> <\tau>&= \frac<1> <6>(\pmb_1 + 2\pmb_2 + 2\pmb_3 + \pmb_4) & & \end $$

Многошаговые методы

В методах Рунге—Кутта в вычислениях участвуют значения приближенного решения только в двух соседних узлах $ \pmb^n $ и $ \pmb^ $ — один шаг по переменной $ t $. Линейный $ m $-шаговый разностный метод записывается в виде $$ \begin \tag <13>\frac<1> <\tau>\sum_^m a_i \pmb^ = \sum_^ b_i \pmb(t_, \pmb^), \quad n = m-1, m, \ldots \end $$ Вариант численного метода определяется заданием коэффициентов $ a_i $, $ b_i $, $ i = 0, 1, \ldots, m $, причем $ a_0 \ne 0 $. Для начала расчетов по рекуррентной формуле (13) необходимо задать $ m $ начальных значений $ \pmb^0 $, $ \pmb^1 $, $ \dots $, $ \pmb^ $ (например, можно использовать для их вычисления метод Эйлера).

Различные варианты многошаговых методов (методы Адамса) решения задачи с начальными условиями для систем обыкновенных дифференциальных уравнений могут быть получены на основе использования квадратурных формул для правой части равенства $$ \begin \tag <14>\pmb(t_) — \pmb(t_n) = \int_^> \pmb(t, \pmb) dt \end $$

Для получения неявного многошагового метода используем для подынтегральной функции интерполяционную формулу по значениям функции $ \pmb^ = \pmb(t_, \pmb^) $, $ \pmb^n $, $ \dots $, $ \pmb^ $, т.е. $$ \begin \tag <15>\frac<\pmb^ — \pmb^n> <\tau>= \sum_^ b_i \pmb(t_, \pmb^) \end $$

Для интерполяционного метода Адамса (15) наивысший порядок аппроксимации равен $ m+1 $.

Для построения явных многошаговых методов можно использовать процедуру экстраполяции подынтегральной функции в правой части (14). В этом случае приближение осуществляется по значениям $ \pmb^n $, $ \pmb^ $, $ \dots $, $ \pmb^ $ и поэтому $$ \begin \tag <16>\frac<\pmb^ — \pmb^n> <\tau>= \sum_^ b_i \pmb(t_, \pmb^) \end $$

Для экстраполяционного метода Адамса (16) погрешность аппроксимации имеет $ m $-ый порядок.

На основе методов Адамса строятся и схемы предиктор–корректор. На этапе предиктор используется явный метод Адамса, на этапе корректора — аналог неявного метода Адамса. Например, при использовании методов третьего порядка аппроксимации в соответствии с (18) для предсказания решения положим $$ \frac<\pmb^ — \pmb^n> <\tau>= \frac<1> <12>(23 \pmb^ -16\pmb^ + 5\pmb^). $$ Для уточнеия решения (см. (17)) используется схема $$ \frac<\pmb^ — \pmb^n> <\tau>= \frac<1> <24>(9\pmb^ + 19\pmb^ — 5\pmb^ + \pmb^). $$ Аналогично строятся и другие классы многошаговых методов.

Жесткие системы ОДУ

При численном решении задачи Коши для систем обыкновенных дифференциальных уравнений (3), (4) могут возникнуть дополнительные трудности, порожденные жесткостью системы. Локальные особенности поведения решения в точке $ u = w $ передаются линейной системой $$ \begin \frac

= \sum_^ \frac<\partial f_i> <\partial u_j>(t, w) v + \bar(t), \quad t > 0. \end $$
Пусть $ \lambda_i(t) $, $ i = 1, 2, \ldots, m $ — собственные числа матрицы $$ \begin A(t) = \< a_(t) \>, \quad a_(t) = \frac<\partial f_i><\partial u_j>(t, w). \end $$ Система уравнений (3) является жесткой, если число $$ \begin S(t) = \frac <\max_<1 \leq i \leq m>|Re \lambda_i(t)|> <\min_<1 \leq i \leq m>|Re \lambda_i(t)|> \end $$ велико. Это означает, что в решении присутствуют составляющие с сильно различающимися масштабами изменения по переменной $ t $.

Для численное решения жестких задач используются вычислительные алгоритмы, которые имеют повышенный запас устойчивости. Необходимо ориентироваться на использование $ A $-устойчивых или $ A(\alpha) $-устойчивых методов.

Метод называется $ A $-устойчивым, если при решении задачи Коши для системы (3) область его устойчивости содержит угол $$ \begin |\arg(-\mu)|

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

Читайте также:

I. Постановка задачи

I. Цели и задачи БЖД. Основные понятия.

II СИТУАЦИОННЫЕ ЗАДАЧИ ПО ЧАСТНОЙ ГИСТОЛОГИИ

II. Основные задачи службы

II. Решение задачи с ограничениями.

II. Цели и задачи курса

III. Объемно-планировочное решение здания.

III. Основные проблемы и их решение

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

XII. Заполните схему. Выпишите задачи физического воспитания учащихся и воспитательные дела, предназначенные для их решения.

Автокорреляция в остатках. Критерий Дарбина-Уотсона в оценке качества уравнений, построенных по временным рядам

Адаптивный подход в сельском хозяйстве и задачи агроэкологии

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

Пусть на отрезке [a, b] требуется найти численное решение задачи Коши y′ = dy(t)/dt = f (y, t)

с соответствующим начальным условием

, где a = t . Как и в предыдущем методе разобьем этот участок на n равных частей и построим последовательность значений t , t ₁ , t ₂ , …, t _n аргумента t искомой функции y(t). Предполагаем существование непрерывных производных функции y(t) до пятого порядка.

Формулу для решения можно записать в виде:

где ∆y _k — приращение искомой функции y(t) на (k+1)-ом шаге интегрирования.

Придадим аргументу t приращение, равное шагу интегрирования h, и разложим функцию y(t + h) в ряд Тейлора в окрестности точки t, сохранив в нем пять членов:

Перенося первое слагаемое в этой сумме в левую часть получим, что

Здесь производные определяются последовательным дифференцированием уравнения y = y (t)

Вместо непосредственных вычислений по формуле (10.6) в методе Рунге-Кутта для каждого значения ∆y _k = ∆ y(x _k ) определяются четыре числа:

Доказывается, что если числа k 1k , k 2k , k 3k , k 4k , последовательно умножить на 1/6, 1/3, 1/3, 1/6 и сложить между собой, то выражение (10.10), соответствующее формуле Рунге-Кутта:

имеет погрешность R _n (h 5 ).

Таким образом, рабочая формула Рунге-Кутта для интегрирования имеет вид:

Метод Рунге-Кутта может быть использован и при решении систем дифференциальных уравнений.

Дата добавления: 2015-04-24 ; Просмотров: 635 ; Нарушение авторских прав? ;

Нам важно ваше мнение! Был ли полезен опубликованный материал? Да | Нет

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

МИНИСТЕРСТВО ОБРАЗОВАНИЯ И НАУКИ УКРАИНЫ
ГОСУДАРСТВЕННОЕ ВЫСШЕЕ УЧЕБНОЕ ЗАВЕДЕНИЕ
УКРАИНСКИЙ ГОСУДАРСТВЕННЫ Й
студент группы 3ИC27 Куделя С.В.
ассистент Ильхман Яна Викторовна
1. ОБЪЕКТНО-ОРИЕНТИРОВАННОЕ ПРОГРАММИРОВАНИЕ
2. ОПИСАНИЕ ПРЕДМЕТНОЙ ОБЛАСТИ (ПО)
2.1 Назначение программного продукта
2.2 Основные задачи
2.3 Входные и выходные данные
3.1 Выделение основных объектов ПО
3.2 Описание полей и методов
3.3 Иерархия классов на основе выделенных объектов
4. ОСНОВНЫЕ ФОРМЫ И КОМПОНЕНТЫ, ИСПОЛЬЗУЕМЫЕ ДЛЯ РЕАЛИЗАЦИИ ПРОГРАММЫ. ОСНОВНЫЕ АЛГОРИТМИЧЕСКИЕ РЕШЕНИЯ РЕЗУЛЬТАТЫ РАБОТЫ ПРОГРАММЫ
4.1 Метод Рунге-Кутта
4.2 Описание программы ” РЕШЕНИЕ ОДУ “
4.3 Назначение элементов графического окна программы
4.4 Реакция программы при возникновении ошибок
4.5 Перечень компонент DELPHI использованных в программе
5. ТЕХНИЧЕСКИЕ ХАРАКТЕРИСТИКИ И ТРЕБОВАНИЯ К ПО
6. ТЕКСТ ПРОГРАММЫ
7. РЕШЕНИЕ ДИФФЕРЕНЦИАЛЬНОГО УРАВНЕНИЯ Y = Y−2X/Y МЕТОДОМ РУНГЕ – КУТТА В СРЕДЕ EXCEL
1. ОБЪЕКТНО-ОРИЕНТИРОВАННОЕ ПРОГРАММИРОВАНИЕ
Delphi является объектно-ориентированной средой программирования. В качестве языка программирования используется язык Object Pascal.
Исторически сложилось так, что программирование возникло и развивалось как процедурное программирование, которое предполагает, что основой программы является алгоритм, процедура обработки данных.
Объектно-ориентированное программирование (ООП) — это методика разработки программ, в основе которой лежит понятие объект. Объект — это некоторая структура, соответствующая объекту реального мира, его поведению. Задача, решаемая с использованием методики ООП, описывается в терминах объектов и операций над ними, а программа при таком подходе представляет собой набор объектов и связей между ними. Объектно-ориентированное программирование позволяет конструировать новые и производные (дочерние) классы на основе существующих классов.
По сравнению с традиционными способами программирования ООП обладает рядом преимуществ. Главное из них заключается в том, что эта концепция в наибольшей степени соответствует внутренней логике функционирования операционной системы (ОС) Windows. Программа, состоящая из отдельных объектов, отлично приспособлена к реагированию на события, происходящие в ОС. К другим преимуществам ООП можно отнести большую надежность кода и возможность повторного использования отработанных объектов.
Основные понятия ООП в языке Delphi: объект, класс, компонент;
Основные механизмы ООП: инкапсуляция, наследование и полиморфизм;
Класс — помимо описания данных, включает описание процедур и функций, которые могут быть выполнены над представителем класса — объектом. Данные класса называются полями, процедуры и функции — методами. Переменные в зависимости от предназначения именуются полями или свойствами.
При создании объекта он наследует структуру (переменные) и поведение (методы) своего класса.
Класс, называемый потомком, производным или дочерним классом (подклассом), также может быть создан на основе другого родительского класса (предка) и при этом наследует его структуру и поведение.
Методы — это процедуры и функции, описанные внутри класса и предназначенные для операций над его полями. В состав класса входит указатель на специальную таблицу, где содержится вся информация, нужная для вызова методов. От обычных процедур и функций методы отличаются тем, что им при вызове передается указатель на тот объект, который их вызвал. Поэтому обрабатываться будут поля именно того объекта, который вызвал метод. Внутри метода указатель на вызвавший его объект доступен под зарезервированным именем self.
Понятие свойства можно определить его как поле, доступное для чтения и записи не напрямую, а через соответствующие методы.
Классы могут быть описаны либо в секции интерфейса модуля, либо на верхнем уровне вложенности секции реализации. Не допускается описание классов внутри процедур и других блоков кода.
Любой компонент (элемент управления) или объект в Delphi всегда является экземпляром класса.
Программно объект представляет собой переменную объектного типа.
Для каждого компонента Delphi существует свой класс, наследуемый от TComponent. Предком всех объектов, включая компоненты, является класс TObject.
Инкапсуляция — это создание защищенных объектов, доступ к свойствам и методам которых разрешен только через определенные разработчиком . Иначе говоря, инкапсуляция — это предоставление разработчику конкретного набора свойств и методов для управления поведением и свойствами объекта, определяемыми внутри класса.
Инкапсуляция обеспечивает скрытие полей объекта с целью обеспечения доступа к ним только посредством методов класса.
В языке Delphi ограничение доступа к полям объекта реализуется при помощи свойств объекта. Свойство объекта характеризуется полем, сохраняющим значение свойства, и двумя методами, обеспечивающими доступ к полю свойства. Метод установки значения свойства называется методом записи свойства (write), а метод получения значения свойства — методом чтения свойства (read).
В описании класса перед именем свойства записывают слово property (свойство). Ключевое слово property обозначает свойство объекта, которое в отличие от обычных полей (переменных) класса имеет спецификаторы доступа обеспечивающие контроль доступа к свойствам объекта.
После имени свойства указывается его тип, затем — имена методов, обеспечивающих доступ к значению свойства. После слова read указывается имя метода, обеспечивающего чтение свойства, после слова write — имя метода, отвечающего за запись свойства.
Пример описания класса TPerson, содержащего два свойства: Name и Address:
TName = string [15]; TAddress = string[35];
TPerson = class // класс
FName: TName; // значениесвойства Name
FAddress: TAddress; // значениесвойства Address
Function GetName: TName;
Function GetAddress: TAddress;
Property Name: Tname // свойство Name
read GetName; // доступно только для чтения
Property Address: TAddress // свойство Address

read GetAddress // доступнодлячтения
write SetAddress; // изаписи
В программе для установки значения свойства записать обычную инструкцию присваивания значения свойству. Например, чтобы присвоить значение свойству Address объекта student, достаточно записать
student.Address := ‘Гвардейский, ул.Зенитная 1, кв.10’;
Внешне применение свойств ничем не отличается от использования полей объекта. Однако между свойством и полем объекта существует принципиальное отличие: при присвоении и чтении значения свойства автоматически вызывается процедура, которая выполняет некоторую работу.
Наследование позволяет определять новые классы в терминах существующих классов.
Новые классы возможно определять посредством добавления полей, свойств и методов к уже существующим классам. Такой механизм получения новых классов называется порождением. При этом новый, порожденный класс (потомок) наследует свойства и методы своего базового, родительского класса.
В объявлении класса-потомка указывается класс родителя. Например, класс TEmployee (сотрудник) может быть порожден от рассмотренного выше класса TPerson путем добавления поля FDepartment (отдел). Объявление класса TEmplioyee в этом случае может выглядеть так:
TEmployee = class (TPerson)
FDepartment: integer; // номеротдела
constructor Create(Name:TName; Dep:integer);
Заключенное в скобки имя класса TPerson показывает, что класс TEmployee является производным от класса TPerson. В свою очередь, класс TPerson является базовым для класса TEmployee.
Класс TEmpioyee должен иметь свой собственный конструктор, обеспечивающий инициализацию класса-родителя и своих полей.
Пример реализации конструктора класса TEmployee:
В приведенном примере директивой inherited вызывается конструктор родительского класса. После этого присваивается значение полю класса-потомка.
После создания объекта производного класса в программе можно использовать поля и методы родительского класса. Ниже приведен фрагмент программы, демонстрирующий эту возможность.
engineer.address := ‘ул.Блохина, д.8, кв.10’;
Первая инструкция создает объект типа TEmployee, вторая — устанавливает значение свойства, которое относится к родительскому классу.
Полиморфизм — это возможность различных объектов реагировать по-разному на одни и те же события, иначе говоря, возможность использовать одинаковые имена для методов, входящих в различные классы. Концепция полиморфизма обеспечивает в случае применения метода к объекту использование именно того метода, который соответствует классу объекта. Синтаксис языка поддерживает общепринятую для объектно-ориентированного программирования нотацию: имя_объекта: свойство для ссылки на свойство объекта или имя_объекта: метод для вызова метода объекта. При этом имя_объекта может содержать как простую, так и квалифицированную ссылку. Квалифицированная ссылка содержит отделенные точкой имена родительских объектов.
Пусть определены три класса, один из которых является базовым для двух других:
// базовыйкласс TPerson = class
fname: string; // имя
function info: string;
// производныйот TPerson TStud = class(TPerson)
constructor Create(name:string; gr:integer);
function info: string; override; end;
// производныйот TPerson TProf = class(TPerson)
function info: string;
В каждом из этих классов определен метод info. В базовом классе при помощи директивы virtual метод info объявлен виртуальным. Объявление метода виртуальным дает возможность дочернему классу произвести замену виртуального метода своим собственным. В каждом дочернем классе определен свой метод info, который замещает соответствующий метод родительского класса (метод порожденного класса, замещающий виртуальный метод родительского класса, помечается директивой override).
Определение метода info для каждого класса:
result := fname + ‘ гp.’ + IntTostr(fgr);
result := fname + ‘ каф.’ + fdep;
Так как оба класса порождены от одного и того же базового, объявить список студентов и преподавателей можно так (следует помнить, что объект — это указатель):
list: array [l. .SZL] of TPerson;
Объявить подобным образом список можно потому, что язык Delphi позволяет указателю на родительский класс присвоить значение указателя на дочерний класс. Поэтому элементами массива list могут быть как объекты класса TStud, так и объекты класса TProf.
Вывести список студентов и преподавателей можно применением метода info к элементам массива. Например, так:
for i:=l to SZL do // SZL — размер массива-списка
then st := st + list[i].Info
+ #13; ShowMessage (st);
Во время работы программы каждый элемент массива может содержать как объект типа xstud, так и объект типа TProf. Концепция полиморфизма обеспечивает применение к объекту именно того метода, который соответствует типу объекта.
Есть еще одна, совершенно особенная разновидность методов — перегружаемые.
Перегрузка методов нужна, чтобы произвести одинаковые или похожие действия с разнотипными данными.
Пример, иллюстрирующий статические методы:
procedure SetData(AValue: Extended);
procedure SetData(AValue: Integer); end;
В этом случае попытка вызова из объекта Т2 методов
вызовет ошибку компиляции на первой из двух строк. Для компилятора внутри Т2 статический метод с параметром типа extended перекрыт, и он его «не признает».
Для выхода из сложившегося положения можно переименовать один из методов, например, создать SetlntegerData и SetExtendedData, но если методов не два, а, например, сто, моментально возникнет путаница. Сделать методы виртуальными нельзя, поскольку тип и количество параметров в одноименных виртуальных методах должны в точности совпадать. Для разрешения этой ситуации существуют перегружаемые методы, объявляемые при помощи директивы overload:
procedure SetData(AValue: Extended);overload;
procedure SetData(AValue: Integer); overload;
Объявив метод SetData перегружаемым, в программе можно использовать обе его реализации одновременно. Это возможно потому, что компилятор определяет тип передаваемого параметра (целый или с плавающей точкой) и в зависимости от этого подставит вызов соответствующего метода: для целочисленных данных — метод объекта T2ndobj, для данных с плавающей точкой — метод объекта Tistobj.
Можно перегрузить и виртуальный (динамический) метод. Надо только в этом случае добавить директиву reintroduce:
procedure SetData(AValue: Extended); overload; virtual;
procedure SetData(AValue: Integer); reintroduce; overload ;
На перегрузку методов накладывается ограничение — нельзя перегружать методы, находящиеся в области видимости published, т. е. те, которые будут использоваться в Инспекторе объектов.
В модели объектов языка Delphi существует механизм доступа к составным частям объекта, определяющий области, где ими можно пользоваться (области видимости). Поля и методы могут относиться к четырем группам (секциям), отличающимся областями видимости. Методы и свойства могут быть общими (секция public ), личными (секция private ), защищенными (секция protected ) и опубликованными (секция published ). Есть еще и пятая группа, automated , она ранее использовалась для создания объектов СОМ; теперь она присутствует в языке только для обратной совместимости с программами на Delphi версий 3—5.
Области видимости, определяемые первыми тремя директивами, таковы.
· Поля, свойства и методы секции public не имеют ограничений на видимость. Они доступны из других функций и методов объектов, как в данном модуле, так и во всех прочих, ссылающихся на него.
· Поля, свойства и методы, находящиеся в секции private , доступны только в методах класса и в функциях, содержащихся в том же модуле, что и описываемый класс. Такая директива позволяет полностью скрыть детали внутренней реализации класса. Свойства и методы из секции private можно изменять, и это не будет сказываться на программах, работающих с объектами этого класса. Единственный способ для кого-то другого обратиться к ним — переписать заново созданный вами модуль (если, конечно, доступны исходные тексты).
· Поля, свойства и методы секции protected также доступны только внутри модуля с описываемым классом. Но — и это главное — они доступны в классах, являющихся потомками данного класса, в том числе и в других модулях. Такие элементы особенно необходимы для разработчиков новых компонентов — потомков уже существующих. Оставляя свободу модернизации класса, они все же скрывают детали реализации от того, кто только пользуется объектами этого класса.
СОЗДАНИЕ НОВОГО КЛАССА
Для того чтобы создать новый класс, в interface-секции кода модуля следует записать:
Каждая форма в проекте, разрабатываемом в Delphi, описывается отдельным модулем (создаваемым автоматически при создании новой формы). Этот модуль описывает новый класс для компонента Form. Первоначально по умолчанию создается класс TForml, наследуемый от класса TForm из VCL-библиотеки. Это автоматически записывается в модуле следующим образом:
type (Объявление класса>
[Объявление private переменных и методов>
Объявление переменных и методов класса

Переменные класса указываются после модификаторов доступа (public, private , protected, published, automated ), определяющих их область видимости. Свойства, указанные после модификатора доступа published, являются общедоступными и отображаются в инспекторе объектов.
После имени переменной или списка имен, разделенных через запятую, указывается символ : и тип переменной. Типом может быть как базовый тип Delphi (например. Integer, Boolean), так и производный тип, в том числе реализуемый, как некоторый класс. Такой тип иногда называется объектным типом.
При объявлении методов класса перед именем метода указывается ключевое слово function или procedure. Для функций также после имени функции через символ указывается тип возвращаемого значения.
Var2, Var3: TVarTypeClass;
function F1: Integer;
Объявление класса содержит только объявление переменных и методов. Реализация методов — функций и процедур — записывается в implementation-секции модуля.
Каждый модуль, создаваемый на основе разрабатываемой формы, представляет собой описание класса. Как правило, производного от класса TForm. Любой компонент, располагаемый в форме, также является экземпляром некоторого класса.
Классы в Delphi образуют иерархическое дерево. Будем называть классы из VCL-библиотеки Delphi базовыми классами. Иерархическое дерево для некоторого класса любого компонента имеет корневым элементом класс TObject. Просмотреть иерархию классов-потомков можно в окне Exploring. Для того чтобы перейти в него, достаточно выполнить команду меню View|Browser или нажать клавиши Shift+CtrL+B.
На рис.1 представлена страница Classes окна Exploring Classes. На ней отображено иерархическое дерево наследования для класса TForm 1. В правой части окна расположена панель, содержащая три страницы — Scope, Inheritance, References. Страница Scope содержит древовидную диаграмму всех объектов, переменных и методов выделенного на левой панели класса. При этом ветвь Inherited содержит имена класса-предка и класса-потомка. Страница Inheritance содержит поддерево иерархии классов, начиная с класса-предка для выделенного на левой панели класса.
На странице References можно узнать имена всех модулей и номера строк, в которых встречается имя данного класса.
Для того чтобы добавить в проект собственные описания производных классов, наиболее целесообразно создать отдельный модуль и в interface-секции модуля записать все объявления классов.
Все классы VCL-библиотеки Delphi разбиты на группы, которые расположены в каталоге Delphi7\Source\VCL . Для того чтобы просмотреть файл библиотеки, достаточно выполнить File | Open и выбрать каталог и имя файла. Справа в окне кода программы (рис.2) будет показан код модуля, а слева — список всех объявленных в нем классов.
СВОЙСТВА / МЕТОДЫ И ОБРАБОТЧИКИ СОБЫТИЙ
Каждый объект обладает набором свойств. Свойства могут быть как наследуемые от родительского класса, так и добавленные индивидуально для создаваемого объекта. Список всех свойств объекта и их значений отображается в диалоговом окне Object Inspector.
Ссылка на свойство в программном модуле записывается как
Метод — это процедура или функция, ассоциируемая с некоторым объектом.
Ссылка на методов программном модуле записывается как
Delphi-приложение выполняется в среде Windows, и как любое Windows-приложение, получает сообщения о возникающих для него событиях. Управление приложением фактически сводится к обработке получаемых сообщений.
Методы, в которых содержится код обработки события, называются обработчиками событий (Event Handler). Delphi автоматически генерирует процедуры обработки событий – обработчикисобытий для любого компонента. При этом имя обработчика событий формируется из имени компонента и названия события (например, EditlClick). Имя обработчика события автоматически квалифицируется именем класса формы.
Например: TForml.ButtonlClick(Sender: TObject);.
Для каждого компонента предусмотрено одно стандартное событие. Например, для командной кнопки, флажка, списка, поля ввода — это событие Click, а для формы — событие FormCreate.
Для того чтобы автоматически добавить в модуль объявление и описание обработчика стандартного события, достаточно выполнить на компоненте формы или самой форме двойной щелчок мышью. Объявление события добавляется в interface-секцию модуля, а пустое описание события — в implementation-секцию модуля. Далее в редакторе кода внутри уже имеющегося блока begin end; следует только ввести код обработчика события.
procedure TForml.ButtonlClick(Sender: TObject);
Для того чтобы добавить обработчик любого события, можно выбрать в инспекторе объектов страницу Events и выполнить двойной щелчок мышью в поле, расположенном справа от имени события. Это поле представляет собой комбинированное окно списка — в него можно вводить новое значение имени обработчика события двойным щелчком мыши или выбирать имя уже существующей процедуры. Это позволяет при необходимости определять одну процедуру обработки событий одновременно для нескольких событий.
2. ОПИСАНИЕ ПРЕДМЕТНОЙ ОБЛАСТИ (ПО)
2.1 Назначение программного продукта
Программа предназначена для решения заданных программно обыкновенных дифференциальных уравнений первого порядка методом Рунге – Кутта, вывода результата решения ОДУ на экран в виде графика в декартовой системе координат.
2.2 Основные задачи
Программа обеспечивает решение следующих задач:
− ввод исходных данных;
− решение ОДУ и вывод результата решения в численном и аналитическом виде;
− обнуление результатов решения ОДУ;
− контроль корректности ввода исходных данных и вывод на экран сообщение о содержании
ошибки с рекомендацией по её устранению;
− контроль возникновения в процессе вычислений ошибки «деление на ‘0’» и вывод на экран
соответствующего сообщения о содержании ошибки с рекомендацией по её устранению.
2.3 Входные и выходные данные
Входными данными для программы являются:
− начальное условие решения ОДУ − у'(х ) = у .,
− начальное и конечное значения отрезка, в пределах которого находится решение ОДУ;
−величина шага дифференцирования,
−формула образцовой функции.
Выходными данными программы являются:
− массив (х₁ ,у₁ ; х₂ ,у₂ ;…; х_i ,у_i ) − решений выбранного дифференциального уравнения на заданном
− график функции, которая, будучи подставленной, в исходное образцовое уравнение, обращает его в
тождество и одновременно удовлетворяет начальному условию.
3.1 Выделение основных объектов ПО
− Объект класса TRKutta (Form1) – главная окно программы.
− Объект класса TRngeKutta (Form2) – окно вывода графика функции-решения ДУ.
− Объект класса TSpravka(Form3) – окно «О программе».
− Объект класса TRungeKutta – координатная плоскость и график функции;
3.2 Описание полей и методов
x1 — значение x1(начало отрезка);
x2 — значение x2(конец отрезка);
yc — начальные значения (Y0) для передачи в графический модуль;
xc — начальные значения (х0)для передачи в графический модуль;
h– зачение величины шага вычислений;
f- значение функции при начальных условиях
zx — массив значений аргумента;
zy — массив значений функции;
line_arr — размерность массивa;
procedureButton1Click – вычисление значений функции – решения ОДУ;
procedureButton2Click – очистка полей ввода/вывода данных;
procedure Button3Click — выводокнаГРАФИК;
procedure Button4Click — выходизпрограммы;
procedure RadioGroup1Click — выборобразцовойфункции;

procedure Button5Click — активациявводаобразцовыхфункций;
procedure Button6Click — деактивациявводаобразцовыхфункций;
procedureN7Click — вывод информации о программе;
Класс TRungeKutta ( родительский классTObject )
x0,y0 — координаты начала координатных осей;
dx,dy- шаг координатной сетки (в пикселях)
y1,y2,xP,xL — засечки по оси Y и X
dlx,dly — шаг меток (оцифровки) линий сетки по X и Y
cross- счетчик неоцифрованных линий сетки
dcross- количество неоцифрованных линий между оцифрованными
razm– размер массивов;
uzx- Динамический массив координат-X
uzy- Динамический массив координат-Y
uxc,uyc — Оцифровка по осям
mx, my- масштабы по осям X и Y;
BaseMasht_X,BaseMasht_Y — МАХ значения элементов в массивах
functionMaxAbsElementArray– определяет MAX элемент массива по модулю;
procedureUstanovkaParametrov – вычисляет исходные параметры необходимые для построения декартовой плоскости;
procedureKoordPloskost – вычерчивает координатную плоскость;
constructorTRungeKutta.CreateGr — создание объекта (график функции, координатная сетка, координатные оси)
RungeKutta — переменная-объект классаTRungeKutta ;
procedureN4Click — вывод графика функции в окне ‘График’;
procedureN5Click — закрытие окна ‘График’;
procedureButton1Click – процедура вывода информации о программе
3.3 Иерархия классов на основе выделенных объектов
4. ОСНОВНЫЕ ФОРМЫ И КОМПОНЕНТЫ, ИСПОЛЬЗУЕМЫЕ ДЛЯ РЕАЛИЗАЦИИ ПРОГРАММЫ. ОСНОВНЫЕ АЛГОРИТМИЧЕСКИЕ РЕШЕНИЯ РЕЗУЛЬТАТЫ РАБОТЫ ПРОГРАММЫ
1. Решение дифференциальных уравнений
Пусть дано дифференциальное уравнение первого порядка:
и начальное условие его решения:
Тогда решить уравнение — это значит найти такую функцию у — φ(х), которая, будучи подставленной, в исходное уравнение, обратит его в тождество и одновременно будет удовлетворено начальное условие. Задача отыскания функции у = φ (х) называется в математике задачей Коши. При решении дифференциального уравнения порядка nзадача Коши формулируется следующим образом.
Дано дифференциальное уравнение порядка n:
у ( n ) = f(x, y, у’ ’ ,…,y n -1 )
Необходимо найти такую функцию у = φ (х), которая, будучи подставленной в исходное уравнение, обратит его в тождество и одновременно будут удовлетворены следующие п начальных условий:
4.1 Метод Рунге-Кутта
Метод Рунге-Кутта обладает более высокой точностью, чем методы Эйлера за счет снижения методических ошибок. Идея метода состоит в следующем.
По методу Эйлера решение дифференциального уравнения первого порядка определяется из соотношения:
Тогда приращение Δ y_i , может быть найдено путем интегрирования:
Вычислим теперь интеграл по методу прямоугольников:
Из полученного выражения видно, что вычисление интеграла по методу прямоугольников приводит к формуле Эйлера.
Вычислим интеграл по формуле трапеций:
Из выражения видно, что оно совпадает с расчетной формулой усовершенствованного метода Эйлера-Коши. Для получения более точного решения дифференциального уравнения следует воспользоваться более точными методами вычисления интеграла. В методе Рунге-Кутта искомый интеграл представляется в виде следующей конечной суммы:
где Pi— некоторые числа, зависящие от q; K_i (h) — функции, зависящие от вида подынтегральной функции f(x,y) и шага интегрирования h, вычисляемые по следующим формулам:
Значения p, α, β получают из соображений высокой точности вычислений. Формулы Рунге-Кутта третьего порядка (q= 3) имеют следующий вид:
Наиболее часто используется метод Рунге-Кутта четвертого порядка, для которого расчетные формулы имеют следующий вид:
Формулы Рунге-Кутта имеют погрешности порядка h q +1 . Погрешность метода Рунге-Кутта четвертого порядка имеет порядок h 5
4.2 Описание программы ” РЕШЕНИЕ ОДУ “
Программа ”Решение ОДУ“ достаточно проста в использовании.
При запуске программы открывается главное окно программы (рис. 4 ), с установленными по умолчанию начальными условиями в полях ввода.
Назначение элементов ввода данных.
1. Поля X 1, X 2, Y ( x 1), H предназначены для ввода начального и конечного значений отрезка, на котором ищется решение дифференциального уравнения, значения функции при аргументе равном Х1 и величины шага дифференцирования;
2. В поле dY выводится формула дифференциального уравнения 1-й степени, выбранная для решения;
3. В поле dY ( x 1, y 1) выводится значение производной в исходной точке.
Назначение элементов управления и контроля.
1. При нажатии кнопки EXAMPLE активируются “радиокнопки” выбора уравнений;
2. Щелчком “радиокнопки” выбирается соответствующее ей уравнение, вид формулы контролируется по её отображению в поле dY ;
3. Щелчком по кнопке ВЫЧИСЛИТЬ находятся приближенные решения выбранного дифференциального уравнения на заданном интервале;
4. Решения дифференциального уравнения в виде пар значений X — Y выводятся в поля X и Y ; (рис. 5.)
По окончании вычислений активируются кнопка ГРАФИК и пункт меню ГРАФИК главного окна системы.
4.3 Назначение элементов графического окна программы
Вход в графическое окно осуществляется с помощью кнопок ГРАФИК наглавной формеили пункт меню ГРАФИК (рис. 6).
С помощью кнопки ВЫЧЕРТИТЬ на координатную плоскость выводится график функции – решение дифференциального уравнения на заданном интервале.
4.4 Реакция программы при возникновении ошибок
При вводе пользователем ошибочных данных (отсутствии начальных условий, некорректных значений переменных) программа выдает сообщение об ошибке (рис.7 а, б) рис.7 а. рис.7 б.
4.5 Перечень компонент DELPHI использованных в программе
В Form 1 использованы компоненты:
— Edit1.text, Edit2.text, Edit3.text, Edit4.text – для ввода начальных условий дифференциального
— Memo4.TMemo – для вывода формулы уравнения;

— Memo1.TMemo, Memo2.TMemo — для вывода результатов вычислений;
— Memo3.TMemo – для вывода значения производной в точке (Х0,Y0)
— ScrollBars ssVertical всвойствахMemo1.TMemo, Memo2.TMemo;
— Button1 “Вычислить”, Button2 “Очистить”, Button3 “График”, Button4 “Выход”,
Button5 “Example”, Button6 “UnExample”;
— Label1.TLabel — Label9.TLabel – дляотображенияназначениякомпонентов Memo иEdit;
— RadioGroup – для выбора вида уравнения;
В Form 2 использованы компоненты:
— MainMenu- для построения графика;
В Form 3 использованы компоненты:
— Panel1.TPanel – для размещения информации о программе;
— Label1.TLabel – Label14.TLabel – для отображения информации о программе;
— Button1.TButton “OK” – для выхода из окна
5. ТЕХНИЧЕСКИЕ ХАРАКТЕРИСТИКИ И ТРЕБОВАНИЯ К ПО
Программа работает в среде операционных систем Windows 9х, NT.
Требования к ПО
Минимальные системные требования
aпроцессор Intel 486 с рабочей частотой 66 MHz и выше;
b) операционная система Windows 95, 98, NT 4.0, 2000, XP;
с) 16 Мбайт оперативной памяти (или более);
d) 3 Мбайт свободного пространства на жёстком диске.
6. ТЕКСТ ПРОГРАММЫ
Windows, Messages, SysUtils, Classes, Graphics, Controls, Forms, Dialogs,
StdCtrls, CheckLst, ComCtrls, ExtCtrls,math, Menus;
Метод рунге кутта решения диф уравнений и их систем
Пример 1. Найти приближённое решение задачи Коши методом Рунге-Кутты 4 порядка на заданном отрезке с шагом h = 0,1.
Решение:
Для начала, найдем точное решение этого линейного уравнения первого порядка
Тогда точное решение имеет вид :
Найдем приблизительное численное решение дифференциального уравнения методом Рунге-Кутты ( Рунге-Кутта ) 4-го порядка.
Формулы для метода Рунге-Кутты:

Название: Программа для решения дифференциальных уравнений первого порядка методом Рунге-Кутта
Раздел: Рефераты по информатике, программированию
Тип: курсовая работа Добавлен 21:04:57 05 апреля 2011 Похожие работы
Просмотров: 3135 Комментариев: 14 Оценило: 4 человек Средний балл: 5 Оценка: неизвестно Скачать

n Точное значение

0,0025 0,002501 0,010018 0,003337

1 0,1 0,000334 0,01002 0,022632 0,022675 0,040633 0,023545 0,000334

2 0,2 0,002688 0,040644 0,064369 0,064592 0,094933 0,065583 0,002688

3 0,3 0,009246 0,09496 0,133009 0,133709 0,181431 0,134971 0,009246

4 0,4 0,022744 0,181492 0,241149 0,242961 0,318567 0,244713 0,022743

5 0,5 0,047215 0,318698 0,414566 0,418916 0,543032 0,421449 0,047215

6 0,6 0,08936 0,543305 0,703721 0,713888 0,926793 0,717553 0,089359

7 0,7 0,161115 0,927332 1,207835 1,231503 1,61371 1,23662 0,161115

8 0,8 0,284777 1,614692 2,127454 2,183025 2,901619 2,189545 0,284779

9 0,9 0,503732 2,903203 3,884036 4,016816 5,434521 4,023238 0,503741

10 1 0,906055 0,906094

Пример 2. Найти приближённое решение задачи Коши методом Рунге-Кутты 4 порядка на заданном отрезке с шагом h = 0,01.
Найдем точное решение этого уравнения:

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

Найдем численное приблизительное решение дифференциального уравнения методом Рунге-Кутты 4-го порядка.
Сведем заменой переменных это уравнение 3 – го порядка к системе диф. уравнений

И начальные условия примут вид:
Для системы вида:

Напомним формулы для этого метода

Соответственно, здесь будет:

Опять вспомним систему:
, откуда

Внимание, в таблице даны не все шаги при h=0,01, обращайте внимание на n.

Программа для решения дифференциальных уравнений первого порядка методом Рунге-Кутта (стр. 4 из 6)

Из полученного выражения видно, что вычисление интеграла по методу прямоугольников приводит к формуле Эйлера.

Вычислим интеграл по формуле трапеций:

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

где Pi— некоторые числа, зависящие от q; K_i (h) — функции, зависящие от вида подынтегральной функции f(x,y) и шага интегрирования h, вычисляемые по следующим формулам:

Значения p, α, β получают из соображений высокой точности вычислений. Формулы Рунге-Кутта третьего порядка (q= 3) имеют следующий вид:

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

Формулы Рунге-Кутта имеют погрешности порядка h q +1 . Погрешность метода Рунге-Кутта четвертого порядка имеет порядок h 5

4.2 Описание программы ” РЕШЕНИЕ ОДУ “

Программа ”Решение ОДУ“ достаточно проста в использовании.

При запуске программы открывается главное окно программы (рис. 4 ), с установленными по умолчанию начальными условиями в полях ввода.

Назначение элементов ввода данных.

1. Поля X 1, X 2, Y ( x 1), H предназначены для ввода начального и конечного значений отрезка, на котором ищется решение дифференциального уравнения, значения функции при аргументе равном Х1 и величины шага дифференцирования;

2. В поле dY выводится формула дифференциального уравнения 1-й степени, выбранная для решения;

3. В поле dY ( x 1, y 1) выводится значение производной в исходной точке.

Назначение элементов управления и контроля.

1. При нажатии кнопки EXAMPLE активируются “радиокнопки” выбора уравнений;

2. Щелчком “радиокнопки” выбирается соответствующее ей уравнение, вид формулы контролируется по её отображению в поле dY ;

3. Щелчком по кнопке ВЫЧИСЛИТЬ находятся приближенные решения выбранного дифференциального уравнения на заданном интервале;

4. Решения дифференциального уравнения в виде пар значений X — Y выводятся в поля X и Y ; (рис. 5.)

По окончании вычислений активируются кнопка ГРАФИК и пункт меню ГРАФИК главного окна системы.

4.3 Назначение элементов графического окна программы

Вход в графическое окно осуществляется с помощью кнопок ГРАФИК наглавной формеили пункт меню ГРАФИК (рис. 6).

С помощью кнопки ВЫЧЕРТИТЬ на координатную плоскость выводится график функции – решение дифференциального уравнения на заданном интервале.

4.4 Реакция программы при возникновении ошибок

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

4.5 Перечень компонент DELPHI использованных в программе

В Form 1 использованы компоненты:

— Edit1.text, Edit2.text, Edit3.text, Edit4.text – для ввода начальных условий дифференциального

— Memo4.TMemo – для вывода формулы уравнения;

— Memo1.TMemo, Memo2.TMemo — для вывода результатов вычислений;

— Memo3.TMemo – для вывода значения производной в точке (Х0,Y0)

— ScrollBars ssVertical всвойствахMemo1.TMemo, Memo2.TMemo;

— Button1 “Вычислить”, Button2 “Очистить”, Button3 “График”, Button4 “Выход”,

Button5 “Example”, Button6 “UnExample”;

— Label1.TLabel — Label9.TLabel – дляотображенияназначениякомпонентов Memo иEdit;

— RadioGroup – для выбора вида уравнения;

В Form 2 использованы компоненты:

— MainMenu- для построения графика;

В Form 3 использованы компоненты:

— Panel1.TPanel – для размещения информации о программе;

— Label1.TLabel – Label14.TLabel – для отображения информации о программе;

— Button1.TButton “OK” – для выхода из окна

5. ТЕХНИЧЕСКИЕ ХАРАКТЕРИСТИКИ И ТРЕБОВАНИЯ К ПО

Программа работает в среде операционных систем Windows 9х, NT.

Требования к ПО

Минимальные системные требования

aпроцессор Intel 486 с рабочей частотой 66 MHz и выше;

b) операционная система Windows 95, 98, NT 4.0, 2000, XP;

с) 16 Мбайт оперативной памяти (или более);

d) 3 Мбайт свободного пространства на жёстком диске.

6. ТЕКСТ ПРОГРАММЫ

Windows, Messages, SysUtils, Classes, Graphics, Controls, Forms, Dialogs,

StdCtrls, CheckLst, ComCtrls, ExtCtrls,math, Menus;

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

1. Постановка задачи

3. Выбор метода реализации программы

6. Идентификация переменных

8. Обсуждение результатов

9. Инструкция к программе

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

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

j ( x, y, y 1 , . y (n) )=0. 1.1

Из теории ОДУ известно, что уравнение (1.1) эквивалентно системе n уравнений первого порядка

j k (x, y 1, y 1 ’ ,y 2 ,y 2 ’ , . ,y n ,y n ’ )=0. 1.2

Уравнение (1.1) и эквивалентная ему система (1.2) имеют бесконечное множество решений. Единственные решения выделяют с помощью дополнительных условий, которым должны удовлетворять искомые решения. В зависимости от вида таких условий рассматривают три типа задач, для которых доказано существование и единственность решений.

Первый тип – это задачи Коши, или задачи с начальными условиями. Для таких задач кроме исходного уравнения (1.1) в некоторой точке xo должны быть заданы начальные условия, т.е. значения функции y(x) и ее производных

y(x 0 )=y 0 ’ , y ’ (x 0 )=y 10 , . , y (n-1) (x 0 )=y n-1 , 0 .

Для системы ОДУ типа (1.2) начальные условия задаются в виде

y 1 (x 0 )=y 10 , y 2 (x 0 )=y 20 , . , y n (x 0 )=y n0 . 1.3

Ко второму типу задач относятся так называемые граничные, или краевые задачи, в которых дополнительные условия задаются в виде функциональных соотношений между искомыми решениями. Количество условий должно совпадать с порядком n уравнения или системы. Если решение задачи определяется в интервале x є [ x 0 ,x k] , то такие условия могут быть заданы как на границах, так и внутри интервала. Минимальный порядок ОДУ, для которых может быть сформулирована граничная задача, равен двум.

Третий тип задач для ОДУ – это задачи на собственные значения. Такие задачи отличаются тем, что кроме искомых функций y(x) и их производных в уравнения входят дополнительно m неизвестных параметров l 1, l 2, ј , хm, которые называются собственными значениями. Для единственности решения на интервале [x0, xk] необходимо задать m+n граничных условий. В качестве примера можно назвать задачи определения собственных частот, коэффициентов диссипации, структуры электромагнитных полей и механических напряжений в колебательных системах, задачи нахождения фазовых коэффициентов, коэффициентов затухания, распределения напряженностей полей волновых процессов и т. д.

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

Большинство методов решения ОДУ основано на задаче Коши, алгоритмы и программы для которой рассматриваются в дальнейшем.

1. Постановка задачи

Многие процессы химической технологии описываются СДУ — начиная от кинетических исследований и заканчивая химическими технологическими процессами. В основу математических способов описания процессов положены СДУ и СЛАУ. Эти уравнения описывают материальные и тепловые балансы объектов химической технологии, а так же структуры потоков технических веществ в этих аппаратах.

Для получения, распределения технологических параметров во времени и в пространстве (в пределах объекта), необходимо произвести СДУ методом, которых дал бы высокую точность решения при минималььных затратах времени на решение, потому что ЭВМ должна работать в режиме реального времени и успевать за ходом технологического процесса. Если время на решение задачи большое, то управляющее воздействие, выработанное на ЭВМ может привести к отрицательным воздействиям. Методов решения существует очень много. В данной работе будет рассмотрен метод решения СДУ методом Рунге-Кутта 4 порядка.

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

Так как коэффициенты K1,K2,K3,K4 являются константами, то можно уравнение записать в следущем виде.

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

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

Методы Рунге-Кутта обладают следующими свойствами:

1. Эти методы являются одноступенчатыми: чтобы найти у m+1, нужна

информация о предыдущей точке xm, ym.

2. Они согласуются с рядом Тейлора вплоть до членов порядка hp, где степень р

различна для различных методов и называется порядковым номером или

3. Они не требуют вычисления производных от f (x, y), а требуют вычисления

Рассмотрим сначала геометрическое построение и выведем некоторые формулы на основе геометрических аналогий. После этого мы подтвердим полученные результаты аналитически.

Предположим, нам известна точка xm, ym на искомой кривой. Тогда мы можем провести прямую линию с тангенсом угла наклона уў m =f(x m, y m ), которая пройдет через точку x m, y m. Это построение показано на рис. 1, где кривая представляет собой точное, но конечно неизвестное решение уравнения, а прямая линия L1 построена так, как это только что описано.

Тогда следующей точкой решения можно считать ту, где прямая L 1 пересечет ординату, проведенную через точку x=x m+1 =x m +h.

Уравнение прямой L 1 выглядит так: y=y m +yў m (x-x m ) так как yў =f(x m, y m ) и кроме того, x m+1 =x m +h тогда уравнение примет вид

y m+1 =y m +h*f(x m, y m ) 1. 1

Ошибка при x=x m+1 показана в виде отрезка е. Очевидно, найденное таким образом приближенное значение согласуется с разложением в ряд Тейлора вплоть до членов порядка h, так что ошибка ограничения равна e t =Кh2

Заметим, что хотя точка на графике 1 была показана на кривой, в действительности y m является приближенным значением и не лежит точно на кривой.

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

Рассмотрим исправленный метод Эйлера и модификационный метод Эйлера. В исправленном методе Эйлера мы находим средний тангенс угла наклона касательной для двух точек: x m, y m и x m +h, y m +hyў m. Последняя точка есть та самая, которая в методе Эйлера обозначалась x m+1, y m+1. Геометрический процесс нахождения точки x m+1, y m+1 можно проследить по рис. 2. С помощью метода Эйлера находится точка x m +h, y m +hyў m, лежащая на прямой L 1. В этой точке снова вычисляется тангенс, дает прямую L . Наконец, через точку x m, y m мы проводим прямую L, параллельную L . Точка, в которой прямая L пересечется с ординатой, восстановленной из x=x m+1 =x m +h, и будет искомой точкой x m+1, y m+1.

Тангенс угла наклона прямой L и прямой L равен

Ф(x m, y m, h)=Ѕ [f(x m, y m )+f(x m +h, y m +yў m h)] 1. 2

где yў m =f(x m, y m ) 1. 3

Уравнение линии L при этом записывается в виде

y=y m +(x-x m )Ф(x m, y m, h),

y m+1 =y m +hФ(x m, y m, h). 1. 4

Соотношения 1. 2, 1. 3, 1. 4 описывают исправленный метод Эйлера.

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

f(x, y)=f(x m, y m )+(x-x m )¶ f/¶ x+(y-y m )¶ f/¶ x+ј 1. 5

где частные производные вычисляются при x=x m и y=y m.

Подставляя в формулу 1. 5 x=x m +h и y=y m +hyў m и используя выражение 1. 3 для yў m, получаем

f(x m +h, y m +hyў m )=f+hf x +hff y +O(h2),

где снова функция f и ее производные вычисляются в точке x m, y m. Подставляя результат в 1. 2 и производя необходимые преобразования, получаем

Ф(x m, y m, h)=f+h/2(f x +ff y )+O(h2).

Подставим полученное выражение в 1. 4 и сравним с рядом Тейлора

y m+1 =y m +hf+h2/2(f x +ff y )+O(h3).

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

Рассмотрим модификационный метод Эйлера. Рассмотрим рис. 3 где первоначальное построение сделано так же, как и на рис. 2. Но на этот раз мы берем точку, лежащую на пересечении этой прямой и ординатой x=x+h/2. На рисунке эта точка образована через Р, а ее ордината равна y=y m +(h/2)yў m. Вычислим тангенс угла наклона касательной в этой точке

Ф(x m, y m, h)=f+(x m +h/2, y m +h/2*yў m ), 1. 6

где yў m =f(x m, y m ) 1. 7

Прямая с таким наклоном, проходящая через Р, обозначена через L *. Вслед за тем, мы проводим через точку xm, ym прямую параллельную L *, и обозначаем ее через L 0. Пересечение этой прямой с ординатой x=x m +h и даст искомую точку x m+1, y m+1. Уравнение прямой можно записать в виде y=y m +(x-x m )Ф(x m, y m, h),

где Ф задается формулой 1. 6. Поэтому

y m+1 =y m +hФ(x m, y m, h) 1. 8

Соотношения 1. 6, 1. 7, 1. 8 описывают так называемый модификационный метод Эйлера и является еще одним методом Рунге-Кутта второго порядка. Обобщим оба метода. Заметим, что оба метода описываются формулами вида

y m+1 =y m +hФ(x m, y m, h) 1. 9

и в обоих случаях Ф имеет вид

Ф(x m, y m, h)=a 1 f(x m, y m )+a 2 f(x m +b 1 h, y m +b 2 hyў m ), 1. 10

где yў m =f(x m, y m ) 1. 11

В частности, для исправленного метода Эйлера

В то время как для модификационного метода Эйлера

Формулы 1. 9, 1. 10, 1. 11 описывают некоторый метод типа Рунге-Кутты. Посмотрим, какого порядка метод можно рассчитывать получить в лучшем случае и каковы допустимые значения параметров a 1, a 2, b 1 и b 2 .

Чтобы получить соответствие ряду Тейлора вплоть до членов степени h, в общем случае достаточно одного параметра. Чтобы получить согласование вплоть до членов степени h2, потребуется еще два параметра, так как необходимо учитывать члены h2f x и h2ff y. Так как у нас имеется всего четыре параметра, три из которых потребуются для создания согласования с рядом Тейлора вплоть до членов порядка h2, то самое лучшее, на что здесь можно рассчитывать — это метод второго порядка.

В разложении f(x, y) в ряд 1. 5 в окрестности точки x m, y m положим x=x m +b 1 h,

Тогда f(x m +b 1 h, y m +b 2 hf)=f+b 1 hf x +b 2 hff y +O(h2), где функция и производные в правой части равенства вычислены в точке x m, y m.

Тогда 1. 9 можно переписать в виде y m+1 =y m +h[a 1 f+a 2 f+h(a 2 b 1 f x +a 2 b 2 ff y )]+O(h3).

Сравнив эту формулу с разложением в ряд Тейлора, можно переписать в виде

y m+1 =y m +h[a 1 f+a 2 f+h(a 2 b 1 f x +a 2 b 2 ff y )]+O(h3).

Если потребовать совпадения членов hf, то a 1 +a 2 =1.

Сравнивая члены, содержащие h2f x, получаем a 2 b 1 =1/2.

Сравнивая члены, содержащие h2ff y, получаем a 2 b 2 =1/2.

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

Положим, например, a 2 =w № 0. тогда a 1 =1-w , b 1 =b 2 =1/2w и соотношения 1. 9, 1. 10, 1. 11 сведутся к

y m+1 =y m +h[(1-w )f(x m, y m )+w f(x m +h/2w , y m +h/2w f(x m, y m ))]+O(h3) 1. 12

Это наиболее общая форма записи метода Рунге-Кутта второго порядка. При w =1/2 мы получаем исправленный метод Эйлера, при w =1 получаем модификационный метод Эйлера. Для всех w , отличных от нуля, ошибка ограничения равна

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

y m+1 =y m +h/6(R 1 +2R 2 +2R 3 +R 4 ) 1. 14

где R 1 =f(x m, y m ), 1. 15

R 2 =f(x m +h/2, y m +hR 1 /2), 1. 16

R 3 =f(x m +h/2, y m +hR 2 /2), 1. 17

R 4 =f(x m +h/2, y m +hR 3 /2). 1. 18

Ошибка ограничения для этого метода равна e t =kh5

так что формулы 1. 14-1. 18 описывают метод четвертого порядка. Заметим, что при использовании этого метода функцию необходимо вычислять четыре раза.

3. Выбор метода реализации программы

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

этот метод является одноступенчатым и одношаговым

требует информацию только об одной точке

имеет небольшую погрешность

значение функции рассчитывается при каждом шаге

Решение систем дифференциальных уравнений методом Рунге-Кутты

Имеется система диффуров

где k1=50.93,k2= 1569.12
Необходимо решить эту систему уравнений методом Рунге-Кутты 4го порядка
Имеются наброски но не совсем рабочие(
06.06.2020, 22:33
Решение дифференциального уравнения методом Рунге-Кутты
Не могу понять почему не правильно выдает значения. #include «stdafx.h» #include .

Решение задачи Коши методом Эйлера и Рунге-Кутты
Доброго времени суток всем:) Писал прогу для решения задачи Коши методом Эйлера и Рунге-Кутты на.

Численное интегрирование системы дифференциальных уравнений методом Рунге — Кутта
Доброго времени суток. Пытаюсь разработать функцию для численного интегрирования систем.

Решение диф. уравнения для колебаний методом Рунге—Кутты 4го порядка.(С++)
Ребят помогите пожалуйста решить такое на C++ уравнение x»=(-w^2)*x где w-const там замену.

Решение систем дифференциальных уравнений методами Эйлера
Здравствуйте, дорогие программисты и просто любители языка С++! В принципе, моя будущая.

Илон Маск рекомендует: Работа с датой и временем в PHP

n	Точное значение
				0,0025	0,002501	0,010018	0,003337
1	0,1	0,000334	0,01002	0,022632	0,022675	0,040633	0,023545	0,000334
2	0,2	0,002688	0,040644	0,064369	0,064592	0,094933	0,065583	0,002688
3	0,3	0,009246	0,09496	0,133009	0,133709	0,181431	0,134971	0,009246
4	0,4	0,022744	0,181492	0,241149	0,242961	0,318567	0,244713	0,022743
5	0,5	0,047215	0,318698	0,414566	0,418916	0,543032	0,421449	0,047215
6	0,6	0,08936	0,543305	0,703721	0,713888	0,926793	0,717553	0,089359
7	0,7	0,161115	0,927332	1,207835	1,231503	1,61371	1,23662	0,161115
8	0,8	0,284777	1,614692	2,127454	2,183025	2,901619	2,189545	0,284779
9	0,9	0,503732	2,903203	3,884036	4,016816	5,434521	4,023238	0,503741
10	1	0,906055	0,906094

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

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

Отправить свою хорошую работу в базу знаний просто. Используйте форму, расположенную ниже

Введение

Цель курсовой работы: изучить метод Рунге — Кутта 4-го порядка для решения обыкновенных дифференциальных уравнений.

Постановка задачи: необходимо составить программу, позволяющую решать обыкновенные дифференциальные уравнения методом Рунге — Кутта 4-го порядка.

Курсовая работа состоит из 3 разделов, содержит 6 рисунков, 3 листинга, 1 приложение и 18 страниц.

1. Теоретическая часть

= f(x, y) (1)

Разобьём отрезок [Х, XN ] на N частей так, что xn+1 — хn = hn ,

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

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

Численные методы решения задачи Коши

Явный метод Эйлера

Программная реализация явного метода Эйлера

Неявный метод Эйлера

Программная реализация неявного метода Эйлера

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

Многошаговые методы

Жесткие системы ОДУ

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

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

Метод рунге кутта решения диф уравнений и их систем

Программа для решения дифференциальных уравнений первого порядка методом Рунге-Кутта (стр. 4 из 6)

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

Решение систем дифференциальных уравнений методом Рунге-Кутты

Разобьём отрезок [Х, X_N ] на N частей так, что x_n₊₁ — х_n = h_n _,