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

Для численного интегрирования систем ОДУ в MathCAD также имеется выбор – либо использовать вычислительный блок Given/Odesolve, либо встроенные функции rkfixed, Rkadapt и Bulstoer.

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

где Y и Y’ – соответствующие неизвестные векторные функции переменной t, а F – вектор правых частей системы уравнений первого порядка. Именно векторное представление используется для ввода системы ОДУ в среде MathCAD.

Если в систему ОДУ входят и уравнения высших порядков, то оно тоже сводится к системе уравнений первого порядка, как было показано выше. При этом количество нулевых условий для вычислительного блока Given/Odesolve, а также размер вектора начальных условий y и размер вектора правых частей F(x,y) для встроенных функций rkfixed, Rkadapt и Bulstoer должны быть равны сумме порядков всех уравнений.

Вначале покажем решение систем ОДУ первого порядка с использованием вычислительного блока Given/Odesolve

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

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

Рассмотрим решение этой же системы ОДУ первого порядка с использованием встроенной функции rkfixed.

Полученное решение полностью соответствует вышеприведенному решению с использованием вычислительного блока Given/Odesolve. Следует отметить, что начальные условия здесь задаются в виде вектора y, а функциям x(t) и y(t) соответствуют элементы этого вектора y₁ и y₂. Вектор начальных условий y и вектор правых частей F имеют размер равный двум, т.к. система состоит из двух уравнений первого порядка. Для системы ОДУ, состоящей из двух уравнений второго порядка, размер этих векторов будет равен четырем

Решение жестких ОДУ и систем ОДУ.

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

Для жестких систем (stiff) не работает обычный метод Рунге-Кутта или Булирша-Штера. Наличие быстро затухающего возмущения приводит к тому, что эти численные методы дают расходящееся решение. Для жестких задач разрабатываются специальные методы. В MathCAD предусмотрены три различные функции для решения жестких задач:

- Radau – метод Radaus для жестких систем. Полностью аналогичен использованию функции odesolve с выбранным в контекстном меню методом Stiff.

- Stiffb – метод Булирша-Штера, адаптированный для жестких систем.

- Stiffr – метод Розенброка.