|
Набор парамтров для этих функций :
bulstoer(y,x1,x2,acc,F,kmax,save)
rkadapt(y,x1,x2,acc,F,kmax,save)
stiffb(y,x1,x2,acc,F,J,kmax,save)
stiffr(y,x1,x2,acc,F,J,kmax,save)
Первые три параметра и пятый (F) этих функций те же, что идля функции Rkadapt. Дополнительные параметры:
acc — параметр, контролирующий точность решения (реком. асс=0.001)
kmax — максимальное число промежуточных точек в которых ищется решение
save — минимально допустимый интервал между точками, в которых ищется решение
Решение граничных задач для обыкновенных дифференциальных уравнений.
Если для дифференциального уравнения n-го порядка k граничных условий заданы в начальной точке х1, а (n-k) граничных условий — в конечной точке х2, то такая задача называется краевой. В MathCAD реализованы две функции, позволяющие численно найти недостающие условия в точках х1 и х2.
Двухточечная краевая задача
Задача решается в два этапа. Сначала с помощью функции sbval находятся недостающие начальные значения, а затем применяется одна из выше описанных функций для решения стандартной задачи Коши на отрезке.
sbval(v,x1,x2,F,load,score) — ищет недостающие начальные условия в точке х1
v — вектор началных приближений для искомых начальных значений в точке х1,
х1,х2 — граничные точки интервала
F(x,y) — вектор-столбец из n элементов, содержит правые части дифференциальных уравнений
load(x1,v) — вектор-столбец из n элементов, содержит начальные значения в точке х1; некоторые из значений- константы, другие неизвестны и будут найдены в процессе решения.
score(x2,y) — вектор-столбец размерности вектора v, содержащий разность между начальным условием в точке х2 и значеием искомого решения в этой точке.
Краевая задача с условиями внутри интервала.
На первом этапе используется функция
balfit(V1,V2,x1,x2,xf,F,load1,load2,score) — ищет недостающие начальные условия в точках х1 и х2, сшивая решения, выходящие из этих точек, в точке xf
V1,V2 — вектора началных приближений для искомых начальных значений в точках х1 и х2
х1,х2 — граничные точки интервала
load1(x1,V1) — вектор-столбец из n элементов, содержит начальные значения в точке х1; некоторые из значений- константы, другие неизвестны и будут найдены в процессе решения
load2(x2,V2) — вектор-столбец из n элементов, содержит начальные значения в точке х2; некоторые из значений- константы, другие неизвестны и будут найдены в процессе решения.
score(xf,y) — вектор-столбец размерности n, содержащий разность между решениями, начинающимися в точках х1 и х2, в точке xf
БлогNot. MathCAD: решаем основные типы дифференциальных уравнений встроенными функциями
MathCAD: решаем основные типы дифференциальных уравнений встроенными функциями
Решать дифференциальные уравнения (далее ДУ) в MathCAD, составляя собственные подпрограммы-функции, не всегда удобно и экономично по времени, хотя и полезно на этапе обучения. Опишем в этой заметке способы решения основных типов ДУ с помощью стандартных средств пакета, ограничимся простыми примерами.
1. ДУ с разделяющимися переменными. Общая постановка задачи: y’=f(x,y)=g(x)*h(y) , y(x0)=y0 . То есть, f(x,y) допускает представление в виде произведения функций от x и от y .
Для решения уравнения достаточно задать его правую часть как пользовательскую функцию MathCAD, определить интервал поиска решения [x0,x1] , начальное условие y0 и применить стандартную функцию Odesolve . Покажем этот процесс на примере уравнения y’=2x-y+x 2 , x∈[0,2] , y(0)=0 с известным решением y(x)=x 2 :
Знак "равно" в записи уравнений, конечно же, жирный (панель Boolean или сочетание клавиш Ctrl+=).
Функция Odesolve вернула именно функцию y , её нужно смотреть от аргумента, например, y(1)= .
И ещё 2 особенности:
Фича MathCAD — у ряда функций решения ДУ есть "версии", которые пишутся и с Большой, и с маленькой буквы. Функции, имена которых начинаются с маленькой буквы, используются в тех случаях, когда важным является решение задачи в конечной точке интервала.
Фича MathCAD — в блоке Given-Odesolve штрих ` является не просто знаком, а оператором взятия производной. И это не штрих ‘ с клавиши русской "Э". это обратный штрих ` с клавиши "Ё", слева от "1".
Зная точное решение, графически сравним с ним найденное решение. Как видно на графике, MathCAD справился с задачей отлично.
2. Неоднородное ДУ первого порядка. В общем виде такое уравнение можно записать как y’=a(x)*y+b(x) . Оно решается аналитически по формуле, которую можно найти в любой книге по решению обыкновенных ДУ:
Здесь С – константа интегрирования. Остаётся применить формулу к конкретному уравнению (возьмём для примера задачу y’+2xy=x*e -x 2 sin(x) , y(0)=1 ) и оценить её символьно:
Здесь мы получаем решение в общем виде. Нижний оператор оценён символьно (см. панель "Символика"), а аргумент t используется, так как x в документе "уже занят" (для корректной работы символьной оценки переменные не должны быть определены заранее).
После подстановки начального условия получим частное решение y(t,Y) , а для проверки решеня будет достаточно подставить полученную функцию в исходное уравнение и упростить его символьной функцией simplify . Полученный результат в нашем случае совпал с заданной в условии правой частью. Также для y(x,Y) , как и для любой функции, можно построить график на нужном интервале изменения x.
3. Неоднородное ДУ второго порядка. В общем виде имеем уравнение y» + p(x)*y’ + g(x)*y = f(x) плюс набор краевых условий, количество которых соответствует порядку задачи, например y(0)=. y'(0)=. или y(0)=. y(1)=.
Возьмём уравнение, которое мы мучили вот здесь, решим его стандартными средствами, сравним с известным точным решением и построим график:
Здесь при вызове Odesolve второй параметр, равный единице — это правая граница интервала, третий параметр, равный 10, задаёт количество интервалов. Точное решение u(t) взяли по ссылке. Как видим, даже на 10 интервалах всё очень хорошо совпадает.
4. Система ДУ. Подход к решению системы ДУ покажем на примере. Пусть задана система дифференциальных уравнений
x’ = a*x — y — (x 2 + y 2 )*x,
y’ = a*y + x — (x 2 + y 2 )*y,
x(0)=0, y(0)=1, a=-0.2
Чтобы решить эту систему стандартной функцией rkfixed , нужно задать для неё вектор начальных значений x = (x, y) и вектор правых частей D(t,x) .
После этого задача решится вызовом rkfixed , второй и третий параметры ( 0, 20 ) задают интервал по времени t , на котором ищется решение, четвёртый параметр 100 означает количество точек на интервале.
Функция вернёт матрицу решений системы, в которой количество строк соответствует количеству точек на интервале, а количество столбцов — количеству уравнений в системе.
Для построения графика достаточно отобразить зависимость столбцов Zi,1 , Zi,2 от Zi,0 , i=0..99 :
Скачать расчёты из этой статьи в архиве .zip с документом .xmcd, Mathcad 15 (46 Кб)
10.11.2015, 17:16; рейтинг: 22557
Глава 5. Решение дифференциальных уравнений
5 .1 Вычислительный блок Given– Odesolve
Применение функции Odesolve требует записи вычислительного блока, состоящего из трех частей:
1) ключевого слова Given (Дано);
2) дифференциального уравнения и начальных или граничных условий к нему;
3) функции о desolve ( x , xk , n ) (решение ОДУ), где x – имя переменной, относительно которой решается уравнение; xk – конец интервала интегрирования (начало интервала интегрирования указано ранее, в начальных условиях); n – необязательный внутренний параметр, определяющий число шагов интегрирования, на которых решается дифференциальное уравнение.
Примеры использования функции Odesolve приведены на рис. 5.2–5.5.
Given
граничные условия можно задавать лишь в двух точках:
в начале и в конце интервала интегрирования
Рис. 5.2 Решение уравнения с граничными условиями
Если вид дифференциального уравнения изменяется с изменением аргумента, например, уравнение содержит ступенчатую разрывную функцию, то такое уравнение можно записать с использованием условия, то есть записать несколько уравнений в виде одного уравнения (рис. 5 . 3 ).
MathCAD позволяет записать условие несколькими способами:
– с использованием булевых операторов:
;
– с использованием функции if :
.
В данном случае невозможно записать условие с использование условного оператора программирования if , так как вертикальной черте (обязательному признаку программирования) должен предшествовать оператор присваивания :=, тогда как в нашем случае использован знак + или знак логического равенства =.
поменяйте местами закрашенные выражения
Given
Рис. 5.3 Объединение двух дифференциальных уравнений в одно
В MathCAD все встроенные функции для решения дифференциальных уравнений требуют записи производной высшего порядка в явном виде. Однако, как показано на рис. 5 . 2 , присутствие коэффициента (постоянного или переменного) перед производной высшего порядка – не помеха при решении ОДУ. Если уравнение есть сложная функция от высшей производной, предварительно надо решить это уравнение алгебраически относительно высшей производной. На рис. 5 .4. исходное уравнение решено символьно относительно у'(х). Для выполнения этой операции надо выделить переменную, затем в главном меню выбрать команду Symbolics → Variable → Solve (Символьные вычисления→Переменная→Решить). Найденный результат подставлен в дифференциальное уравнение, которое решено с помощью функции Odesolve.
В MathCAD 2001 i возможности функции Odesolve были расширены. Теперь она может решать и системы дифференциальных уравнений. При этом несколько изменяется обращение к ней.
Дано уравнение
решим его относительно у'(х)
has solution(s)
уравнение в стандартном виде подставим в блок Given — Odesolve
два корня уравнения дают решение ОДУ
Given
Рис. 5.4 Приведение дифференциального уравнения к стандртному виду и его решение
Обращение к функции для решения одного уравнения выглядит так:
о desolve ( x , xk , n ),
для решения системы дифференциальных уравнений–
о desolve ((вектор имен неизвестных), x , xk , n ).
Н а рис. 5 .5 приведено решение системы из трех дифференциальных уравнений первого порядка. соответственно, вектор имен неизвестных содержит имена трех неизвестных, х, у и z . Граничные условия заданы в точке t = –1. конечная точка интервала интегрирования t =3. Следовательно, интервал интегрирования от –1 до 3. Конечное значение всегда должно быть больше начального.
Если график строится для значений, выходящих за пределы интервала интегрирования, MathCAD экстраполирует значения функции за пределы интервала интегрирования, получая чаще всего неправдоподобные результаты.
Если в MathCAD 11 функция Odesolve решает системы уравнений только с начальными параметрами, то в MathCAD 12, 13, 14 – и с начальными, и с граничными условиями.
поменяйте вид уравнений и граничные условия
Рис. 5.5 Решение системы ОДУ функцией Odesolve
В MathCAD 2001 i у функции Odesolve появилась возможность принимать в качестве ограничений алгебраические уравнения (рис. 5.6 и 5.7), но только в задачах с начальными условиями.
Рассмотрим механизм работы блока Odesolve: блок Odesolve преобразует заданные дифференциальные уравнения к стандартному виду, используемому функциями rkfixed и Rkadapt . С помощью указанных функций он решает систему ОДУ (или одно ОДУ) с начальными условиями. В результате расчета получают массив значений функции и ее производные.
Далее используются функции интерполяции lspline и interp , преобразующие массив решений системы ОДУ в функцию, которую затем можно дифференцировать или интегрировать. Необязательный параметр steps в обращении к функции Odesolve как раз и задает число точек интерполяции. Чем больше steps , тем выше точность интерполяции, но тем больше время решения задачи.
При решении системы ОДУ есть возможность выбрать метод решения. Для этого надо установить курсор на слове Odesolve и нажать правую кнопку мыши. В открывшемся контекстном меню выбрать Fixed (Решение с фиксированным шагом функцией rkfixed ), Adaptive (Решение с переменным шагом функцией Rkadapt ) или Stiff (Решение жесткой системы ОДУ функцией Radau ) , как показано на рис. 5.7 .
дифференциальные
уравнения
алгебраическое уравнение
требуется 6 начальных условий:
Рис. 5. 6 Решение системы ОДУ с алгебраическими ограничениями функцией Odesolve
Given
Рис. 5.7 В ыбор метода решения дифференциальных уравнений