Самым простым случаем является дифференциальное уравнение 1 - го порядка, разрешенное относительно производной. Оно имеет вид y' = f(t,y). Если переменная y является вектором, то речь идет о системе дифференциальных уравнений. Дифференциальные уравнения и системы порядка 2 и выше обычно можно свести к системам уравнений 1 - го порядка, введя такое количество вспомогательных функций, которое соответствует порядку уравнения. Например, чтобы свести к системе уравнений 1 - го порядка уравнение y'' = F(t,y,y'), нужно использовать подстановку: y1 = y, y2 = y'. В результате получим систему дифференциальных уравнений
Для численного решения обыкновенных дифференциальных уравнений произвольного порядка и систем с начальными условиями, т. е. задачи Коши, в MATLAB имеютя необходимые решатели (солверы) – ode23, ode45, ode113, ode15s, ode23s, ode23t, ode23tb. Справку по ним можно получить с помощью команды doc ode45. Солверы с суффиксом s предназначены для решения так называемых жестких систем. К ним относится уравнение Ван-дер-Поля y'' - μ(1 - y2)y'+y = 0, решение которого с помощью ode15s приводится во встроенной в MATLAB справочной системе Help (приложение 1, рис. П.4).
Для всех остальных систем наиболее употребительным является ode45, реализующий алгоритм Рунге-Кутта 4 - 5 порядков (разные порядки точности используются для контроля шага интегрирования).
Все вышеупомянутые солверы позволяют задать дополнительный параметр options, контролирующий вычислительный процесс. Пример контроля точности вычислений рассмотрен в разделе 6.4.
В самом простом случае синтаксис вызова перечисленных выше солверов
[t,Y]=solver('fun', [t0, tn], Y0).
Здесь приняты следующие обозначения:
solver – название соответствующего солвера;
'fun' – имя файл-функции, которая вычисляет вектор-столбец правых частей системы дифференциальных уравнений;
t – вектор-столбец, содержащий значения независимой переменной (обычно это значения времени);
Y – матрица значений неизвестных функций в соответствующие моменты времени;
[t0, tn] – вектор с начальным и конечным временем наблюдения;
Y0 – скаляр или вектор-столбец, в котором задаются начальные условия.
Схема решения состоит из следующих этапов:
приведение дифференциального уравнения к системе дифференциальных уравнений первого порядка;
оформление в виде файл - функции правой части системы уравнений;
вызов подходящего солвера;
визуализация результата.
Пример:
Решить задачу Коши для дифференциального уравнения 2-го порядка
y''+4y'+3y = cost, y(0) = 1, y'(0) = 0.
Решение:
Выполнив подстановку y1 = y, y2 = y', получим систему дифференциальных уравнений
Теперь составим файл-функцию для вычисления правых частей системы дифференциальных уравнений. Она должна содержать два входных аргумента (переменную t, по которой производится дифференцирование и вектор, размер которого равен числу неизвестных функций системы) и один выходной аргумент (вектор правой части системы). Для нашего примера текст файл-функции будет таким:
function F=osc(t,y)
F=[y(2);-4*y(2)-3*y(1)+cos(t)];
Для вычисления решения системы на интервале [0;15] используем солвер ode45. С учетом начальных условий обращение к функции ode45 будет иметь такой вид:
>> [T,Y]=ode45('osc',[0,15],[1;0]);
В силу проделанных замен y1 = y, y2 = y', первый столбец матрицы Y содержит как раз значения неизвестной функции y(t), входящей в исходное дифференциальное уравнение, а второй столбец – значения ее производной.
Отобразим на рис. 6.5 график решения и график производной от этого решения (т. е. графики изменения координаты точки и ее скорости в зависимости от времени) с помощью команды
>> plot(T,Y(:,1),T,Y(:,2),'k--')
В разделе 7.14 найдено точное (аналитическое) решение
y = cost+sint - e -3t+e -t
данного дифференциального уравнения. Выведем на рис. 6.5 график точного решения в с маркерами виде кружков:
Как видно из рис. 6.5, график приближенного решения совпадает с графиком точного решения.
Приведенный здесь пример касается лишь задачи Коши для систем обыкновенных дифференциальных уравнений. Однако система MATLAB предлагает солверы для граничных задач, а также солверы для дифференциальных уравнений с частными производными.
Об этих солверах можно узнать во встроенной в пакет MATLAB справочной системе Help (приложение 1, рис. П.5).