ode45 and euler not working for random signal

2 次查看(过去 30 天)
When trying to simulate a first order system subject to a random input that changes value at fixed intervals, neither my ode45 nor euler integration attempts make sense to me. Neither of them follow the inut signal. Please help!
(This code reqs 2016b or later I think)
clc;clear;close all;format compact
w_filt = 10*2*pi; % first order filter
t = linspace(0,1,1000);
dt = t(2)-t(1);
inits = [0];
%% Try ODE 45
Eqns = @(t,s) eqns(t,s,dt,w_filt);
[t,x] = ode45(Eqns,t,inits);
for i = 1:length(t)
[dx(i,:),ext(i)] = Eqns(t(i),x(i));
end
in = ext;
du = dx;
u = x;
%% Try Euler
X = inits;
for i = 1:length(t)-1
dX(i) = w_filt*(in(i) - X(i));
X(i+1) = du(i)*dt + X(i);
end
dX(end+1) = 0;
%% Plots
figure(1);
title('ODE45');
plot(t,u,t,du,t,in)
legend('u','du','input')
figure(2);
title('Euler Integration')
plot(t,X,t,dX,t,in)
legend('u','du','input')
%% EOM
function [dx,ext] = eqns(t,x,dt,w_relax)
global randval
u = x;
in = randval;
if rem(t,100*dt) == 0 || t == 0
in = 0.1*pi/180*(2*rand-1);
randval = in;
end
du = w_relax*(in - u);
dx = du;
ext = in;
end

回答(1 个)

James Tursa
James Tursa 2020-6-5
编辑:James Tursa 2020-6-5
You can't use random inputs with ode45( ). ode45( ) relies on the ability to call the derivative function at arbitrary times to control esimated errors. It may even call your derivative function in non-increasing times because of this. It needs the derivatives to be consistent in order for this to work. Having randomness in your derivative function makes the derivatives inconsistent from call to call. ode45( ) will either fail with an error message because it can't get consistent estimated error results, or worse it will give you a garbage answer. You need to use another method besides ode45( ) to solve your stochastic problem.
  3 个评论
James Tursa
James Tursa 2020-6-5
And ...? Are you saying that your Euler scheme is not working, or doesn't give the expected results?

请先登录,再进行评论。

类别

Help CenterFile Exchange 中查找有关 Programming 的更多信息

标签

产品

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by