Смекни!
smekni.com

Решение дифференциальных уравнений в среде MathCAD (стр. 7 из 8)

Если не заданы начальные условия, тодифференциальное уравнение n – го порядкаимеет бесконечное множество решений, при задании начальных условий y(x0)= y0, y’(x0)= y0,1, y’’(x0)= y0,2, …, y(n-1)(x0)= y0,n-1 решение становится единственным (задача Коши).

Задача Коши для дифференциальногоуравнения n – го порядка может быть сведена к задаче Коши длянормальной системы n дифференциальных уравнений 1го порядка, которая в векторной форме имеет вид

Y’ = F(x,Y), Y(x0) = Y0

Y(x0) = Y0 – вектор начальных условий;

  Y’=(y’1, y’2, …, y’n) – вектор первых производных;

F(x, Y)= (y2, y3, …, yn, f(x,y1, … , yn) – вектор правых частей;

Y = (y2, y3, …, yn) – вектор искомого решения.

Эта система получается в результатеследующей замены:

   
,где
 

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

Вычислительный блок Given/Odesolve

Ниже приведены два примера для решениядифференциальных уравнений первого и второго порядка с использованиемвычислительного блока решения Given/Odesolve.


Вычислительный блок для решения одногоОДУ состоит из трех частей:

-   ключевое слово given;

-   ОДУ и начальные условия, записанные с помощьюлогического равенства;

-   встроенная функция Odesolve(x,b) относительно независимой переменной xна интервале [a, b]; b –верхняя граница отрезкаинтегрирования. Допустимо и даже предпочтительнее задание функции Odesolve(a, b, step) с тремя параметрами, где step –внутренний параметр численного метода, определяющий количество шагов; чембольше step, тем с лучшей точностью будет получен результат, нотем больше времени будет затрачено на его поиск.

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

Функция Odesolveиспользует для решения дифференциальных уравнений наиболее популярный алгоритмРунге-Кутта четвертого порядка, описанный в большинстве книг по методамвычислений. Он обеспечивает малую погрешность для широкого класса систем ОДУ заисключением жестких систем. Если щелчком правой кнопки мыши на блоке формул сфункцией Odesolve вызвать контекстное меню, то можно изменить методвычисления решения, выбрав один из трех вариантов: Fixed –метод Рунге-Кутта с фиксированным шагом интегрирования (этот метод используетсяпо умолчанию), Adaptive – также метод Рунге-Кутта, но с переменным шагом,изменяемым в зависимости от скорости изменения функции решения, Stiff – метод, адаптированный для решения жестких уравнений и систем(используется так называемый метод PADAUS).

Альтернативный метод решения ОДУ заключается в использовании одной извстроенных функций: rkfixed, Rkadapt, или Bulstoer. Все они решают задачуКоши для системы дифференциальных уравнений первого порядка, нокаждая из них использует для этого свой метод. Для простых систем неиграет большой роли, какой метод использовать – все равно получите решениедостаточно быстро и с высокой точностью. Но для сложных или специфическихсистем бывает, что некоторые методы вообще не могут датьудовлетворительного решения за приемлемое время. Именно для таких сложных, ноне редких случаев в MathCAD и введено несколько различных методов решениясистем ДУ.

-   rkfixed – метод Рунге-Кутта с фиксированным шагом интегрирования. Самыйпростой и быстрый метод, но далеко не всегда самый точный. Полностью аналогичениспользованию функции Odesolve с выбранным в контекстном меню методом Fixed.

-   Rkadapt – метод Рунге-Кутта с переменным шагом интегрирования. Величина шагаадаптируется к скорости изменения функции решения. Данный метод позволяетэффективно находить решения уравнений, в случае если оно содержит как плавные,так и быстро меняющиеся участки. Там, где решение меняется слабо, шагивыбираются более редкими, а в областях его сильных изменений – частыми. Врезультате для достижения одинаковой точности требуется меньшее число шагов,чем для rkfixed. Полностью аналогичен использованию функции Odesolve с выбранным в контекстном меню методом Adaptive.

-   Bulstoer – метод Булирша – Штера. Этот метод более эффективен, чем методРунге-Кутта, в случае если решение является плавной функцией.

Имена функций Rkadapt и Bulstoer начинаются с прописной буквы. В MathCAD длянекоторых имен функций неважно, с какой буквы они записаны, но дляперечисленных функций это принципиально, т.к. в MathCADтакже существуют функции с такими же именами, только записанные с маленькойбуквы – rkadap, bulstoer. Эти функции используются в тех случаях, когдаважным является решение задачи в конечной точке интервала интегрирования.


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

Применение встроенных функций вдокументах MathCAD выглядит сходным образом, т.е. функции Rkadapt и Bulstoer имеют тот же синтаксис, что и выше приведеннаяфункция rkfixed. Назначение аргументов в этих встроенных функцияхследующее:

-    y – вектор начальных значений неизвестных функций,входящих в систему. В случае одного уравнения и одной неизвестной функции – этопросто число.

-    а – начало отрезка, на которомищется решение системы (отрезка интегрирования). Именно в этой точке значениянеизвестных функций принимаются равными элементам вектора y.

-    b – конец отрезка интегрирования.

-    n – количество частей, на которые разбивается отрезок[a, b] при решении системы. Чем больше это число, темточнее получается решение, но расчет занимает больше времени.

-    F(x,y) – векторная функция, элементы которой содержатправые части уравнений системы в нормальной форме (когда левые части –первые производные от соответствующих функций, а в правых частях производныеотсутствуют). Аргументами этой функции являются вектор y, элементыкоторого соответствуют различным неизвестным функциям системы, и скалярныйаргумент x , соответствующий независимой переменной в системе.В случае одного уравнения функция F может быть скалярнойфункцией, зависящей от двух скалярных переменных x и y.

Возвращаемым значением всехвышеперечисленных встроенных функций является матрица. Первый столбец этойматрицы – это точки, на которые разбивается отрезок [a, b],а остальные столбцы – это значения функций системы в этих точках. Если варгументе функции rkfixed было указано количество частей n= 100, то матрица будет содержать 101 строку вместе с начальной.

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

Для численного интегрирования систем ОДУв 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 наглядность теряется, и для визуализациифазового портрета приходится строить его различные проекции.