Решение дифференциальных уравнений в matlab
Решение ОДУ в Matlab
Доброго времени суток! Сегодня мы поговорим о решении ОДУ (обыкновенных дифференциальных уравнений) в Matlab. Перед тем как мы начнём обсуждать данную тему, советую вам ознакомиться с темой: Численное дифференцирование в Matlab, чтобы лучше понимать теоретическую составляющую решения ОДУ.
Обыкновенные дифференциальные уравнения
С помощью дифференциальных уравнений можно описать разные задачи: движения системы, взаимодействующих материальных точек, химической кинетики и т.д. Различают три типа задач для систем диф. уравнений:
- Задача Коши
- Краевая задача
- Задача на собственные значения
Кратко расскажу о их сути:
Задача Коши предполагает дополнительные условия в виде значения функции в определённой точке.
Краевая задача подразумевает поиск решения на заданном отрезке с краевыми (граничными) условиями в концах интервала или на границе области.
Задача на собственные значения — помимо искомых функций и их производных, в уравнение входят дополнительное несколько неизвестных параметров, которые являются собственными значениями.
Методы решения дифференциальных уравнений
Решение ОДУ в Matlab и не только, в первую очередь, сводится к выбору порядка численного метода решения. Порядок численного метода не связан с порядком дифференциального уравнения. Высокий порядок у численного метода означает его скорость сходимости.
В случае большого интервала, с помощью алгоритмов с низким порядком сжимают интервал с решениями и находят приблизительные корни, а затем уже уточняют корни с помощью методов с высоким порядком.
Решение обыкновенных дифференциальных уравнений в Matlab можно реализовать «своими ручками», прописав алгоритм по разным схемам. Но также в Matlab есть встроенные функции, выполняющие все стандартные задачи.
Метод Рунге-Кутта первого порядка
Методы Рунге-Кутта представляют собой разложения в ряд Тейлора и от количества использованных элементов ряда зависит порядок этого метода. Следовательно, помимо Рунге-Кутта первого порядка, вы сможете увидеть методы других порядков. Иногда их называют другими именами.
Например, Метод Рунге-Кутта первого порядка, также известен как Метод Эйлера или Метод ломаных. Информацию о его математическом и графическом представлении советую поискать в гугл. Мы же поговорим о том, как Метод Рунге-Кутта первого порядка реализуется в Matlab для решения ОДУ. Например:
Решить и привести график ошибки уравнения y’ = y*x методом Рунге-Кутта первого порядка (Методом Эйлера, Методом ломаных).
Погрешность Метода Рунге-Кутта 1 порядка
» data-medium-file=»https://i1.wp.com/codetown.ru/wp-content/uploads/2017/02/Рунге-1-погрешность.png?fit=300%2C236&ssl=1″ data-large-file=»https://i1.wp.com/codetown.ru/wp-content/uploads/2017/02/Рунге-1-погрешность.png?fit=622%2C489&ssl=1″ src=»https://i1.wp.com/codetown.ru/wp-content/uploads/2017/02/%D0%A0%D1%83%D0%BD%D0%B3%D0%B5-1-%D0%BF%D0%BE%D0%B3%D1%80%D0%B5%D1%88%D0%BD%D0%BE%D1%81%D1%82%D1%8C.png?resize=622%2C489″ alt=»Погрешность метода 1 порядка» w />
На данном графике показана зависимость величины ошибки от шага.
Метод Рунге-Кутта второго порядка
Также известен как Метод Эйлера-Коши. Как видите, во второй части уравнения происходит обращения к следующему шагу. Но как тогда быть, если нам ещё не известен следующий шаг? Всё просто. Метод Рунге-Кутта второго порядка — это всё тот же метод первого порядка, однако, на половине шага происходит нахождение «первичного» решения, а затем происходит его уточнение. Это позволяет поднять порядок скорости сходимости до двух.
Решить и привести график ошибки уравнения u’ = u*x методом Рунге-Кутта второго порядка.
По сравнению с Рунге-Куттом первого порядка изначальная ошибка уже гораздо меньше.
Мы не будем говорить о третьем порядке, потому что задачи на третий порядок встречаются редко, но если будет необходимо, пишите в комментариях, выложу.
Метод Рунге-Кутта четвёртого порядка
Метод Рунге-Кутта четвёртого порядка считается самым распространённым. Тем не менее, работает он аналогично второму и третьему порядку.
Решить и привести график ошибки уравнения u’ = u*x методом Рунге-Кутта четвёртого порядка.
Как видите, на последней картинке размерность ошибки на столько мала, что пришлось воспользоваться loglog() для лучшей видимости.
Решение ОДУ в Matlab стандартными средствами
Стоит отметить, что мы с вами разобрали только один самый известный метод решения ОДУ с разными порядками. Однако, методов очень много.
Для решения дифференциальных уравнений и систем в MATLAB предусмотрены следующие функции:
ode45 (f, interval, X0, [options])
ode23 (f, interval, X0, [options])
ode113 (f, interval, X0, [options])
ode15s (f, interval, X0, [options])
ode23s (f, interval, X0, [options])
ode23t (f, interval, X0, [options])
ode23tb (f, interval, X0, [options])
Входными параметрами этих функций являются:
- f — вектор-функция для вычисления правой части уравнения системы уравнений;
- interval — массив из двух чисел, определяющий интервал интегрирования дифференциального уравнения или системы;
- Х0 — вектор начальных условий системы дифференциальных уравнений;
- options — параметры управления ходом решения дифференциального уравнения или системы.
Все функции возвращают:
- массив Т — координаты узлов сетки, в которых ищется решение;
- матрицу X, i-й столбец которой является значением вектор-функции решения в узле Тi.
В функции ode45 реализован метод Рунге-Кутта 4-5 порядка точности, в функции ode23 также реализован метод Рунге-Кутта, но 2-3 порядка, а функция ode113 реализует метод Адамса.
Для решения жёстких систем предназначены функция ode15s, в которой реализован метод Гира, и функция ode23s, реализующая метод Розенброка. Для получения более точного решения жёсткой системы лучше использовать функцию ode15s. Для решения системы с небольшим числом жёсткости можно использовать функцию ode23t, а для грубой оценки подобных систем служит функция ode23tb.
Символьное решение обыкновенных дифференциальных уравнений произвольного порядка осуществляет функция dsolve r = dsolve(‘eq1,eq2,…’, ‘cond1,cond2,…‘, ‘v’)
Пример использования:
На этом мы закончим. Если остались вопросы, задавайте их в комментариях. Также вы можете скачать исходники чтобы лучше понять тему: «Решение ОДУ в Matlab».
Решение дифференциальных уравнений в среде MatLab
MatLab предоставляет возможности для решения обыкновенных дифференциальных уравнений произвольного порядка и систем с начальными условиями, т.е. задачи Коши.
Решатели реализуют следующие методы решения систем дифференциальных уравнений, причем для решения жестких систем уравнений рекомендуется использовать только специальные решатели ode15s , ode23s, ode23t, ode23tb:
· ode45 — одношаговые явные методы Рунге-Кутта 4-го и 5-го порядка. Это классический метод, рекомендуемый для начальной пробы решения. Во многих случаях он дает хорошие результаты;
· ode23 — одношаговые явные методы Рунге-Кутта 2-го и 4-го порядка. При умеренной жесткости системы ОДУ и низких требованиях к точности этот мето;. может дать выигрыш в скорости решения;
· ode113 — многошаговый метод Адамса-Башворта-Мултона переменного порядка Это адаптивный метод, который может обеспечить высокую точность решения
· ode23tb — неявный метод Рунге-Кутта в начале решения и метод, использующий формулы обратного дифференцирования 2-го порядка в последующем
Несмотря на сравнительно низкую точность, этот метод может оказаться более эффективным, чем ode15s;
· ode15s — многошаговый метод переменного порядка (от 1 до 5, по умолчанию 5), использующий формулы численного дифференцирования. Это адаптивный метод, его стоит применять, если решатель ode45 не обеспечивает решения;
· ode23s — одношаговый метод, использующий модифицированную формулу Розенброка 2-го порядка. Может обеспечить высокую скорость вычислений при низкой точности решения жесткой системы дифференциальных уравнений;
· ode23t — метод трапеций с интерполяцией. Этот метод дает хорошие результаты при решении задач, описывающих колебательные системы с почти гармоническим выходным сигналом;
Все решатели могут решать системы уравнений явного вида у’ = F(t, y). Решатели ode15s и ode23t способны найти корни дифференциально-алгебраических уравнений M(t)y’ = F(t, у), где М называется матрицей массы. Решатели ode15s, ode23s, ode23t и ode23tb могут решать уравнения неявного вида M(t,y) у’ = F(t, у).
· ode23tb, ode23s служат для решения жестких дифференциальных уравнений.
· ode15s — жестких дифференциальных и дифференциально-алгебраических уравнений,
· ode23t — умеренно жестких дифференциальных и дифференциально-алгебраических уравнений.
В описанных далее функциях для решения систем дифференциальных уравнений приняты следующие обозначения и правила:
· options — аргумент, создаваемый функцией odeset — позволяет вывести параметры, установленные по умолчанию;
· tspan — вектор, определяющий интервал интегрирования [t tfinal]. Для получения решений в конкретные моменты времени t, tl. tfinal (расположенные в порядке уменьшения или увеличения) нужно использовать tspan = [t tl . tfinal];
· у0 — вектор начальных условий;
· Т, Y — матрица решений Y, где каждая строка соответствует времени, возвращенном в векторе-столбце Т.
Этапы решения ОДУ:
1 Приведение дифференциального уравнения к системе дифференциальных уравнений первого порядка. Для этого вводится столько дополнительных функций, каков порядок уравнения.
2 Написание специальной файл-функции для системы уравнений. Файл-функция содержит два входных аргумента: переменную t, по которой производится дифференцирование, даже если она входит в уравнение неявно, и вектор, размер которого равен числу неизвестных функций системы.
3 Вызов подходящего солвера (встроенной функции). Входными аргументами солвера, в простом случае, являются имя файл-функции в апострофах, вектор с начальным и конечным значениями переменной
Покажем применение решателя ОДУ на ставшем классическом примере — решении уравнения Ван-дер-Поля, записанного в виде системы из двух дифференциальных уравнений:
при начальных условиях
Перед решением нужно записать систему дифференциальных уравнений в виде ode-функции. Для этого в главном меню выберем File > New > M-File и введем
function dydt = vanderpol(t, y)
% функция правых частей дифф.уравнения
Сохраним m-файл-функцию. Тогда решение решателем ode45 и сопровождающий его график можно получить, используя следующие команды:
tspan = [0, 20]; % интервал интегрирования
y0 = [2; 0]; % начальные условия
[t,y] = ode45(@vanderpol, tspan, y0, []);
% Результат в виде графика
title(‘van der Pol Equation, mu = 1’)
Опора деревянной одностоечной и способы укрепление угловых опор: Опоры ВЛ — конструкции, предназначенные для поддерживания проводов на необходимой высоте над землей, водой.
Папиллярные узоры пальцев рук — маркер спортивных способностей: дерматоглифические признаки формируются на 3-5 месяце беременности, не изменяются в течение жизни.
Организация стока поверхностных вод: Наибольшее количество влаги на земном шаре испаряется с поверхности морей и океанов (88‰).
Поперечные профили набережных и береговой полосы: На городских территориях берегоукрепление проектируют с учетом технических и экономических требований, но особое значение придают эстетическим.