Skip to content

Instantly share code, notes, and snippets.

Show Gist options
  • Save anonymous/183e075a3ef33d8530f3524eba3c9127 to your computer and use it in GitHub Desktop.
Save anonymous/183e075a3ef33d8530f3524eba3c9127 to your computer and use it in GitHub Desktop.
Метод эйлера решения дифференциальных уравнений

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



Метод Эйлера Метод Гюна Методы Рунге-Кутты. Метод Адамса Метод прогноза и коррекции. При решении научных и инженерно-технических задач часто бывает необходимо математически описать какую-либо динамическую систему. Лучше всего это делать в виде дифференциальных уравнений ДУ или системы дифференциальных уравнений. Наиболее часто они такая задача возникает при решении проблем, связанных с моделированием кинетики химических реакций и различных явлений переноса тепла, массы, импульса — теплообмена, перемешивания, сушки, адсорбции, при описании движения макро- и микрочастиц. Обыкновенным дифференциальным уравнением ОДУ n-го порядка называется следующее уравнение, которое содержит одну или несколько производных от искомой функции y x:. В ряде случаев дифференциальное уравнение можно преобразовать к виду, в котором старшая производная выражена в явном виде. Такая форма записи называется уравнением, разрешенным относительно старшей производной при этом в правой части уравнения старшая производная отсутствует:. Именно такая форма записи принята в качестве стандартной при рассмотрении численных методов решения ОДУ. Линейным дифференциальным уравнением называется уравнение, линейное относительно функции y x и всех ее производных. Решением обыкновенного дифференциального уравнения называется такая функция y x , которая при любых х удовлетворяет этому уравнению в определенном конечном или бесконечном интервале. Процесс решения дифференциального уравнения называют интегрированием дифференциального уравнения. Общее решение ОДУ n -го порядка содержит n произвольных констант C 1 , C 2 , …, C n. Это очевидно следует из того, что неопределенный интеграл равен первообразной подынтегрального выражения плюс константа интегрирования. Так как для решения ДУ n -го порядка необходимо провести n интегрирований, то в общем решении появляется n констант интегрирования. Частное решение ОДУ получается из общего, если константам интегрирования придать некоторые значения, определив некоторые дополнительные условия, количество которых позволяет вычислить все неопределенные константы интегрирования. Точное аналитическое решение общее или частное дифференциального уравнения подразумевает получение искомого решения функции y x в виде выражения от элементарных функций. Это возможно далеко не всегда даже для уравнений первого порядка. Численное решение ДУ частное заключается в вычислении функции y x и ее производных в некоторых заданных точках , лежащих на определенном отрезке. То есть, фактически, решение ДУ n -го порядка вида получается в виде следующей таблицы чисел столбец значений старшей производной вычисляется подстановкой значений в уравнение:. Например, для дифференциального уравнения первого порядка таблица решения будет представлять собой два столбца — x и y. Сами координаты при этом называют узлами сетки. Чаще всего, для удобства, используются равномерные сетки , в которых разница между соседними узлами постоянна и называется шагом сетки или шагом интегрирования дифференциального уравнения. Для определения частного решения необходимо задать дополнительные условия, которые позволят вычислить константы интегрирования. Причем таких условий должно быть ровно n. Для уравнений первого порядка — одно, для второго - 2 и т. В зависимости от способа их задания при решении дифференциальных уравнений существуют три типа задач:. Необходимо найти такое частное решение дифференциального уравнения, которое удовлетворяет определенным начальными условиям, заданным в одной точке:. Эта точка х 0 называется начальной. Например, если решается ДУ 1-го порядка, то начальные условия выражаются в виде пары чисел x 0 , y 0. Такого рода задача встречается при решении ОДУ , которые описывают, например, кинетику химических реакций. В качестве примера можно так же привести задачу о теплопереносе или массопереносе диффузии , уравнение движения материальной точки под действием сил и т. В этом случае известны значения функции и или ее производных в более чем одной точке, например, в начальный и конечный момент времени, и необходимо найти частное решение дифференциального уравнения между этими точками. Сами дополнительные условия в этом случае называются краевыми граничными условиями. Естественно, что краевая задача может решаться для ОДУ не ниже 2-го порядка. Ниже приведен пример ОДУ второго порядка с граничными условиями заданы значения функции в двух различных точках:. Задачи этого типа похожи на краевую задачу. При их решении необходимо найти, при каких значениях какого-либо параметра решение ДУ удовлетворяет краевым условиям собственные значения и функции, которые являются решением ДУ при каждом значении параметра собственные функции. Например, многие задачи квантовой механики являются задачами на собственные значения. Рассмотрим некоторые численные методы решения задачи Коши начальной задачи обыкновенных дифференциальных уравнений первого порядка. Запишем данное уравнение в общем виде, разрешенном относительно производной правая часть уравнения не зависит от первой производной:. Сложность, однако, заключается в том, что интеграл в правой части есть интеграл от неявно заданной функции, нахождение которого в аналитическом виде в общем случае невозможно. Из множества разработанных для решения ОДУ первого порядка методов рассмотрим методы Эйлера , Рунге-Кутта и прогноза и коррекции. Они достаточно просты и дают начальное представление о подходах к решению данной задачи в рамках численного решения задачи Коши. Исторически первым и наиболее простым способом численного решения задачи Коши для ОДУ первого порядка является метод Эйлера. В его основе лежит аппроксимация производной отношением конечных приращений зависимой y и независимой x переменных между узлами равномерной сетки:. Сравнивая формулу Эйлера с общим выражением, полученным ранее, видно, что для приближенного вычисления интеграла в 6. Графическая интерпретация метода Эйлера также не представляет затруднений см. Действительно, исходя из вида решаемого уравнения 6. Если искомая функция сильно отличается от линейной на отрезке интегрирования, то погрешность вычисления будет значительной. Ошибка метода Эйлера прямо пропорциональна шагу интегрирования:. Процесс вычислений строится следующим образом. При заданных начальных условиях x 0 и y 0 можно вычислить. Таким образом, строится таблица значений функции y x с определенным шагом h по x на отрезке [x 0 , x N ]. Ошибка в определении значения y x i при этом будет тем меньше, чем меньше выбрана длина шага h что определяется точностью формулы интегрирования. При больших h метод Эйлера весьма неточен. Он дает все более точное приближение при уменьшении шага интегрирования. Данное уравнение уже записано в стандартном виде, резрешенном относительно производной искомой функции. Для первых четырех узлов сетки вычисления будут следующими:. Во второй колонке таблицы для сравнения приведены значения, вычисленные по аналитическому решению данного уравнения. Во второй части таблицы приведена относительная погрешность полученных решений. Таблица 1 Решение уравнения 6. Относительные погрешности вычисленных значений функции при различных h. В таблице 1 приведены для сравнения вычисления для некоторых других значений N , вплоть до Поэтому метод Эйлера практически не используется в вычислительной практике. Точность метода Эйлера можно повысить, если воспользоваться для аппроксимации интеграла более точной формулой интегрирования — формулой трапеций. Таким образом получается метод Гюна или метод Эйлера с пересчетом. Для каждого узла интегрирования производится следующая цепочка вычислений. Благодаря более точной формуле интегрирования, погрешность метода Гюна пропорциональна уже квадрату шага интегрирования. Подход, использованный в методе Гюна, используется для построения так называемых методов прогноза и коррекции , которые будут рассмотрены позже. Что намного точнее значения, полученного методом Эйлера при том же шаге интегрирования. Отметим существенное увеличение точности вычислений метода Гюна по сравнению с методом Эйлера. Такая же точность вычислений по формуле Эйлера достигается при числе отрезков интегрирования N примерно Следовательно, при использовании метода Гюна при одинаковой точности вычислений понадобится примерно в 15 раз меньше времени ЭВМ, чем при использовании метода Эйлера. Решение ОДУ в некоторой точке x i называется устойчивым, если найденное в этой точке значение функции y i мало изменяется при уменьшении шага интегрирования. Для проверки устойчивости, таким образом, надо провести два расчета значения y i — с шагом интегрирования h и при уменьшенной например, двое величине шага. Такая проверка может осуществляться и для всех решений на всем интервале значений x. Если условие не выполняется, то шаг снова делится пополам и находится новое решение и т. Дальнейшее улучшение точности решения ОДУ первого порядка возможно за счет увеличения точности приближенного вычисления интеграла в выражении 6. Мы уже видели, какое преимущество дает переход от интегрирования по формуле прямоугольников метод Эйлера к использованию формулы трапеций метод Гюна при аппроксимации этого интеграла. Чтобы воспользоваться этой формулой, надо использовать некоторое приближение для вычисления этих значений. При использовании различных методов приближенного вычисления этих величин, получаются выражения для методов Рунге-Кутты различного порядка точности. Алгоритм Рунге-Кутты третьего порядка - РК3 погрешность порядка h Алгоритм Рунге-Кутты четвертого порядка- РК4 погрешность порядка h Алгоритмы третьего и четвертого порядков требуют на каждом шаге трех и четырех вычислений функции соответственно, но являются весьма точными. Используя алгоритм Рунге-Кутты третьего 6. Как видно, точность решения, полученного методом Рунге-Кутты четвертого порядка, намного превышает точность решения, полученного методом Гюна и методом Рунге-Кутты 3-го порядка сравнте с данными Таблицы 2. Высокая точность, вместе с достаточной простотой реализации делает метод Рунге-Кутты четвертого порядка одним из весьма распространенных численных методов решения задачи Коши ОДУ и систем ОДУ первого порядка. Рассмотренные ранее методы Эйлера, Гюна, Рунге-Кутты используют значение функции на одном предшествующем шаге, поэтому они относятся к так называемым одношаговым методам. Если используются значения в k предыдущих узлах, то говорят о k-шаговом методе интегрирования уравнения. Одним из способов построения многошаговых методов заключается в следующем. По значениям функции, вычисленным в k предшествующих узлах, строится интерполяционный полином степени k-1 - , который используется при интегрировании дифференциального уравнения по выражению 6. Интеграл при этом выражается через квадратурную формулу:. Значения квадратурных коэффициентов для k от 2 до 4 приведены в таблице. Полученное таким образом семейство формул называется явной k -шаговой схемой Адамса методы Адамса-Башфорта. Однако обычно это уравнение не решается, а значение в правой части заменяется на рассчитанное по какой-либо явной формуле - например, формуле Адамса-Башфорта. Такой подход лежит в основе методов "прогноза-коррекции". Достоинством многошаговых методов Адамса при решении ОДУ заключается в том, что в каждом узле рассчитывается только одно значение правой части ОДУ - функции F x,y. К недостаткам можно отнести невозможность старта многошагового метода из единственной начальной точки, так как для вычислений по k -шаговой формуле необходимо знание значения функции в k узлах. Поэтому приходится k-1 решение в первых узлах x 1 , x 2 , …, x k-1 получать с помощью какого-либо одношагового метода, например метода Рунге-Кутты 4—го порядка. Другой проблемой является невозможность изменения шага в процессе решения, что легко реализуется в одношаговых методах. Несколько иной подход используется в многошаговых методах прогноза и коррекции. В качестве иллюстрирующего примера рассмотрим 2-х шаговый метод прогноза и коррекции. Затем это значение корректируется по более точной формуле, в данном случае — по формуле трапеций неявная формула Адамса второго порядка. Обычно значение в узле x 1 определяется каким-либо одношаговым методом методом Рунге-Кутты или Гюна. На каждом шаге построения решения методом прогноза и коррекции требуется вычислить всего одно значение функции, а одно берется из предыдущего узла сетки. Поэтому он весьма экономичен по затратам времени вычислений при достаточной точности. Аналогичные схемы прогноза-коррекции могут быть получены сочетанием явных прогноз и неявных коррекция формул Адамса для различных k. Решение задачи Коши для систем дифференциальных уравнений 1-го порядка. Системой M дифференциальных уравнений первого порядка в общем случае можно назвать следующую совокупность Обыкновенных дифференциальных уравнений:. Начальными условиями при решения задачи Коши для такой системы будут являться значение независимой переменной и значения всех M функций при этом значении:. Все описанные ранее методы решения задачи Коши для уравнений легко обобщаются на случай решения систем ДУ первого порядка. Формулы выбранного метода применяются последовательно к каждому уравнению системы уравнений для определения значения соответствующей функции. Из первого уравнения определяется значение y 1 i , из второго — y 2 i , …, из M -го - y Mi. В качестве примера рассмотрим применение метода Рунге-Кутты 4-го порядка для решения системы двух ОДУ первого порядка. Адаптируем формулу Рунге-Кутты 4-го порядка для данной системы уравнений. Из первого уравнения будем вычислять значения функции u x , а из второго — функции v x это функции, чьи производные стоят в левой части соотетствующих уравнений:. С учетом вышесказанного, для расчета коэффициентов k u0 - k u3 используем правую часть первого уравнения F 1 , а для коэффициентов k v0 - k v 3 - второго F 2. Кроме этого, для расчета приращения функции u используем коэффициенты k u , а для расчета приращения функции v - k v. Таким образом, коэффициенты рассчитываются по следующим формулам:. Данная модель довольно широко применяется при описании временной зависимости объема популяций в биологических системах, при моделировании экономических и физических процессов. Задача формулируется следующим образом. Пусть в системе в некоторый момент времени t имеются хищники например, волки в количестве v t и жертвы например, зайцы в количестве u t. Модель "хищник — жертва" утверждает, что u t и v t удовлетворяют системе ОДУ первого порядка:. Где A , B , C и D — некоторые числовые константы. Действительно, если зайцы имеют достаточно травы для питания, то скорость роста популяции будет прямо пропорциональна их числу первое слагаемое в первом уравнении. Второе слагаемое описывает гибель зайцев при встрече с хищниками, так как вероятность их встречи равна произведению. Второе уравнение описывает изменение популяции хищников. Скорости роста популяции способствует их хорошее питание первое слагаемое второго уравнения пропорционально вероятности встречи хищника и жертвы - , а избыток хищников приводит к их гибели за счет голода второе слагаемое. Применим метод Рунге-Кутты 4-го порядка для решения полученной системы уравнений. Значения функций u и v находятся по уже известным формулам:. Решение задачи Коши для дифференциальных уравнений второго и более высоких порядков. Аналогично, ОДУ порядка n сведется к системе из n дифференциальных уравнений первого порядка. Движение тела под действием пружины Рассмотрим некоторое материальное тело массой m , которое движется по горизонтальной поверхности в общем случае — с трением под действием пружины. Сила трения направлена всегда против движения тела и пропорциональна его скорости:. Баланс сил, действующих на тело в каждый момент времени можно записать так:. С учетом того, что координата тела есть функция от времени, а скорость и ускорение — это первая и вторая производная координаты во времени, соответственно, получим:. Таким образом, изменение координаты тела от времени описывается ОДУ 2-го порядка, которое в стандартном виде записывается так:. Как было показано, ОДУ 2-го порядка сводится к системе двух уравнений 1-го порядка подстановкой:. Правые части уравненией имеют вид:. Из первого уравнения рассчитываем значения функции x t , из второго — v t. Теперь запишем соответствующие выражения для расчета коэффициентов k x0 — k x3 и k v0 - k v3. Решение других проблем, связанных с дифференциальными уравнениями — задачи с граничными условиями и задачи на собственные значения и функции в данном курсе не рассматриваются. Обыкновенным дифференциальным уравнением ОДУ n-го порядка называется следующее уравнение, которое содержит одну или несколько производных от искомой функции y x: Такая форма записи называется уравнением, разрешенным относительно старшей производной при этом в правой части уравнения старшая производная отсутствует: Общее решение ОДУ n -го порядка содержит n произвольных констант C 1 , C 2 , …, C n Это очевидно следует из того, что неопределенный интеграл равен первообразной подынтегрального выражения плюс константа интегрирования Так как для решения ДУ n -го порядка необходимо провести n интегрирований, то в общем решении появляется n констант интегрирования. То есть, фактически, решение ДУ n -го порядка вида получается в виде следующей таблицы чисел столбец значений старшей производной вычисляется подстановкой значений в уравнение: В зависимости от способа их задания при решении дифференциальных уравнений существуют три типа задач: Необходимо найти такое частное решение дифференциального уравнения, которое удовлетворяет определенным начальными условиям, заданным в одной точке: Ниже приведен пример ОДУ второго порядка с граничными условиями заданы значения функции в двух различных точках: Численные методы решения задачи Коши ОДУ первого порядка Рассмотрим некоторые численные методы решения задачи Коши начальной задачи обыкновенных дифференциальных уравнений первого порядка. Запишем данное уравнение в общем виде, разрешенном относительно производной правая часть уравнения не зависит от первой производной: Метод Эйлера Исторически первым и наиболее простым способом численного решения задачи Коши для ОДУ первого порядка является метод Эйлера. В его основе лежит аппроксимация производной отношением конечных приращений зависимой y и независимой x переменных между узлами равномерной сетки: Ошибка метода Эйлера прямо пропорциональна шагу интегрирования: При заданных начальных условиях x 0 и y 0 можно вычислить Таким образом, строится таблица значений функции y x с определенным шагом h по x на отрезке [x 0 , x N ]. Используя метод Эйлера, построить приближенное решение для следующей задачи Коши: Для первых четырех узлов сетки вычисления будут следующими: И так далее, до конца отрезка. Проведем вычисления для уравения 6. Таблица 2 Решение уравнения методами Эйлера и Гюна x Точное Метод Гюна Метод Эйлера y отн. Проверка устойчивости решения Решение ОДУ в некоторой точке x i называется устойчивым, если найденное в этой точке значение функции y i мало изменяется при уменьшении шага интегрирования. Методы Рунге-Кутты Дальнейшее улучшение точности решения ОДУ первого порядка возможно за счет увеличения точности приближенного вычисления интеграла в выражении 6. Алгоритм Рунге-Кутты третьего порядка - РК3 погрешность порядка h 3: Многошаговые методы Метод Адамса Рассмотренные ранее методы Эйлера, Гюна, Рунге-Кутты используют значение функции на одном предшествующем шаге, поэтому они относятся к так называемым одношаговым методам. Интеграл при этом выражается через квадратурную формулу:


РЕШЕНИЕ ОБЫКНОВЕННЫХ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ


Обыкновенные дифференциальные уравнения встречаются достаточно часто в различных прикладных задачах. Ими описываются задачи движения систем материальных точек, электрических цепей и др. Порядок старшей производной, входящей в это уравнение, называется порядком дифференциального уравнения. Поиск решения такого уравнения называется его интегрированием, а полученное решение — интегралом этого уравнения. Частное решение определяется конкретным значением С. Задача отыскания решения обыкновенного дифференциального уравнения при начальных условиях называется задачей Коши для обыкновенного дифференциального уравнения. Их можно свести к системе n обыкновенных дифференциальных уравнений первого порядка заменой на неизвестную функцию Р 1 х , на Р 2 х и т. Численное решение состоит в построении таблицы приближенных значений y 1 ,y 2 ,…,y n -1 точного решения y x уравнения. Мы не знаем вида функции, не можем сразу вычислить значение функции. Первый шаг — разбиваем отрезок на конечное число узловых точек узлы сетки. Второй шаг — зная начальные условия вычисляем значение первой производной в точке Х 0. Геометрическая интерпритация метода состоит в замене интегральной кривой на отрезке касательной к ней в точке x i ,y i. На каждом шаге заново определяется касательная, и, следовательно, соответствующая приближенному решению кривая будет ломаной линией. Поэтому метод называют еще методом ломаных. Это схема вычислений для уравнения 1 порядка. Метод Эйлера для уравнений 2 порядка надо уравнение с начальными условиями свести к системе 2 уравнений 1-го порядка:. На практике вычисления проводятся снизу вверх: Для уравнения второго порядка приходится дважды применять интегрирование. Восстанавливаем значение производной, а затем значение искомой функции, то есть решаем вначале одно уравнение, а затем второе. На отрезке [a,b] составить таблицу значений приближенного решения дифференциального уравнения , удовлетворяющего начальным условиям. Функция y зависит от аргумента х. Точность метода или ошибкой обрыва называют ошибку, которую делают при переходе от предыдущего к последующему Х, если заменяют дифференциальное уравнение конечностным выражением. Принято считать, если расчетные формулы численного метода согласуются с разложением в ряд Тейлора до членов порядка h p , то число р называют порядком точности метода. В методе Эйлера сравнение формулы вычисления с разложением в ряд Тейлора согласуется до первого порядка по h, то есть метод имеет первый порядок точности с локальной погрешностью O h 2. При решении задач на ЭВМ важное значение имеет наглядность и удобство быстрого восприятия изучаемых явлений. После освоения средств графического режима Турбо Паскаля можно строить графики, а также получать различные геометрические фигуры. Какие бы изображения не выводились на экран, все они строятся из точек, группы которых в свою очередь образуют отрезки и кривые. Вывод точки осуществляется процедурой PutPixel X, Y, Color , где Х и Y-координаты экрана, где будет расположена точка. Color — ее цвет. Здесь введены имена переменных: Шаг изменения аргумента X определяется выражением:. В организованном цикле табулирования функции вычисляются значения Y x , используя подпрограмму — функцию с именем F. Индекс каждого значения Х и соответствующее ему значение Y x , будут помещаться в определенные элементы массива У. Это достигается программным путем применением операторов. В том случае, если в окрестности точки разрыва функция получала бы положительное значение, следовало бы осуществлять ограничение сверху введением операторов. Для определения масштаба вычерчиваемого графика необходимо найти минимальное и максимальное значения элементов массива Y. Это достигается с помощью операторов. В результате выполнения этого фрагмента программы минимальное значение функции будет храниться в переменной yMin, а максимальное значение — в переменной yMax. Эти значения используются для задания масштаба, чтобы минимальное и максимальное значения функции изображались точками, отстоящими друг от друга не более чем на ширину экрана. Кроме того, они используются при формировании вертикальной и горизонтальной сетки, а также при задании оси абсцисс. Зависимости, по которым определяются эти параметры, приводятся в программе с подробными комментариями. Для того, чтобы установить графический режим, необходимо использовать стандартную процедуру InitGraph модуля Graph, которая устанавливает один из возможных графических режимов. В модуле Graph имеется более 70 подпрограмм. SetTextStule Font, Direct, CharSize ; —процедура вывода текста с указанием вида шрифта, наклона и его размера. Вид шрифта параметр Font в модуле Graph определяется константами:. Для определения направления вывода текста используются константы: Третий параметр определяет размер шрифта. Максимальное значение константы равно При формировании оси X, горизонтальной и вертикальной сетки можно использовать следующие процедуры:. Для установки курсора в указанную координату используется процедура MoveTo, X, Y , а для вычерчивания линии от текущей линии до указанной координаты — процедура LineTo X, Y. При работе с графическими изображениями необходимо использовать систему координат. В зависимости от выбранного графического режима устанавливается разрешающая способность экрана в х, х или другое число точек пикселей. Программным путем можно установить максимальное число пикселей в строке с помощью функции GetMaxX, а число строк— с помощью функции GetMaxY. Начало координат имеет значение 0, 0 , центр экрана — GetMaxX div 2, GetMaxY div 2 , а правый нижний угол — GetMaxX, GetMaxY. При выполнении графических операций всегда указывается, что должен делать курсор, однако сам курсор в графическом режиме невидимый. Сдача сессии и защита диплома - страшная бессонница, которая потом кажется страшным сном. Cуть методов численного интегрирования Http: Методика изучения ценностных ориентаций М. Организационные методы подходы II Основные теоретико-методологические понятия. Дисконтные методы анализа эффективности инвестиционных проектов. Метод начальных параметров II. Но предоставляет возможность бесплатного использования. Есть нарушение авторского права? Метод Эйлера Обыкновенные дифференциальные уравнения встречаются достаточно часто в различных прикладных задачах. В неявном виде дифференциальное уравнение записывается так: Существуют общее и частное решения дифференциального уравнения. Мы рассматриваем уравнения разрешенными относительно старшей производной. Для решения мы будем рассматривать методы Эйлера и Рунге-Кутта. Решение обыкновенных дифференциальных уравнений.


https://gist.github.com/258de9677d49512029e5a3c835ba8250
https://gist.github.com/eced23886d10ff0efc7672bab210b758
https://gist.github.com/56fe3ec191c455cc7691c5058a34d404
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment