Plotting a System of Two Second-Order Differential Equations

I'm trying to reduce a system of two second-order differential equations into a system of first-order equations, solve them, and plot the result. The differential equations are x1'' = -(k1*x1 - k2*(x1 - x2))/m and x2'' = -(k2(x2 - x1))/m2.
This is what I have so far:
[t, y] = ode45(@func,[0 100],[3, 0])
figure (9)
plot(t,y(:,1),t,y(:,3))
function dy = func(t,y)
m1 = 2;
m2 = 4;
k1 = 3;
k2 = 5;
dy = zeros(4,1)
dy(1) = y(2);
dy(2) = -(-(k1 + k2)*y(1)+k2*y(3))./m1;
dy(3) = y(4);
dy(4) = (k1*y(1) - k2*y(3))./m2;
end
It says that there is an error in line 1 and errors in ode45 and odearguments. There is also an error in line 11 and that "index exceeds the number of array elements."

 채택된 답변

James Tursa
James Tursa 2019년 12월 19일
편집: James Tursa 2019년 12월 19일
You've got a 4th order system, so your initial state must contain four elements including the x1' and x2', not two. E.g.,
[t, y] = ode45(@func,[0 100],[3, 0, 0, 0]); % x1, x1', x2, x2' in initial conditions
Also, your dy(2) and dy(4) don't match your posted derivative equations, so you should double check and correct those.

댓글 수: 2

Thanks! I noticed that the equation for dy(2) also did not match what I wrote in the description, so I changed it to:
dy(2) = (-(y(1) - y(3))*k2-k1*y(1))/m1;
Is that correct?
Based on what you wrote, this equation
x1'' = -(k1*x1 - k2*(x1 - x2))/m
translates to this code
dy(2) = -(k1*y(1) - k2*(y(1) - y(3))) / m1;
and this equation
x2'' = -(k2(x2 - x1))/m2
translates to this code
dy(4) = -k2*(y(3) - y(1)) / m2;
I don't know the point of rearranging the terms in your code since it just makes it more difficult to compare directly to the derivative equation. And in fact it looks like you still don't have the dy(2) code correct, so you might consider just doing the straghtforward translation and using my supplied code.

댓글을 달려면 로그인하십시오.

추가 답변 (0개)

카테고리

도움말 센터File Exchange에서 Numerical Integration and Differential Equations에 대해 자세히 알아보기

질문:

2019년 12월 19일

편집:

2019년 12월 19일

Community Treasure Hunt

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

Start Hunting!

Translated by