Решение дифференциальных уравнений и систем в прикладной программе Mathcad.
Mathcad имеет ряд встроенных функций для решения обыкновенных дифференциальных уравнений. В результате решения получается матрица, содержащая значения функции, вычисленные на некотором множестве точек (на некоторой сетке значений).
Каждая из встроенных функций, предназначенных для решения ОДУ, требует, чтобы было задано следующее:
начальные условия;
множество точек, в которых необходимо найти решение;
cамо дифференциальное уравнение;
Наиболее употребляемой для решения ОДУ является функция rkfixed, которая использует для поиска решения метод Рунге-Кутта четвертого порядка. В результате решения получается матрица, имеющая два столбца. Первый столбец – точки, в которых ищется решение ОДУ. Второй – это столбец значений найденного решения в соответствующих точках. Функция rkfixed (y, x1, x2, npoints, D) имеет следующие аргументы:
у – вектор начальных условий размерности n, где n – порядок дифференциального уравнения или число уравнений в системе. Для уравнения первого порядка этот вектор вырождается в точку.
x1, х2 – граничные точки интервала, на котором ищется решение.
Npoints – число точек (не считая начальной), в которых ищется приближенное решение. Число строк в возвращаемой матрице решения определяется как (1+npoints).
D(x, y) – функция, которая возвращает значение в виде вектора из n элементов, содержащих первые производные неизвестных функций.
Дифференциальные уравнения первого порядка.
Рассмотрим на примере, как решаются дифференциальные уравнения первого порядка.
Пример. Решить уравнениеу¢ = 2у/x + х с начальным условиему(1) = 0 на отрезке[1; 1,5]
Для решения дифференциального уравнения первого порядка
У0: = 0 – задание начального значения
D(x,
y):
=
– определить функцию, задающую
производную
z:= rkfixed (y,1,1.5,50,D). – находится решение в 50 точках отрезка [1, 1,5]
Дифференциальные уравнения второго порядка.
Для решения дифференциального уравнения второго порядка необходимо задать:
вектор начальных условий состоит из двух элементов: значений функции и ее первой производной в начальной точке интеграла.
функция D(х, у) – вектор с двумя элементами:
D(x,
y) =
Матрица решения содержит три столбца: первый значения х, в которых ищется решение, второй – у(х) и третий – у¢(x)
Пример.
Решить
дифференциальное уравнение второго
порядка
со следующими начальными условиямиу(0) = 1, у¢(0)
= 3 на
отрезке [0;
0,5].
Системы дифференциальных уравнений.
Для решения системы обыкновенных дифференциальных уравнений первого порядка необходимо:
определить вектор, содержащий начальные значения для каждой неизвестной функции;
определить функцию, которая возвращает значение в виде вектора из n элементов, содержащих первые производные каждой из неизвестных функций;
выбрать точки для поиска приближенного решения.
задать функцию rkfixed.
Функция вернет матрицу, первый столбец которой – значения аргумента, а остальные столбцы содержат значения найденных приближенных решений в соответствующих точках.
Пример.
Решить систему двух ОДУ
с начальным условиями х0(0) = 0, хi(0) = 1
Примечание. Mathcad имеет ещё две встроенные функции для решения задачи Коши. Это функции Rkadapt и Bulstoer. Эти функции имеют те же самые аргументы и возвращают решения в такой же форме, что и функция rkfixed. Первая из этих функций использует метод Рунге-Кутты с переменным шагом, что позволяет повысить точность вычислений и сократить их объём, если искомое решение имеет области, где её значения меняются быстро и области плавного изменения. Функция Rkadapt будет варьировать величину шага в зависимости от скорости изменения решения.
Функция Bulstoer реализует иной численный метод – Метод Булирша-Штёра. Её следует применять, если известно, что решение является гладкой функцией.
Решение краевой задачи для дифференциального уравнения.
Пусть имеется дифференциальное уравнение второго порядка
.
Решение
уравнения требуется найти на интервале
[0; 1], причём известно, что у(0)=у0,
у(1)=у1.
Понятно, что произвольный интервал [a,
b]
заменой переменных
может быть сведен к единичному.
Для решения краевой задачи обычно применяют метод стрельб, позволяющий свести решение краевой задачи к задаче Коши. Пусть k некоторый параметр, такой что y(0)=k. Для некоторого пробного значения k может быть решена задача Коши, например, методом Рунге-Кутты. Мы хотим найти такое значение параметра, чтобы выполнялось условие y(1,k)=y1. Фактически мы свели исходную задачу к задаче решения трансцендентного уравнения с таблично заданной функцией. Если найдены такие значения параметра k1 и k2, что y(1,k1)>y1 и y(1,k2)<y1, то дальнейшее уточнение значения параметра можно проводить методом деления отрезка пополам.
Пусть
требуется найти решение дифференциального
уравнения
при
условиях
и
.
При
значениях параметров 
.
Для решения краевой задачи имеется встроенная функция sbval, реализующая метод стрельб и позволяющая свести краевую задачу к задаче Коши. Функция sbval имеет следующие параметры:
v– вектор, содержащий начальные приближения для недостающих начальных условий,
xmin, xmax – границы интервала, на котором ищется решение,
D(x,y) – вектор–функция, содержащий правые части системы дифференциальных уравнений первого порядка, эквивалентной исходному уравнению, размер вектора n совпадает со степенью старшей производной дифференциального уравнения
load(xmin,v) – вектор–функция, элементы которой соответствуют n значениям функций на левой границе интервала. Часть этих значений известна, а для части заданы начальные приближения в векторе v. Их уточненные значения будут найдены в процессе вычисления
score(xmax,y) – вектор–функция, имеющая то же число элементов, что и v. Каждое значение является разностью между начальными значениями в конечной точке интервала и соответствующей оценки для решения. Этот вектор показывает, на сколько близко найденное решение к истинному.
Наша задача сводится к системе двух дифференциальных уравнений первого порядка:
Поэтому
функция D имеет вид 
Задаем
граничные условия: 
Задача
Коши для дифференциального уравнения
второго порядка содержит два начальных
условия. Нам известно только одно.
Начальное приближение для недостающего
значения задаем в векторе v,
который в нашем случае состоит только
из одного элемента. Несмотря на это,
индекс о должен быть обязательно указан,
чтобы подчеркнуть векторный характер
этой величины: 
На
левой границе интервала нам известно
значение
и
задано начальное приближение для
.
Это значение записано в
.
Задаем вектор-функциюload.
Ее нулевой элемент – начальное
значение для
,
первый – для
.
Теперь, когда нам стало известно недостающее начальное условие в задаче Коши, можно воспользоваться, например, функцией rkfixed
