Info

Diese Frage ist geschlossen. Öffnen Sie sie erneut, um sie zu bearbeiten oder zu beantworten.

which differential equation solver should I use?

1 Ansicht (letzte 30 Tage)
February
February am 15 Feb. 2015
Geschlossen: MATLAB Answer Bot am 20 Aug. 2021
I have a differential equation of the form following:
dy(t)/dt = ay(t) + bx(t) + cx(t)/dt
y(t) is the output I am looking for, when x(t) is a given input (discrete samples from collected data). a, b, c are constants.
I followed ode45 instructions for differential equation of time dependent variables. To test out the code, I set x(t) = sin(t). Accordingly I set x(t)/dt = cos(t). this probably is not the way I will test out the real data since they are noisy discrete samples. I need to figure out how to find calculate x(t)/dt for those but that's another question.
I got a straight linear line with positive slope for output y(t). I was expected oscillating function like modified sin(t) and cos(t). I tried on wolfram alpha and that displayed oscillating (damping) graph for y(t), so I think I did not solve this equation correctly on MATLAB. I am pretty sure I implemented ode45 correctly, so I am wondering I may be using a wrong differential equation solver.
Could you give me any feedback? Thanks.
--------------------------------------
t = linspace(0,pi*20,pi*20/0.01);
Pa = sin(t)';
dPa = cos(t)';
f = 4 * -2.6021e-07 * ones(size(Pa,1),1);
g = 4 * 2.6021e-07 * Pa + .6 * dPa;
Tspan = [0 pi*20];
IC = 5 * ones(size(Pa,1),1);
t = t';
Tspan = Tspan';
[T Pic] = ode45(@(t,Pic) f1(t,Pic,t,f,t,g),Tspan, IC');
plot(t,Pa,'b') % input
hold on;
plot(T,Pic(:,45)) % one possible output
%
and function f1 is defined as following:
function dydt = f1(t,y,ft,f,gt,g)
dydt = f.*y + g;
end
  3 Kommentare
Torsten
Torsten am 17 Feb. 2015
You have to supply g at the times indicated by the variable "t" in dydt. Thus your function routine must look like
function dydt = f1(t,y,ft,f,gt,g)
f = 4*(-2.6021e-07);
Pa = sin(t);
dPa = cos(t);
g = 4*2.6021e-07*Pa + .6*dPa;
dydt=f*y+g;
end
And you only have to supply one value as initial condition for y, not a vector of length 5*size(Pa).
Best wishes
Torsten.
February
February am 20 Feb. 2015
Thank you, Torsten. I finally figured out.

Antworten (0)

Diese Frage ist geschlossen.

Tags

Community Treasure Hunt

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

Start Hunting!

Translated by