2

我想对我的 ode45 运动方程计算有两个限制:位置和时间。我已经让时间事件开始工作,但我不确定是否以及如何添加另一个事件来限制位置。编辑:我在一个 ODE 方程中也有许多不同的粒子耦合在一起,并且一旦它们到达“屋顶”就需要它们单独停止,因为它们都以不同的速度行进......我能通过事件来实现这一点吗?我有一个想法,我将如何做到这一点,但它非常复杂,可能会很慢......

4

1 回答 1

0

我不确定你是否可以做你想做的事,但可以用事件做很多事情。首先,这听起来像是某种首次通过时间(又名首次击球时间)的数值计算。如果这些“粒子”是随机的,请停止并且不要使用ode45适合 SDE 的方法。

据我所知,您可以拥有多少个事件函数 - 或者更确切地说是事件函数的维度(类似于您的 ODE 函数的维度) - 并且它们的数量与您拥有多少个 ODE 方程无关. events 函数接收当前时间和当前状态向量。您可以使用其中的任何或所有元素来创建每个事件。您是正确的,更多的事件函数和更复杂的事件会减慢集成速度。性能还取决于检测事件的频率。如果您的每个粒子都到达您所说的“屋顶”,并且只触发一个事件,那么这不会太糟糕。

在实现方面,这里有一个基于 Matlab 示例的简单示例ballode,仅在垂直方向上模拟 N 个弹道粒子。有 N 个非终止事件来捕捉每个粒子通过 y = 0 时的时间和速度。添加一个额外的终止事件来检查是否所有粒子都通过了 y = 0(如果我们知道这是哪个粒子就像我们在这里所做的那样,我们可以使该事件成为终止事件)。

function eventsdemo

% Initial conditions for n balls
n = 10;
y0(2*n,1) = 0;
y0(n+1:end) = linspace(20,40,n);

% Specify events function
options = odeset('Events',@(t,y)efun(t,y,n));

% Integrate
[t,y,te,ye,ie] = ode45(@(t,y)f(t,y,n),[0 10],y0,options);

figure;
plot(t,y(:,1:n),'b',te(1:n),ye(sub2ind(size(ye),ie(1:n),(1:n).')),'r.');

function dydt = f(t,y,n)
% Differential equations for ballistic motion
dydt = [y(n+1:end);zeros(n,1)-9.8];

function [value,isterminal,direction] = efun(t,y,n)
% Last event checks that all balls have hit ground and terminates integration
yn = y(1:n);
value = [yn;all(yn < 0)];
zn = zeros(n,1);
isterminal = [zn;1];
direction = [zn-1;1];

在某些方面,这有点低效,因为我们继续模拟所有 N 个系统,即使它们中的一些已经通过 y = 0。但是,它很简单,并且输出数组是矩形的而不是参差不齐的。

我不清楚您所说的“耦合在一起”并需要粒子“停止”是什么意思。如果您需要做的不仅仅是记录事件数据,例如更改系统参数或以其他方式更改微分方程,那么您需要在每个事件后终止并重新开始积分。查看ballode示例(edit ballode在 Matlab 命令窗口中键入)以查看一些建议以使其更有效。

于 2013-05-22T00:47:00.307 回答