Решение нестандартных ОДУ

Эта страница содержит два примера решения нежестких обыкновенных дифференциальных уравнений с помощью ode45. MATLAB® имеет три решателя для нежестких ОДУ.

  • ode45

  • ode23

  • ode113

Для большинства нежестких проблем, ode45 выполняет лучше всего. Однако ode23 рекомендуется для проблем, которые разрешают немного более грубый ошибочный допуск или в присутствии умеренной жесткости. Аналогично, ode113 может быть более эффективным, чем ode45 для проблем со строгими ошибочными допусками.

Если нежесткие решатели занимают много времени, чтобы решить задачу или последовательно приводить интегрирование к сбою, то проблема может быть жесткой. Смотрите Решают Жесткие ОДУ для получения дополнительной информации.

Пример: Нежесткое уравнение Ван дер Поля

Уравнение Ван дер Поля является ОДУ второго порядка

$$y''_1 - \mu \left( 1 - y_1^2\right) y'_1+y_1=0,$$

где$\mu > 0$ скалярный параметр. Перепишите это уравнение как систему ОДУ первого порядка путем создания замены$y'_1 = y_2$. Получившаяся система ОДУ первого порядка

$$
\begin{array}{cl}
y'_1 &= y_2\\
y'_2 &= \mu (1-y_1^2) y_2 - y_1.\end{array}
$$

Система ОДУ должна быть закодирована в файл функции, который может использовать решатель ОДУ. Общая функциональная подпись функции ОДУ

  dydt = odefun(t,y)

Таким образом, функция должна принять оба t и y как входные параметры, даже если это не использует t для любых расчетов.

Файл функции vdp1.m кодирует использование уравнения Ван дер Поля$\mu = 1$. Переменные$y_1$ и$y_2$ представлены y(1) и y(2), и двухэлементный вектор-столбец dydt содержит выражения для$y'_1$ и$y'_2$.

function dydt = vdp1(t,y)
%VDP1  Evaluate the van der Pol ODEs for mu = 1
%
%   See also ODE113, ODE23, ODE45.

%   Jacek Kierzenka and Lawrence F. Shampine
%   Copyright 1984-2014 The MathWorks, Inc.

dydt = [y(2); (1-y(1)^2)*y(2)-y(1)];

Решите ОДУ с помощью ode45 функция на временном интервале [0 20] с начальными значениями [2 0]. Выход является вектор-столбцом моментов времени t и массив решения y. Каждая строка в y соответствует времени, возвращенному в соответствующей строке t. Первый столбец y соответствует$y_1$, и второй столбец к$y_2$.

[t,y] = ode45(@vdp1,[0 20],[2; 0]);

Постройте решения для$y_1$ и$y_2$ против t.

plot(t,y(:,1),'-o',t,y(:,2),'-o')
title('Solution of van der Pol Equation (\mu = 1) using ODE45');
xlabel('Time t');
ylabel('Solution y');
legend('y_1','y_2')

vdpode функция решает ту же задачу, но это принимает пользовательское заданное значение для$\mu$. Уравнения Ван дер Поля становятся жесткими как$\mu$ увеличения. Например, со значением$\mu = 1000$ необходимо использовать жесткий решатель, такой как ode15s решить систему.

Пример: Нежесткие эйлеровы уравнения

Уравнения Эйлера для твердого тела без внешних сил являются стандартной тестовой задачей для решателей ОДУ, предназначенных для нежестких проблем.

Уравнения

$$ \begin{array}{cl} y'_1 &= y_2y_3 \\ y'_2 &= -y_1y_3 \\ y'_3 &=
-0.51y_1y_2. \end{array}$$

Файл функции rigidode задает и решает эту систему уравнений первого порядка по временному интервалу [0 12], использование вектора начальных условий [0; 1; 1] соответствие начальным значениям$y_1$$y_2$, и$y_3$. Локальная функция f(t,y) кодирует систему уравнений.

rigidode вызовы ode45 без выходных аргументов, таким образом, решатель использует выходную функцию по умолчанию odeplot автоматически построить точки решения после каждого шага.

function rigidode
%RIGIDODE  Euler equations of a rigid body without external forces.
%   A standard test problem for non-stiff solvers proposed by Krogh.  The
%   analytical solutions are Jacobian elliptic functions, accessible in
%   MATLAB.  The interval here is about 1.5 periods; it is that for which
%   solutions are plotted on p. 243 of Shampine and Gordon.
%
%   L. F. Shampine and M. K. Gordon, Computer Solution of Ordinary
%   Differential Equations, W.H. Freeman & Co., 1975.
%
%   See also ODE45, ODE23, ODE113, FUNCTION_HANDLE.

%   Mark W. Reichelt and Lawrence F. Shampine, 3-23-94, 4-19-94
%   Copyright 1984-2014 The MathWorks, Inc.

tspan = [0 12];
y0 = [0; 1; 1];

% solve the problem using ODE45
figure;
ode45(@f,tspan,y0);

% --------------------------------------------------------------------------

function dydt = f(t,y)
dydt = [    y(2)*y(3)
   -y(1)*y(3)
   -0.51*y(1)*y(2) ];


Решите нежесткие уравнения Эйлера путем вызова rigidode функция.

rigidode
title('Solution of Rigid Body w/o External Forces using ODE45')
legend('y_1','y_2','y_3','Location','Best')

Смотрите также

| |

Похожие темы