Понятие о численном решении дифференциальных уравнений метод эйлера

Дифференциальным уравнением первого порядка называется уравнение вида F(x,y,у’)=0 или у’=f(x,y). Функция y(x), при подстановке которой уравнение обращается в тождество, называется решением дифференциального уравнения.

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

Рассмотрим два варианта вывода расчетных формул

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

Описание: Недостатки метода Эйлера 4. Идея метода Эйлера очень проста. В результате приходим к приближённому уравнению: Поскольку по определению у= окончательно имеем следующее уравнение являющееся основой метода Эйлера: 8 Конечно это уравнение является лишь приближённым и мы надеемся что чем меньше величина шага h тем оно будет более точным уменьшается локальная погрешность метода то есть погрешность на одном его шаге.

Дата добавления: 2015-08-14

Размер файла: 153.35 KB

Работу скачали: 23 чел.

Поделитесь работой в социальных сетях

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

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

Содержание

1. Метод Эйлера

2. О решении ОДУ высших степеней и их систем

3. Недостатки метода Эйлера

4. Четырёхточечный метод Рунге-Кутты

5. Вычислительный эксперимент

Литература

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

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

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

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

Уравнение (1) является уравнением первого порядка и поэтому его общее решение зависит лишь от одной произвольной постоянной. В связи с этим для выделения частного решения достаточно задать лишь одно начальное условие:

Это условие означает, что при фиксированном значении аргумента искомая функция у(х) должна иметь некоторое известное значение .

Таким образом, перед нами стоит вопрос о решении простейшей задачи Коши, которая определяется заданием дифференциального уравнения и некоторого начального условия (более подробно смотри далее):

Приведённую задачу Коши мы собираемся решать численно. Что это означает? Всем хорошо известны из школьного курса физики «Четырёхзначные математические таблицы» Брадиса. В них разные функции, в частности тригонометрические, задаются в табличной форме: в первой колонке указаны дискретные значения аргумента, например =0.1, 0.2, 0.3, 0.4,…, а в соседней колонке – соответствующие им значения табулируемой функции, например, синуса: sin (0.1), sin (0.2), sin (0.3), sin (0.4),….

Аналогичным образом численное решение рассматриваемой нами задачи Коши будет представлено в форме таблицы значений аргумента

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

Шаг табулирования h предполагается достаточно малым. Величина его, очевидно, зависит от решаемой нами задачи. Например, при описании движения траектории Земли, возможно, вполне достаточно будет в качестве временного шага выбрать один день (55 точек на одном обороте Земли вокруг Солнца), но вряд ли кому-нибудь потребуется при решении этой задачи выбрать в качестве величины этого шага одну микросекунду.

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

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

Первая строка этой таблицы содержит известные нам величины, входящие в начальное условие (3 b ) – и . Во второй строке – аргумент =+ h известен по построению, а (то есть значение функции у()) является неизвестным.

Идея метода Эйлера очень проста. По определению производной,

В случае производной от функции у(х) в точке х=имеем

поскольку в нашем случае . Если шаг h достаточно мал , то приближённое значение производной () можно найти, опуская предел в формуле (5), т.е. полагая

С другой стороны, уравнение (3а), по определению, справедливо при любом значении аргумента х и , стало быть, оно верно и в точках :

Подставим теперь в это точное равенство приближённое выражение для первой производной (6). В результате приходим к приближённому уравнению:

Поскольку, по определению у()=, окончательно имеем следующее уравнение, являющееся основой метода Эйлера:

Конечно, это уравнение является лишь приближённым, и мы надеемся, что, чем меньше величина шага h , тем оно будет более точным (уменьшается локальная погрешность метода, то есть погрешность на одном его шаге).

Заметим, что в численном анализе не принято писать знак приближённого равенства, вместо него используется знак точного равенства (по умолчанию предполагается, что все формулы численного анализа являются приближёнными).

Полагая k =0 из уравнения (8) имеем

Все величины в правой части этого уравнения нам известны (см. начальные условия (3 b )), что даёт возможность вычислить неизвестное значение и заполнить тем самым вторую строку таблицы 1.

Полагая далее k =1, из уравнения (8) имеем

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

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

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

Из рассмотренной выше вычислительной схемы метода Эйлера со всей очевидностью следует, что она одинаково применима для любого дифференциального уравнения первого порядка: как линейного, так и нелинейного. 2

2. О решении ОДУ высших степеней и их систем

Мы разобрали технику применения метода Эйлера для численного решения уравнения первого порядка. С другой стороны, в обсуждавшиеся ранее математические модели входили уравнения второго порядка, что, фактически, связано с применением второго закона Ньютона. В связи с этим возникает вопрос: как решать ОДУ высших порядков и их системы?

Оказывается, что метод Эйлера, как, впрочем, и другие численные методы, пригодные для решения задачи Коши (3), можно достаточно легко обобщить на случай решения систем уравнений первого порядка в канонической форме (см. далее), а уравнения высших степеней сводятся к таким системам уравнений. Сведение дифференциального уравнения высшего порядка к системе уравнений первого порядка достигается с помощью введения дополнительных переменных, равных последовательным производным неизвестной функции, входящей в это уравнение. В случае уравнения второго порядка типа уравнений гармонического осциллятора или математического маятника достаточно ввести только одну новую переменную, причем имеющую ясный физический смысл. Она является скоростью (линейной – в случае уравнения гармонического осциллятора или угловой – в случае математического маятника).

При рассмотрении уравнения математического маятника положим . Тогда и уравнение можно записать в форме следующих двух уравнений первого порядка

Рассматривая эти уравнения для момента времени имеем

Согласно основной идее метода Эйлера заменим точные значения производных, входящих в уравнения (9), их допредельными образами:

где – шаг по времени, с которым решается рассматриваемое ОДУ. Здесь Подстановка соотношений (11) в уравнения (10) приводит нас к следующим уравнениям метода Эйлера

Полагая k =0,1,2,3 и т. д., мы будем получать последовательные значения угла и соответствующие ему значения угловой скорости u ( t ) в последовательные моменты времени

Совершенно аналогичным образом систему двух ОДУ второго порядка, которая описывает движения планет вокруг Солнца (31), можно свести к системе из четырёх уравнений первого порядка, если ввести две дополнительные переменные, представляющие собой скорости движения тела m вдоль координатных осей х и у соответственно:

Тогда исследуемая система принимает вид:

Заменяя в системе (12) каждую производную её допредельным образом, получаем явные формулы для решения этой системы методом Эйлера.

3. Недостатки метода Эйлера

Существует простая геометрическая интерпретация метода Эйлера. Рассмотрим снова задачу Коши (3) для одного ОДУ первой степени и соответствующее ему в методе Эйлера рекуррентное соотношение (8)

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

Пусть решению рассматриваемой нами задачи Коши на рис. 1 отвечает жирная интегральная кривая. Рекуррентное уравнение метода Эйлера (8) можно записать в форме

Поскольку, согласно решаемому нами ОДУ, . Если считать параметр h непрерывной переменной, то легко понять, что уравнение (13) представляет собой уравнение касательной, проведённой в точке к интегральной кривой у(х), проходящей через эту точку.

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

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

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

Может показаться, что если выбрать шаг интегрирования h более маленьким, то можно существенным образом уменьшить эту погрешность. В общем случае это, однако, не так. Действительно, если мы должны получить решение исходной задачи Коши (3) на некотором заданном макроинтервале [,], то уменьшение шага h влечёт за собой увеличение числа шагов интегрирования N , поскольку . «Локальная погрешность» (погрешность на одном шаге) уменьшается при уменьшении h , но увеличение числа шагов может привести к росту «глобальной погрешности» на заданном интервале [,]. Вышеуказанная ситуация характерна для неустойчивых вычислительных процессов, использование которых на практике может привести к катастрофическим последствиям. В силу этого, необходимо рассматривать такие методы численного решения ОДУ, которые порождают достаточно устойчивые численные алгоритмы. 4 Весьма распространёнными и хорошо зарекомендовавшими себя на практике для решения ОДУ являются методы Рунге-Кутты. Это целый класс методов и мы, в качестве примера, рассмотрим так называемый четырёхточечный метод Рунге-Кутты

4. Четырёхточечный метод Рунге-Кутты

Ниже кратко описано применение четырёхточечного метода Рунге-Кутты для решения задачи Коши (11) для дифференциального уравнения первого порядка, разрешённого относительно производной. Таким образом, мы будем рассматривать ту же самую задачу Коши, решение которой ранее рассматривалось методом Эйлера.

Заметим, прежде всего, что решение дифференциального уравнения (3а) фактически определяет зависимость первой производной от двух независимых переменных – х и у. Это очень хорошо видно из рис. 1: фиксируя х, мы имеем бесконечное множество значение производной , поскольку каждая из пересекаемых вертикальной линией х= const интегральных кривых имеет своё направление касательной.

Описываемый метод Рунге-Кутты, как и метод Эйлера, состоит из последовательности шагов величиной h , но, в отличие от последнего метода, на каждом шаге h находится не одно значение производной (в методе Эйлера находилось лишь ), а несколько. Они соответствуют разным значениям аргументов функции f (х,у). В четырёхточечном методе Рунге-Кутты (откуда и происходит его название) находится четыре различных значения и делается некоторое специфическое их усреднение (то есть берётся не простое среднее арифметическое значение этих производных!). После этого делается перемещение из точки в точку по прямой в направлении тангенса угла наклона, которое определяется этим усреднённым значением производной. Заметим, что разные варианты методов Рунге-Кутты отличаются друг от друга стратегией выбора точек, в которых находятся производные на микроинтервале h и формулой усреднения значения этих производных.

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

Этап I . Находим производную к интегральной кривой в точке 1 с координатами , подставляя эти координаты в правую часть дифференциального уравнения, т.е. вычисляя значение . Из точки

1: перемещаемся на полшага вперёд по прямой, направление которой задаётся этим значением производной (то есть по прямой с тангенсом угла наклона к оси абсцисс, равным значению (этот этап метода Рунге-Кутты полностью аналогичен шагу метода Эйлера с шагом h /2). В результате, в плоскости (х,у) мы переходим в точку

Этап II . Через найденную таким образом точку 2 проходит своя интегральная кривая и мы находим направление касательной к ней, то есть вычисляем значение производной

Далее делается полшага вперёд с найденным значением производной (), но снова из начальной точки микроинтервала [,]. Таким образом, мы переходим в плоскости (х,у) в точку

Этап III . В этой точке находим значение производной , подставляя её координаты в правую часть дифференциального уравнения:

Эта производная определяет направление касательной к интегральной кривой, проходящей через точку 3.

Этап IV . Из начальной точки 1: делаем на сей раз полный шаг вперёд (на величину h по оси х) по прямой, в направлении, которое определяется значением производной . В результате мы переходим в точку

Находим производную в этой точке подстановкой её координат в правую часть дифференциального уравнения:

В результате четырёх описанных выше этапов мы нашли четыре значения производных. Производим их усреднение по формуле

Таким образом, является некоторым средневзвешенным значением найдённых четырёх производных: двум «внутренним» значениям производном соответствуют весовые множители 2, а двум крайним – множители 1 [деление в формуле (14) производится на сумму этих четырёх весовых множителей: 6=1+2+2+1].

Далее мы перемещаемся по прямой из начальной точки 1: в направлении, тангенс угла наклона которого к оси абсцисс определяется средним значением производной из формулы (14). Таким образом, из начальной точки с координатами переходим в точку плоскости (х,у) с координатами

Иными словами, полный шаг метода Рунге-Кутты определяется формулами

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

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

где . Здесь есть пятая производная от искомого решения дифференциального уравнения (3а) в некоторой точке на микроинтервале . Таким образом, локальная погрешность метода Рунге-Кутты (то есть погрешность на одном шаге h ) пропорциональна пятой степени шага h и пятой производной искомого решения дифференциального уравнения.

Оценки точности типа (16) позволяют грубо оценить величину шага интегрирования шага h , необходимого для достижения требуемой точности решения исходного дифференциального уравнения. Для метода Эйлера аналогичная погрешность на одном шаге определяется формулой

Таким образом, четырёхточечный метод Рунге-Кутты на три порядка по шагу точнее метода Эйлера (например, при h =0.01 точность метода Рунге-Кутты в миллион раз выше точности метода Эйлера).

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

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

Эти уравнения и представляют собой систему ОДУ в канонической форме.

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

Численное (а по возможности, и аналитическое) решение ОДУ на языке Maple осуществляется с помощью оператора dsolve , с разными спецификациями, которые, в частности, позволяют выбрать необходимый метод численного интегрирования. По умолчанию используется некоторая модификация чеитырёхточечного метода Рунге-Кутты, которая получила название метода Рунге-Кутты-Фельдберга. Она осуществляет решение ОДУ с переменным шагом, величина которого подбирается в зависимости от скорости изменения искомого решения (то есть от крутизны соответствующей интегральной кривой).

5. Вычислительный эксперимент

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

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

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

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

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

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

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

Касьянов В.И.: Руководство к решению задач по высшей математике. — М.: Юрайт, 2011

НИУ БелГУ ; гл. ред. Л.Я. Дятченко: Научные ведомости Белгородского государственного университета. — Белгород: ИПК НИУ «БелГУ», 2011

Чеканов Н.А.: Применение дифференциальных уравнений в курсе общей физики. — Белгород: НИУ БелГУ, 2011

Абрамочкин Е.Г.: Современная оптика гауссовых пучков. — М.: ФИЗМАТЛИТ, 2010

Алексеев Г.В.: Оптимизация в стационарных задачах тепломассопереноса и магнитной гидродинамики. — М.: Научный мир, 2010

Бирман М.Ш.: Избранные труды. — Ижевск: Ижевский институт компьютерных исследований, 2010

Козлов В.В.: Избранные работы по математике, механике и математической физике. — Ижевск: НИЦ «Регулярная и хаотическая динамика» ; , 2010

Лакс П.Д.: Гиперболические дифференциальные уравнения в частых производных. — Ижевск: Ижевский ин-т компьютерных исследований, 2010

Николаевский В.Н.: Собрание трудов. Геомеханика. — Ижевск: Ижевский институт компьютерных исследований, 2010

НИУ БелГУ ; гл. ред. Л.Я. Дятченко: Научные ведомости Белгородского государственного университета. — Белгород: БелГУ, 2010

Новосадов Б.К.: Методы математической физики молекулярных систем. — М.: ЛИБРОКОМ, 2010

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

2 В нелинейном ОДУ функция f(x,y) нелинейным образом зависит от функции у(х).

3 Интегральные кривые не пересекаются. В противном случае, в данной точке (х,у) можно было бы указать две разные касательные =f(x,у), что противоречит предположению об однозначности f(x,y).

4 Подробное рассмотрение вопросов устойчивости различных алгоритмов решения ОДУ можно найти, например, в [5]

5 Впрочем, авторы статьи [6] пишут, что «в настоящее время численное моделирование в физике конденсированного состояния достигло такого уровня доверия, что иногда возникают сомнения в результатах лабораторного эксперимента, если результаты выполненных экспериментов ему противоречат.»

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

Министерство образования Республики Беларусь

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

ИНФОРМАТИКИ И РАДИОЭЛЕКТРОНИКИ

Кафедра вычислительной математики и программирования

к курсовому проекту

« Решение дифференциальных уравнений по методу Эйлера »

2. Математическое объяснение метода

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

2.2 Исправленный метод Эйлера

2.3 Модифицированный метод Эйлера

3. Блок-схема алгоритма программы

4. Описание программы

Список использованной литературы

Приложение 1 (Текст программы)

Приложение 2 (Результаты работы программы)

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

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

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

2. Математическое объяснение метода

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

Решить дифференциальное уравнение у / =f(x,y) численным методом — это значит для заданной последовательности аргументов х₀ , х₁ …, х_n и числа у₀ , не определяя функцию у=F(x), найти такие значения у₁ , у₂ ,…, у_n , что у_i =F(x_i )(i=1,2,…, n) и F(x₀ )=y₀ .

Таким образом, численные методы позволяют вместо нахождения функции

У=F(x) получить таблицу значений этой функции для заданной последовательности аргументов. Величина h=x_k -x_k-1 называется шагом интегрирования.

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

Рассмотрим дифференциальное уравнение первого порядка

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

Требуется найти решение уравнения (1) на отрезке [а,b].

Разобьем отрезок [a, b] на n равных частей и получим последовательность х₀ , х₁ , х₂ ,…, х_n , где x_i =x₀ +ih (i=0,1,…, n), а h=(b-a)/n-шаг интегрирования.

В методе Эйлера приближенные значения у(х_i )»y_i вычисляются последовательно по формулам у_i +hf(x_i , y_i ) (i=0,1,2…).

При этом искомая интегральная кривая у=у(х), проходящая через точку М₀ (х₀ , у₀ ), заменяется ломаной М₀ М₁ М₂ … с вершинами М_i (x_i , y_i ) (i=0,1,2,…); каждое звено М_i M_i+1 этой ломаной, называемой ломаной Эйлера, имеет направление, совпадающее с направлением той интегральной кривой уравнения (1), которая проходит через точку М_i , смотри рисунок 1.

Если правая часть уравнения (1) в некотором прямоугольнике R<|x-x₀ |£a, |y-y₀ |£b>удовлетворяет условиям:

|df/dx|=|df/dx+f(df/dy)| £ M (M=const),

то имеет место следующая оценка погрешности:

где у(х_n )-значение точного решения уравнения(1) при х=х_n , а у_n — приближенное значение, полученное на n-ом шаге.

Формула (3) имеет в основном теоретическое применение. На практике иногда оказывается более удобным двойной просчет : сначала расчет ведется с шагом h, затем шаг дробят и повторный расчет ведется с шагомh/2. Погрешность более точного значения у_n * оценивается формулой

Метод Эйлера легко распространяется на системы дифференциальных уравнений и на дифференциальные уравнения высших порядков. Последние должны быть предварительно приведены к системе дифференциальных уравнений первого порядка.

2.2 Исправленный метод Эйлера

В исправленном методе Эйлера мы находим средний тангенс наклона касательнй для двух точек: 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₂ . Усреднение двух тангенсов дает прямую L’₃ . Наконец, через точку x_m , y_m мы проводим прямую L₃ параллельную L’₃ . Точка, в которой прямая L₃ пересечется с ординатой, восстановленной из x= x_m+1 =x_m +h, и будет искомой точкой y= y_m+1 = y_m +hy’_m . Тангенс угла наклона L₃ равен:

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

Соотношения 5, 6, 7 и 8 описывают исправленный метод Эйлера. (рис. 2)

2.3 Модифицированный метод Эйлера

Этот метод более точен. Рассмотрим дифференциальное уравнение (1) с начальным условием y(x₀ )=y₀ . Разобьем наш участок интегрирования на n равных частей. На малом участке [x₀ ,x₀ +h] интегральную кривую заменимпрямой линией. Получаем точкуМ_к (х_к ,у_к ). (рис. 3)

N_k / _y=y(x)

М_к М_к /

Y_k+1

Y_k

х_к х_к1/2 x_k+h =x_k1 X

Получаем точку N_k / . В этой точке строим следующую касательную:

Из точки М_к проводим прямую с угловым коэффициентом α_к и определяем точку пересечения этой прямой с прямой x_к1 . Получаем точку М_к / . В качестве у_к+1 принимаем ординату точки М_к / . Тогда:

Эти формулы называются рекуррентными формулами метода Эйлера.

Сначала вычисляют вспомогательные значения искомой функции y_к+1/2 в точках x_к+1/2 , затем находят значение правой части уравнения (1) в средней точке

3. Блок-схема алгоритма

где A — начальное значение x, B — конечное значение x, F(x) — значение функции в точке x_n , N — количество промежутков, st – выбор операции, C1,C2,C3 – константы для формул, nom — сохраняет номер используемой функции.

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

Подсчитывая каждый раз новое значение уравнения F(x), получаем последовательность значений x_n y_n , n=1,2,…

По этим значениям строим график.

4. Описание программы

Программа весьма проста. В ней много предусмотрено моментов неправильного ввода данных, о которых программа предупреждает пользователя и сразу же просит повторно ввести данные.

С самого начала программа предоставляет пользователю меню выполняемых функций, которые выделяются при помощи стрелок ↑ и ↓ выбор клавишей Enter:

После запуска программы нужно выбрать Formyla -> Enter, эта опция позволит из предложенного списка формул выбрать одну, по которой компьютер будет производить расчет и строить график. Все предложенные формулы имеют номерацию; чтобы выбрать интересующий вас пример нажмите на цифру равную номеру примера, и сразу же появится новое окно, в котором сверху будет записан ваш пример. Также в окне будет этот же пример но с нулями на месте констант. Под примером будет высвечена большая буква С, это используется для ввода констант. Для этого вам нужно нажать номер константы, он появится, и после знака равно запишите чему она равна (вводятся целые и вещественные значения). По окончании набора нажать Enter. Операцию повторять пока не будут введены все числа. По окончании нажать Esc. После появится строчка «уточните границы изменения Х, от A= до B= » здесь нужно занести данные на каком промежутке абсциссы будет рассматриваться функция. Следующая строчка попросит ввести начальные данные y(A)=. Следующей строчкой будет вопрос: «сохранить данные в файле? Да/Нет» ответить на этот вопрос с помощью клавиш Д и Н (рус), после чего программа вернется в первоначальное меню. Если данные были сохранены (в папке с программой появляется файл form.txt), то в следующий раз чтобы не набирать снова выберите в меню опцию Formyla -> Open in fails и на экране появятся введенные данные с пометкой снизу, сообщая что данные были прочитаны из файла.

Следующая опция Reshenie. После нажатия в окне просят ввести N(целое число) – число промежутков, на которые разделится рассматриваемый участок (ось ОХ). После появится таблица рассчитанных данных (номер точки, значение абсциссы, значение ординаты). При нажатии любой клавиши произойдет переход в меню.

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

Опция Exit применяется для выхода из программы.

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

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

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

Список использованной литературы

1. Д. Мак-Кракен, У. Дорн. Численные методы и программирование на фортране. –М.: Мир,1977.-389,396-408 с.

2. А.А. Самарский. Введение в численные методы. – М.:Наука,1987.-176 с.

3. Алгоритмы вычислительной математики: Лабораторный практикум по курсу «Программирование» для студентов 1 — 2-го курсов всех специальностей БГУИР/А.К. Синицын, А.А. Навроцкий.- Мн.: БГУИР, 2002.- 65-69 с.

4. ГОСТ 2.105-95. Общие требования к текстовым документам.

5. ГОСТ 7.32-91. Система стандартов по информации, библиотечному и издательскому делу. Отчет о НИР. Структура и правила оформления.

Приложение 1. Текст программы.

void formyl(int p)

if(p==1) printf(«\n 1. C1*y’ = C2*y + C3*x + C4*x*y»);

else if(p==2) printf(«\n 2. y’/(C1-100) = C2*y + C3*x + (C4+x)*y»);

else if(p==3) printf(«\n 3. pow(e,C1)*y’ = C2*y + C3*cos(x) + (C4+x+y)»);

else if(p==4) printf(«\n 4. C1*sin(x)*y’ = e*C2*y + C3*arcsin(x) + C4*y/x»);

else if(p==5) printf(«\n 5. C1*y’ = sin(C2)*y + tg(C3*x) + C4*ln(x)*y»);

else if(p==6) printf(«\n 6. C1*y’ = y*C2 + C3*sin(x) + C4*cos(x)*y»);

else if(p==7) printf(«\n 7. (C1+C2+C3+C4)*y’ = C2*y + (C3-x) + lg(C4*x)*y»);

else if(p==8) printf(«\n 8. y’/C1 = y/C2 + C3*sin(x) + C4*x*y»);

else if(p==9) printf(«\n 9. sin(C1)*y’ = C2*y + |C3|*x + x*y/C4»);

void formyl2(int p,double C1,double C2,double C3,double C4)

else if(p==2) printf(«y’/(%.2f-100)=%.2f*y+%.2f*x+(%.2f+x)*y»,C1,C2,C3,C4);else if(p==3) printf(«pow(e,%.2f)*y’=%.2f*y+%.2f*cos(x)+(%.2f+x+y)»,C1,C2,C3,C4);

else if(p==4) printf(«%.2f*sin(x)*y’=e*%.2f*y+%.2f*arcsin(x)+%.2f*y/x»,C1,C2,C3,C4);

else if(p==5) printf(«%.2f*y’=sin(%.2f)*y+tg(%.2f*x)+%.2f*ln(x)*y»,C1,C2,C3,C4);

else if(p==6) printf(«%.2f*y’=y*%.2f+%.2f*sin(x)+%.2f*cos(x)*y»,C1,C2,C3,C4);

else if(p==8) printf(«y’/%.2f=y/%.2f+%.2f*sin(x)+%.2f*x*y»,C1,C2,C3,C4);

else if(p==9) printf(«sin(%.2f)*y’=%.2f*y+|%.2f|*x+x*y/%.2f»,C1,C2,C3,C4);

double formyl3(int p,double h,double x,double y,double C1,double C2,double C3,double C4)

else if(p==2) Y=h*(C1-100)*(y*C2+C3*x+(C4+x)*y)+y;

else if(p==3) Y=h*(C2*y+C3*cos(x)+C4+x+y)/exp(C1)+y;

else if(p==4) Y=h*(exp(1)*C2*y+C3*asin(x)+C4*y/x)/(C1*sin(x))+y;

else if(p==5) Y=h*(sin(C2)*y+tan(C3*x)+C4*log10(x)*y)/C1+y;

else if(p==6) Y=h*(y*C2+C3*sin(x)+C4*cos(x)*y)/C1+y;

else if(p==7) Y=h*(C2*y+(C3-x)+log10(C4*x)*y)/(C1+C2+C3+C4)+y;

else if(p==8) Y=h*(y/C2+C3*sin(x)+C4*x*y)*C1+y;

else if(p==9) Y=h*(C2*y+abs(C3)*x+x*y/C4)/sin(C1)+y;

int vv=0,vv1=0; // руководит операциями

int N=0,W; // кол промежутков

int i,j,k; // используются во всех «for»

int nom; // номер примера

int st=4,vst=0; // строчка в меню

double C1,C2,C3,C4; // константы

double M; // масштаб

double xtoch,ytoch; // считает y(x) по графику

double A=0,B=0,ii,jj,kk; // пределы интегрирования

double x[102],y[102]; // главные переменные x,y

int g_driver=9,g_mode=2, g_error;

printf(«\n error=%d, reason=%s\n», g_error, grapherrormsg(g_error));