Please tell me how to use ode45 code
2 views (last 30 days)
Show older comments
num=134;
den=[1 16 134];
x=tf(num,den)
subplot(1,2,1)
step(x)
title('(a)')
I want to express this code as ode45, how do I do it?
2 Comments
Ameer Hamza
on 4 Jun 2020
Original Question:
Please tell me how to use ode45 code.
num=134;
den=[1 16 134];
x=tf(num,den)
subplot(1,2,1)
step(x)
title('(a)')
I want to express this code as ode45, how do I do it?
Accepted Answer
Ameer Hamza
on 4 Jun 2020
Edited: Ameer Hamza
on 4 Jun 2020
Following is one way of using ode45 to solve this
[t, y] = ode45(@odeFun, [0 1], [0; 0]);
plot(t, y(:,1), 'o-')
function dydt = odeFun(t, y)
% transfer function is equivalent to following ODE
% y'' = -16y'-124y+134u
u = 1; % step input
dydt = zeros(2, 1);
dydt(1) = y(2);
dydt(2) = -16*y(2)-124*y(1)+134*u;
end
3 Comments
Ameer Hamza
on 4 Jun 2020
Edited: Ameer Hamza
on 4 Jun 2020
Check my answer on your other question. I paste the code here for reference
[t, y] = ode45(@odeFun, [0 1], [0; 0]);
plot(t, y(:,1), 'o-')
function dxdt = odeFun(t, x)
% transfer function is equivalent to following ODE
% x'' = 4/9(-16x'+ea(t))
ea = (1-x(1))*303;
if ea > 100
ea = 100;
elseif ea < -100
ea = -100;
end
dxdt = zeros(2, 1);
dxdt(1) = x(2);
dxdt(2) = 4/9*(-16*x(2)+ea);
end
More Answers (1)
Stephan
on 4 Jun 2020
You need to perform a inverse laplace transformation:
num=134;
den=[1 16 134];
x=tf(num,den)
subplot(1,2,1)
step(x)
title('(a)')
syms s
ode = matlabFunction(ilaplace(134/(s^2 + 16*s + 134)),'Vars',{'t','y'})
[t, y] = ode45(ode,[0 1],0);
subplot(1,2,2)
plot(t,y)
title('(b) - with ode45')
xlabel('Time (seconds)')
ylabel('Amplitude')
3 Comments
Ameer Hamza
on 4 Jun 2020
Edited: Ameer Hamza
on 4 Jun 2020
ode45 is not needed here. It worked in this specific case because of the step input, but for any other input (impulse, ramp), this method will not work. Following is the general method for inverse Laplace
syms s
u = 1/s; % laplace of step signal
ode = matlabFunction(ilaplace(134/(s^2 + 16*s + 134)*u),'Vars','t');
subplot(1,2,2)
fplot(ode, [0 1])
title('(b) - with ode45')
xlabel('Time (seconds)')
ylabel('Amplitude')
See Also
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!