Решение нежестких ОДУ

Эта страница содержит два примера решения нежестких обыкновенных дифференциальных уравнений с помощью 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}
$$

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

  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')

The 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')

См. также

| |

Похожие темы