matlab/ Octave 随机事件ode45

Question

我在理解如何在 Octave/Matlab 中实现事件以及求解微分方程时遇到了一些困难。

例如，考虑这个简单的代码来求解微分方程 y' = -y:

function dy = odefun(t,y)
    dy = -y;
endfunction

options = odeset('RelTol',1e-4,'AbsTol',[1e-4]);
[T,Y] = ode45(@odefun,[0 12],[1],options);

现在我想介绍一个随机事件。例如，在某个固定的时间步长，我想随机改变 y 的值，然后根据微分方程继续演化。我怎样才能做到这一点？

Answer 1

做到这一点的方法是分段集成系统并将每次运行的结果输出附加在一起:

% t = 0 to t = 6
tspan = [0 6];
y0 = 1;
options = odeset('RelTol',1e-4,'AbsTol',1e-4);
[T,Y] = ode45(@odefun,tspan,y0,options);

% t = 6 to t = 12
tspan = [6 12];
y0 = Y(end,:)+randn; % Pertubation
[t,y] = ode45(@odefun,tspan,y0,options);
T = [T;t(2:end)];  % Remove first value as it will be same as last of previous run
Y = [Y;y(2:end,:)];

Matlab 编辑器可能会提示数组 T 和 Y 没有被预分配和/或增长，但在这种情况下没关系，因为它们以大块的形式增长只有几次。

您不想尝试做的(但很多尝试)是将您的扰动添加到积分函数中。首先，如果您希望它恰好在特定时间或满足特定条件时发生，则无法在 ODE 函数内完成此操作。其次，在 ODE 中插入较大的不连续性可能会导致精度降低和计算时间延长(尤其是使用 ode45 - ode15s 可能是更好的选择，或者至少确保绝对和相对容差合适)。你会有效地产生一个非常stiff system .

如果您确实希望在满足特定条件时(而不是在特定时间)发生扰动，那么您需要使用事件函数(例如，请参阅 this question )。当满足条件时，您需要使用事件函数终止积分，然后应用您的扰动，如上面所示。