Nickolay.info. Алгоритмы. Решение задачи Коши методом Эйлера
1. Типы задач для обыкновенных дифференциальных уравнений
Обыкновенные дифференциальные уравнения (ОДУ) широко используются для математического моделирования процессов и явлений в различных областях науки и техники. Переходные процессы в радиотехнических цепях, кинетика химических реакций, динамика биологических популяций, движение космических объектов, модели экономического развития исследуются с помощью ОДУ.
В дифференциальное уравнение n-го порядка в качестве неизвестных величин входят функция у(х) и её первые n производных по аргументу х.
(1)
Из теории ОДУ известно, что уравнение (1) эквивалентно системе n уравнений первого порядка
(2)
Уравнение (1) и эквивалентная ему система (2) имеют бесконечное множество решений. Единственные решения выделяют с помощью дополнительных условий, которым должны удовлетворять искомые решения. В зависимости от вида таких условий рассматривают три типа задач, для которых доказано существование и единственность решений.
Первый тип — это задачи Коши, или задачи с начальными условиями. Для таких задач кроме исходного уравнения (1) в некоторой точке х0 должны быть заданы начальные условия, т.е. значения функции у(х) и ее производных
Для системы ОДУ типа (2) начальные условия задаются в виде
(3)
Ко второму типу задач относятся так называемые граничные, или краевые задачи, в которых дополнительные условия задаются в виде функциональных соотношений между искомыми решениями. Количество условий должно совпадать с порядком n уравнения или системы. Если решение задачи определяется в интервале x принадлежит [х0, хk], то такие условия могут быть заданы как на границах, так и внутри интервала. Минимальный порядок ОДУ, для которых может быть сформулирована граничная задача, равен двум.
Третий тип задач для ОДУ — это задачи на собственные значения. Такие задачи отличаются тем, что кроме искомых функций у(х) и их производных в уравнения входят дополнительно m неизвестных параметров λ1, λ2, …, λm, которые называются собственными значениями. Для единственности решения на интервале [х0, хk] необходимо задать n+m граничных условий. В качестве примера можно назвать задачи определения собственных частот, коэффициентов диссипации, структуры электромагнитных полей и механических напряжений в колебательных системах, задачи нахождения фазовых коэффициентов, коэффициентов затухания, распределения напряженностей полей волновых процессов и т. д.
К численному решению ОДУ приходится обращаться, когда не удается построить аналитическое решение задачи через известные функции. Хотя для некоторых задач численные методы оказываются более эффективными даже при наличии аналитических решений.
Большинство методов решения ОДУ основано на задаче Коши, алгоритм решения которой рассматривается далее.
Систему ОДУ (2) часто удается представить в каноническом виде, в так называемой форме Коши
(4)
При формулировке задачи Коши система (4) дополняется начальными условиями (3). Для простоты рассмотрим задачу Коши для одного уравнения типа (4), а затем полученные алгоритмы обобщим на систему n уравнений
(5)
В окрестности точки х0 функцию у(х) разложим в ряд Тейлора
(6)
который можно применить для приближенного определения искомой функции у(х). В точке х0 + h при малых значениях h можно ограничиться двумя членами ряда (6), тогда
(7)
где O(h 2 ) — бесконечно малая величина порядка h 2 . Заменим производную у'(x0), входящую в формулу (7), на правую часть уравнения (5):
(8)
Теперь приближенное решение в точке х1 = х0 + h можно вновь рассматривать как начальное условие и по формуле (8) найти значение искомой функции в следующей точке х2 = x1 + h1. В результате получен простейший алгоритм решения задачи Коши, который называется методом Эйлера, или методом ломаных. Последнее название связано с геометрической интерпретацией процесса (см. рис.); искомую функцию у(х) мы заменяем ломаной линией, представляющей собой отрезки касательных к этой функции в узлах
Рис. Метод Эйлера
Формула (8) может быть получена из других соображений. Заменим производную в левой части уравнения (5) приближенным конечно-разностным отношением
Нетрудно видеть эквивалентность последнего выражения с алгоритмом Эйлера (8).
На каждом шаге метода Эйлера решение у(х) определяется с погрешностью за счет отбрасывания членов ряда Тейлора, пропорциональных h в степени выше первой. Это означает, что метод Эйлера имеет второй порядок локальной погрешности. Доказано, что глобальная погрешность метода имеет первый порядок; и при постоянном шаге h для оценки погрешности применима первая формула Рунге
, (9)
где yh(x) — приближенное решение дифференциального уравнения в точке х, полученное с шагом h; уkh(х) — приближенное решение того же уравнения с шагом kh; р — порядок метода.
Формула (9) позволяет опытным путем определить шаг h, обеспечивающий требуемую точность решения у(х). Так же, как и при вычислении определенных интегралов, можно осуществлять автоматическое изменение шага в процессе интегрирования дифференциального уравнения.
Для уточнения решения применима вторая формула Рунге
(10)
Формула Эйлера (8) обобщается для систем ОДУ, записанных в форме Коши (4) с начальными условиями (3)
(11)
Представленная ниже программа реализует метод Эйлера решения задачи Коши для системы дифференциальных уравнений. Функция уравнения задаётся подпрограммой f(x), точное решение – подпрограммой ft(x). Пользователь вводит интервал поиска решения [A,B], число шагов N, начальное значение Y(0). Программа выводит найденное решение и оценивает его максимальную погрешность.
В паскале решение дифференциального уравнения
Возможно, ваша ваша проблема уже имеет решение на нашем форуме
04.05.2011, 13:05
#2 (permalink)
Извините — ничего не понял! Уравнение вполне себе легко решается аналитически, и выглядит это так:
1. Записываем в виде y’ — y = 1 + x. 2. Решаем однородное уравнение y’ — y = 0. Его решением, очевидно, является C*exp(x). 3. Частное решение ищем методом вариации постоянной, т.е. в виде C(x)*exp(x). Подставляя в исходное уравнение и интегрируя, находим C(x)=C1-exp(-x)*(x+2). 4. Таким образом, окончательно y=C1*exp(x)-x-2.
И при чем тут Паскаль, программирование, точность и всё такое прочее? Кроме того, полагаю, что Ваше предположение о том, что у=1 имеет какое-то отношение к таинственному «шагу», есть полный абсурд. Скорее всего, у=1 — это значение искомой функции при некотором значении аргумента х (например, при х=0), которое Вы пропустили, и которое требуется для определения константы С1.
Метод Рунге-Кутта решения диф. уравнений и их систем.
Метод позволяет решать системы обыкновенных дифференциальных уравнений (ОДУ) первого порядка следующего вида:
которые имеют решение:
где t — независимая переменная (например, время); X, Y и т.д. — искомые функции (зависимые от t переменные). Функции f, g и т.д. — заданы. Также предполагаются заданными и начальные условия, т.е. значения искомых функций в начальный момент.
Одно диф. уравнение — частный случай системы с одним элементом. Поэтому, далее речь пойдет для определенности о системе уравнений.
Метод может быть полезен и для решения диф. уравнений высшего (второго и т.д.) порядка, т.к. они могут быть представлены системой диф. уравнений первого порядка.
Метод Рунге-Кутта заключается в рекурентном применении следующих формул:
Реализация Метода Рунге-Кутта на Delphi может выглядеть так (привожу полностью модуль):
Модуль полностью работоспособен. Возвращаемое функцией Runge_Kutt значение — код ошибки. Вы можете дополнить список ошибок по своему усмотрению. Рассчитанные функции системы помещаются в массив Res. Чтобы не загромождать код, в модуле опущены проверки (типа блоков try). Рекомендую их добавить по своему усмотрению.
Ниже приводится описание функции Runge_Kutt и типов, использующихся в модуле.
FunArray — вектор функций (правых частей уравнений системы);
First, Last — начальная и конечная точки расчетного интервала;
Steps — число шагов по расчетному интервалу;
InitArray — вектор начальных значений
var Res — матрица результатов включая независимую переменную.
В модуле описаны типы:
Функция возвращает коды ошибок:
0 — нет ошибок;
100 — число уравнений не равно числу начальных условий.
Решение содержится в переменной-матрице Res. Первый индекс матрицы относится к переменной (0 — независимая переменная, 1 — первая зависимая и т.д.), второй — к номеру расчетной точки (0 — начальная точка).
Рассмотрим один пример использования модуля. Создадим новое приложение и подключим к нему модуль. На форме приложения разместим кнопку Button1 и область текста Memo1. Поместим в приложение две функции и обработчик нажатия кнопки:
Нажатие кнопки приведет к расчету точек системы, которые будут выведены в текстовую область.
Модуль с примером и справкой можно скачать бесплатно по адресу RK.zip (ZIP, 15,3Kb) (русский вариант). Английский вариант (условно-бесплатный) можно скачать по адресу RK_Eng.zip (ZIP, 23.4Kb)
Скачала по Вашей ссылке русский вариант, изменила для своей системы диф. уравнений, но при запуске выдаёт ошибку : Project Ex.exe raised exception class EOverflow with message ‘ Floating point overflow ‘ Помогите, пожалуйста .
Вот изменённый мною модуль:
unit Unit1; interface uses SysUtils, Forms, StdCtrls, Controls, Classes, Dialogs, Math; type TForm1 = class(TForm) Memo1: TMemo; rk_But: TButton; procedure rk_ButClick(Sender: TObject); private < Private declarations > public < Public declarations > end; var Form1: TForm1; pn,k,ro,Pzv: Extended;
(1)
(2)
(3)
(4)
(5)
(6)
(7)
(8)
Рис. Метод Эйлера
, (9)
(10)
(11)
