Тема 7. Решение дифференциальных уравнений и систем в MathCad
Краткие теоретические сведения
Для решения дифференциальных уравнений с начальными условиями система Mathcad имеет ряд встроенных функций:
rkfixed – функция для решения ОДУ и систем ОДУ методом Рунге–Кутта четвертого порядка с постоянным шагом;
Rkadapt – функция решения ОДУ и систем ОДУ методом Рунге–Кутта с переменным шагом;
Odesolve – функция, решающая ОДУ блочным методом.
Ниже приведено описание стандартной функции rkfixed с указанием параметров функции.
y – вектор начальных условий из k элементов ( k – количество уравнений в системе);
x1 и x2 – левая и правая границы интервала, на котором ищется решение ОДУ или системы ОДУ;
p – число точек внутри интервала (x1, x2), в которых ищется решение;
D – вектор, состоящий из k-элементов, который содержит первую производную искомой функции или первые производные искомых функций, если речь идет о решении системы.
Результатом работы функции является матрица из p +1 строк, первый столбец которой содержит точки, в которых получено решение, а остальные столбцы – сами решения.
На рисунке 2.7.1 приведены конкретные примеры решения различных дифференциальных уравнений и систем ОДУ в MathCAD .
Рисунок 2.7.1 – Примеры решения дифференциальных уравнений и систем
При решении дифференциального уравнения первого порядка нужно создать вектор начальных условий из одного элемента Y 1 , который затем используется при формировании вектора-функции правой части дифференциального уравнения. При обращении к функции rkfixed указывается имя вектора Y , границы интервала, на котором ищется решение уравнения, например, (0 ; 2), количество точек, в которых ищется решение – 100, вектор-функция, описывающая правую часть дифференциального уравнения – D . В результате получается матрица z , в первом столбце которой содержатся значения аргумента искомой функции, во втором – значения самой результирующей функции. При построении графика функции первый столбец полученной матрицы указывается как аргумент, второй столбец – как функция.
При решении системы дифференциальных уравнений нужно создать вектор начальных условий из двух элементов, например, вектор v , который затем используется при формировании вектора-функции правой части дифференциального уравнения. При обращении к функции rkfixed указывается имя вектора v , и границы интервала, на котором ищется решение уравнения, например, (0 ; 5), количество точек, в которых ищется решение – 100, вектор-функция, описывающая правую часть дифференциального уравнения – D . В результате получается матрица s , в первом столбце которой содержатся значения аргумента искомых функций, во втором и третьем столбцах – значения самих функций при соответствующем значении аргумента. При построении графика можно воспользоваться первым столбцом полученной матрицы как аргументом, а вторым и третьим столбцами – как функциями.
На рисунке 2.7.2 приведен пример решения дифференциального уравнения второго порядка с использованием функции rkfixed . Необходимо решить дифференциальное уравнение второго порядка с заданными начальными условиями вида:
Рисунок 2.7.2 – Пример решения дифференциальных уравнений второго порядка с помощью rkfixed
Для решения уравнения с помощью функции rkfixed нужно выполнить замену переменных и привести дифференциальное уравнение второго порядка к двум дифференциальным уравнениям первого порядка. Вид этих уравнений приведен ниже.
Документ формируется точно так же, как и при решении системы ОДУ.
На рисунке 2.7.2 показана возможность вычисления вектора второй производной найденной функции – вектора а, построены графики исходной функции, функций первой и второй производных.
Практическая часть темы 7
7.1 Решение дифференциальных уравнений первого порядка
Последовательность действий для р ешения дифференциального уравнения первого порядка такова:
q сформировать вектор начальных условий из одного элемента, присвоив начальное значение искомой функции переменной с индексом, например: 
или 
(в зависимости от значения переменной ORIGIN );
q определить вектор-функцию из одного элемента, которая содержит первую производную неизвестной функции:
· набрать имя функции с двумя параметрами: первый параметр – аргумент искомой функции (независимая переменная), второй – имя вектора, содержащего искомую функцию (можно использовать имя вектора начальных условий), например, D ( x , Y );
· набрать оператор «:=» и выражение для первой производной (выразить из дифференциального уравнения), в котором вместо имени искомой функции подставлен первый элемент вектора-параметра, например, для уравнения 
вектор-функция будет определятся следующим образом: 
( если ORIGIN = 0 , подставлять 
);
q присвоить некоторой переменной значение функции rkfixed , указав в скобках следующие параметры:
· первый – имя вектора начальных условий,
· второй – левая граница интервала, на котором ищется решение, в виде числовой константы,
· третий – правая граница интервала, на котором ищется решение, в виде числовой константы,
· четвертый – количество точек, в которых ищется решение,
· пятый – имя вектора-функции, описывающего первую производную, без параметров;
например: 
,
(в результате получится матрица Z , в первом столбце которой содержатся значения аргумента искомой функции, во втором – значения самой функции);
q вывести матрицу, содержащую решение ДУ с помощь оператора «=», например: Z = ;
q построить график найденной функции ( см. тему 5 ), указав в качестве аргумента по оси абсцисс столбец 
, а в качестве значения функции по оси ординат – столбец 
( если ORIGIN = 0 , набирать соответственно 
и ).
Пример 7.1 Найти численное решение дифференциального уравнения первого порядка на интервале от 0.2 до 5 в 1000 точках, при начальном условии y (0)=0.1.
Выполнить графическую интерпретацию результатов.
7.2 Решение систем дифференциальных уравнений
Последовательность действий для р ешения системы дифференциальных уравнений первого порядка такова (описана для значения ORIGIN =0 ):
q перейти в исходной системе уравнений к однотипным обозначениям функций и выразить первые производные,
например, систему 
можно преобразовать в 
;
q в документе MathCad сформировать вектор начальных условий, количество элементов которого равно количеству уравнений системы, присвоив его некоторой переменной (см. тему 2);
например, 
;
q определить вектор-функцию, которая содержит первые производные искомых функций:
· набрать имя функции с двумя параметрами: первый параметр – аргумент искомых функций (независимая переменная), второй – имя вектора, содержащего искомые функции (можно использовать имя вектора начальных условий), например, D ( t , V );
(Замечание: если независимая переменная явно не присутствует в системе, то в качестве ее имени можно выбрать любую переменную)
· набрать оператор «:=» и вставить шаблон вектора, количество элементов которого равно количеству уравнений системы (см. тему 2)
· набрать в качестве элементов вектора правые части системы уравнений, в которых искомые функции представлены соответствующими элементами вектора-параметра, например,
;
q присвоить некоторой переменной значение функции rkfixed , указав в скобках следующие параметры:
· первый – имя вектора начальных условий,
· второй – левая граница интервала, на котором ищется решение, в виде числовой константы,
· третий – правая граница интервала, на котором ищется решение, в виде числовой константы,
· четвертый – количество точек, в которых ищется решение,
· пятый – имя вектора-функции, описывающего первые производные, без параметров;
например: 
,
(в результате получится матрица Z , в первом столбце которой содержатся значения аргумента искомых функций, во втором – значения первой функции, в третьем – значения второй функции и т. д.);
q вывести матрицу, содержащую решение системы ДУ с помощь оператора «=», например: Z = ;
q построить графики найденных функций ( см. тему 5 ), указав в качестве аргумента по оси абсцисс первый столбец матрицы решений, например, , а в качестве значений функций по оси ординат – остальные столбцы матрицы через запятую, например, , 
и т. д.
Пример 7.2 Найти решение системы дифференциальных уравнений
на интервале от 0 до 0.5 в 1000 точках, при следующих начальных условиях: x (0)=0.1 и y (0)=1.
Выполнить графическую интерпретацию результатов.
Составить дифференциальное уравнение в маткаде
Mathcad для студентов
Mathcad для начинающих
Скачать программы бесплатно
Функции предназначенные для решения обыкновенных дифференциальных уравнений в Mathcad
Для решения обыкновенных дифференциальных уравнений и систем обыкновенных дифференциальных уравнений в Mathcad введен ряд функций. Рассмотрим их:
- odesolve(x,b,step) — используется для решения обыкновенного дифференциального уравнения, заданного как в виде задачи Коши, так и в виде краевой задачи. Начальные условия и дифференциальное уравнение должны быть определены в блоке given. Параметры функции: х –переменная, по которой производится интегрирование; b — конечное значение промежутка решения; step – величина шага численного метода (параметр необязательный).
- rkfixed(u,a,b,N,D) – реализует в Mathcad численное решение задачи Коши по методу Рунге – Кутта с фиксированным шагом. Имеет следующие преимущества перед odesolve(x,b,step): может быть использована в программных модулях и позволяет оперативно пересчитывать результаты при изменении параметров. Параметры функции: u-вектор начальных условий; a и b – граничные значения отрезка решения задачи; N – число интервалов разбиения отрезка [a,b]; D(x,y) –вектор-функция, содержащая правые части первых производных, записанные в символьном виде.
- Rkadapt(u, a,b, N, D) — возвращает матрицу в Mathcad, содержащую таблицу значений решения задачи Коши на интервале от a до b для уравнения или системы обыкновенных дифференциальных уравнений, вычисленную методом Рунге-Кутта с переменным шагом и начальными условиями в векторе u, D(x,y) –вектор функция, содержащая правые части первых производных, записанная в символьном виде, n — число шагов.
- Функция Rkadapt() вследствие автоматического подбора шага, как правило, дает более точный результат по сравнению с другими функциями в Mathcad.
Метод конечных разностей в Mathcad
В случае краевых задач для линейных дифференциальных уравнений в Mathcad применяются формулы для аппроксимации производных соответствующими конечно – разностными отношениями. Это позволяет свести решение дифференциальных уравнений к решению системы линейных уравнений. Результаты получают в дискретных i – ых точках интервала решения задачи. При этом отрезок [a,b] разбивается на n частей с шагом h =(b-a)/n. Для аппроксимации соответствующих производных в Mathcad используют следующие формулы:
Таким образом, сделав соответствующую замену, получаем систему линейных уравнений, решение которой средствами Mathcad не представляет сложностей. Решение задачи методом конечных разностей приведено на листинге
Задача Коши в Mathcad
Задача Коши в Mathcad для дифференциальных уравнений n-го порядка с одной неизвестной (обыкновенное дифференциальное уравнение — ОДУ) формулируется следующим образом. Найти решение дифференциального уравнения
в виде функции y=y(x), которая удовлетворяет заданным начальным условиям
где 
— заданные значение. Решение задачи Коши для обыкновенных дифференциальных уравнений второго и более высоких порядков можно свести к системе уравнений. Решение задачи Коши для ОДУ первого порядка в Mathcad с использованием различных функций приведено на листинге.
Наибольшее распространение для решения задачи Коши в Mathcad получил метод Рунге – Кутта. Суть метода состоит в последовательном отыскании искомого значения функции yi+1 по формуле
За h принимается достаточно малый шаг, с помощью которого весь интервал задачи Коши разбивается на дискретные точки, в которых и ищется решение. Погрешность результатов пропорциональна пятой степени шага (h5).
Геометрический смысл метода Рунге – Кутта состоит в следующем. Из очередной точки (xi,yi) выбирается направление (угол) 
, для которого tg()=f(xi,yi). На этом направлении вычисляется точка с координатами 
Затем из точки (xi,yi) выбирается направление (угол) 
, для которого
tg()=f
На этом направлении вычисляется точка с координатами Далее из точки (xi,yi) выбирается направление (угол) 
, для которого 
На этом направлении в Mathcad вычисляется точка с координатами 
После чего из точки (xi,yi) выбирается направление (угол) 
, для которого 
. Все четыре полученных направления усредняются в соответствии с формулой для расчета 
. На этом результирующем направлении и строится расчетная точка с координатами
Метод Рунге – Кутта благодаря высокой точности широко используется при численном решении дифференциальных уравнений и в частности в Mathcad. Существует несколько разновидностей данного метода, которые нашли свое отражение в рассмотренных выше функциях. На листинге можно не только сравнить результаты, полученные на основе различных функций, но и оценить эти результаты с позиций точности расчетов.
Путем сравнения результатов решения задачи, можно сделать вывод о точности решения задачи. Наиболее точный результат позволяет получить функция Rkadapt.
Краевые задачи в Mathcad
Краевые задачи в Mathcad отличаются от задачи Коши состоит тем, что в краевой задаче начальные условия задаются на концах интервала поиска решения. Для решения подобных задач в системе Mathcad используется метод пристрелки, который начальное условие в правой точке интервала преобразует в дополнительное начальное условие для левой точки интервала. После чего краевая задача трансформируется в задачу Коши, методы решения которой были рассмотрены в предыдущем разделе. Для реализации метода пристрелки в Mathcad существует функция sbval. Данная функция определяет недостающие условия в начальной точке для двухточечных краевых задач. Функция имеет следующий синтаксис sbval(z,a,b,D,load,score), где z – вектор приближений недостающих начальных условий на левой границе; a,b – левая и правая граница интервала решений; D(x,y) – вектор-функция, содержащая правые части первых производных, записанная в символьном виде; load(a,z) – вектор-функция, описывающая начальные условия на левой границе интервала; score(b,y) – вектор-функция для задания правых граничных условий. Пример решения краевой задачи приведен на листинге.
Дифференциальное уравнение в Mathcad
Функции предназначенные для решения обыкновенных дифференциальных уравнений
Для решения дифференциальных уравнений и систем дифференциальных уравнений в Mathсad введен ряд функций.
Решение задачи Коши
Задача Коши для дифференциальных уравнений n-го порядка с одной неизвестной.
Краевые задачи
Разница краевой задачи и задачи Коши состоит в том, где задается интервала поиска решения.
Метод конечных разностей
В Mathcad в краевых задачах для уравнений применяются формулы для аппроксимации производных соответствующими отношениями.
Составить дифференциальное уравнение в маткаде
Электронный курс по MathCAD
5.2 Решение дифференциальных уравнений и систем.(Задача Коши и граничные задачи).
Решение одиночного дифференциального уравнения.
Для численного решения одиночного дифференциального уравнения в MathCAD имеется функция Odesolve, с помощью которой может быть решена как задача Коши для обыкновенного дифференциального уравнения, так и граничная задача. Эта функция входит в состав блока решения и сявляется его заключительным ключевым словом.
Odesolve(x,b,[step]) — Возвращает функцию, которая является решением дифференциального уравнения. Используется в блоке с оператором Given.
x — переменная интегрирования, действительное число
b — конечная точка отрезка интегрирования
step — величина шага по переменной интегрирования (необязательный аргумент)
Замечания:
- Уравнение должно быть линейным относительно старшей производной.
- Число заданных начальных или граничных условий внутри блока должно быть равно порядку уравнения.
- При записи уравнения для обозначения производных функции используйте специальные кнопки с панели Math или ‘ (штрих) — [Ctrl+F7], для знака равенства = [Ctrl+=] (в том числе и для дополнительных условий).
- Конечная точка должна быть больше начальной.
- Не допускаются начальные и граничные условия смешанного типа (f ‘(a)+f(a)=5).
- Искомая функция в блоке дложна быть обязательно с аргументом ( f(x))
 |  |
Численное решение задачи Коши для дифференциальных уравнений и систем.
Для численного решения задачи Коши для дифференциальных уравнений и систем могут быть использованы функции:
rkfixed(y,x1,x2,n,F) — возвращает матрицу решений системы уравнений методом Рунге-Кутта 4-го порядка при фиксированном шаге по x
rkadapt(y,x1,x2,n,F) — ищет решение с переменным шагом ( там, где решение меняется медленнее, шаг увеличивается, а в области быстрого изменения решения шаг функции уменьшается). Возвращается решение с равным шагом. Функция работает быстрее, чем rkfixed
Bulstoer(y,x1,x2,n,F) — дает более точное решение (методом Bulirsch-Stoer)
Агрумкнты вышеуказанных функций:
y — вектор начальных условий
x1,x2 — границы интервала для поиска решения
n — количество точек на интервале
F(x,y) — вектор-функция первых производных
При решении дифференциальных уравнений порядка выше первого (или систем уравнений, выше первого порядка) исходное уравнение (систему) необходимо преобразовать к системе дифференциальных уравнений первого порядка.
В результате работы укзанных функций рассчитывается матрица, количество стобцов которой равно порядку уравнения +1(или сумме порядков уравнений в системе +1), а количество строк равно параметру n. Первый столбец содержит значения независимой переменной, второй — значение функции, третий — для диф. уравнений 2-го порядка — значение производной искомой функции (если решается система двух уравнений 1-го порядка, то третий столбец будет содержать значения второй функции). Для выделения решений (функций или их производных) можно воспользоваться стандартным оператором вывода столбцов матрицы M < >
Если матрица правых частей дифференциальных уравнений почти вырождена, то такие системы называются жесткими. В этом случае решения, возвращаемые функцией rkfixed будет неустойчивым и для решения таких систем необходимо применять функции Stiffb , Stiffr
Stiffb(y,x1,x2,n,F,J) — ищет решение диф. уравнения или системы дифференциальных уравнений методом Bulirsch-Stoer
Stiffr(y,x1,x2,n,F,J) — ищет решение диф. уравнения или системы дифференциальных уравнений методом Rosenbrock
 |
 |
 |
