Please tell me how to use ode45 code
조회 수: 2 (최근 30일)
이전 댓글 표시
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
Ameer Hamza
2020년 6월 4일
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?
채택된 답변
Ameer Hamza
2020년 6월 4일
편집: Ameer Hamza
2020년 6월 4일
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
Ameer Hamza
2020년 6월 4일
편집: Ameer Hamza
2020년 6월 4일
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
추가 답변 (1개)
Stephan
2020년 6월 4일
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
Ameer Hamza
2020년 6월 4일
편집: Ameer Hamza
2020년 6월 4일
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')
참고 항목
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!