3

我正在 Matlab 中试验 ode45。我已经学会了如何将参数传递给 ode 函数,但我仍然有一个问题。假设我想计算汽车的轨迹(速度曲线),并且我有一个函数,例如getAcceleration,它可以为我提供汽车的加速度以及正确的档位:[acceleration, gear] = getAcceleration(speed,modelStructure)其中modelStructure代表汽车的型号。

ode 函数将是:

function [dy] = car(t,y,modelStructure)

dy           = zeros(2,1);
dy(1)        = y(2);
[dy(2),gear] = getAcceleration(y(1),modelStructure);

然后我以这种方式调用 Ode45 积分器:

tInit = 0;
tEnd  = 5,
[t,y] = ode45(@car,[tInit tEnd], [speedInitial,accelerationInitial],options,modelStructure);

问题是:如何获得矢量存储齿轮?我应该有类似的东西[t,y,gear]=ode45(....)还是应该geary向量内?


我一直在处理我的代码并使用事件函数,我现在能够获得汽车“齿轮”的变化(作为事件)。现在我有一个与相同代码相关的新问题。想象一下,当我评估 de 'dy' 向量时,我能够得到一个进一步的值 Z,这让我可以大幅加快调用加速度计算 (getAcceleration):

function [dy] = car(t,y,modelStructure)

dy           = zeros(2,1);
dy(1)        = y(2);
[dy(2),Z(t)] = getAcceleration(y(1),modelStructure,Z(t-1)); 

并假设我也能够在初始条件下计算 Z。问题是我无法计算 Z 导数。

有没有办法通过 Z 值抛出步进而不集成它?

多谢你们。

4

1 回答 1

5

首先:为什么微分方程的初始值是初始速度(speedInitial)和初始加速度(accelerationInitial)?这意味着微分方程car将在每个时间计算加速度 ( y(1)) 和加加速度 ( y(2)),即加速度的时间导数t。这似乎不正确......我会说初始值应该是初始位置(positionInitial)和初始速度(speedInitial)。但是,我不知道你的型号,我可能是错的。

现在,直接进入gear解决方案:你不能,不是没有黑客攻击ode45。这也是合乎逻辑的;你也不能dy一直直接得到,可以吗?只是不是这样ode45设置的。

我在这里看到两种方法:

全局变量

免责声明不要使用这种方法。它只是在这里展示大多数人会做的第一次尝试。

您可以存储gear在全局变量中。这可能是最少的编码,但也是最不方便的结果:

global ts gear ii

ii    = 1;
tInit = 0;
tEnd  = 5,
[t,y] = ode45(...
    @(t,y) car(t,y,modelStructure), ...
    [tInit tEnd], ...
    [speedInitial, accelerationInitial], options);

...

function [dy] = car(t,y,modelStructure)
global ts gear ii

dy    = zeros(2,1);
dy(1) = y(2);
[dy(2),gear(ii)] = getAcceleration(y(1),modelStructure);

ts(ii) = t;
ii = ii + 1;

但是,由于 的性质ode45,这将为您提供一系列时间ts和相关联gear,其中包含中间点和/或被拒绝的点ode45。所以,你必须在之后过滤那些:

ts( ~ismember(ts, t) ) = [];

我再说一遍:这不是我推荐的方法。仅在测试或做一些快速而肮脏的事情时使用全局变量,但总是很快地转向其他解决方案。此外,全局变量在 的每次(子)迭代时都会增长ode45,这是不可接受的性能损失。

最好使用下一种方法:

解决后调用

这对你的情况也不是太难,我建议你走的路。首先,修改微分方程如下,并正常求解:

tInit = 0;
tEnd  = 5,
[t,y] = ode45(...
    @(t,y) car(t,y,modelStructure), ...
    [tInit tEnd], ...
    [speedInitial, accelerationInitial], options);

...

function [dy, gear] = car(t,y,modelStructure)    

dy    = [0;0];
dy(1) = y(2);
[dy(2),gear] = getAcceleration(y(1),modelStructure);

然后ode45完成后,执行以下操作:

gear = zeros(size(t));
for ii = 1:numel(t)
    [~, gear(ii)] = car(t(ii), y(ii,:).', modelStructure); 
end

这将使您获得汽车有时会配备的所有齿轮t

car我在这里可以看到的唯一缺点是,您将拥有比单独使用更多的功能评估ode45。但这只是一个真正的问题,如果每次评估都car需要几秒钟或更长时间,我怀疑您的设置不是这种情况。

于 2012-11-27T08:13:43.137 回答