How to solve the coupled second order differential equations by using "ODE45"?

4 views (last 30 days)
Hello,
I try to solve the coupled second order differential equations with ODE45.
But my answer was weird.
My system is two defree of freedom system including spring and damper.
So, the equations are
m1*x1'' + (c1+c2)*x1' + (k1+k2)*x1 - c2*x2' - k2*x2 = 0;
m2*x2'' + x2*x2' + k2*x2 - c2*x1' - k2*x1 = F(t);
And I changed the above equations,
x1'' = v1' = ...
x2'' = v2' = ...
Here are my code..
% mx1'' = -(k1+k2)x1 - (c1+c2)x1' + k2 x2 + c2 x2'
% mx2'' = -k2 x2 - c2 x2' + k2 x1 + c2 x1' + F
tspan = linspace(0,10*pi,2000);
h = tspan(2) - tspan(1);
xini = [0; 0; 0; 0];
[t,x] = ode45(@T2_func, tspan, xini);
------ Function ------
function dxdt = T2_func(t,x)
%%Constants
m1 = 1; % [kg]
m2 = 1; % [kg]
k1 = 10; % [N/m]
k2 = 10; % [N/m]
c1 = 1; % [Ns/m]
c2 = 1; % [Ns/m]
%%Force
p = 10 * ones(1,length(t)); % [N]
%%Matrix
F = zeros(2, length(t));
F(2,:) = p/m2;
K11 = -(k1+k2)/m1;
K12 = k2/m1;
K21 = -k2/m2;
K22 = k2/m2;
Kmat = [K11 K12; K21 K22];
C11 = -(c1+c2)/m1;
C12 = c2/m1;
C21 = -c2/m2;
C22 = c2/m2;
Cmat = [C11 C12; C21 C22];
%%Function
dxdt = zeros(4,1);
dxdt(1) = x(2);
dxdt(2) = - Kmat(1,1)*x(1) - Cmat(1,1)*x(2) + Kmat(1,2)*x(3) + Cmat(1,2)*x(4) + F(1);
dxdt(3) = x(4);
dxdt(4) = - Kmat(2,1)*x(3) - Cmat(2,1)*x(4) + Kmat(2,2)*x(1) + Cmat(2,2)*x(2) + F(2);
Thank you.

Accepted Answer

Torsten
Torsten on 23 Apr 2018
You already included the correct signs of the coefficients in the definition of Kmat and Cmat. Thus all terms in dxdt(2) and dxdt(4) must have a "+" sign.
Furthermore, in your original set of equations, K11 = -k1/m1 instead of K11 = -(k1+k2)/m1.
Additionally, t in T2_func will always be a scalar. Thus
%%Force
p = 10 ; % [N]
%%Matrix
F(1) = 0;
F(2) = p/m2;
Best wishes
Torsten.
  2 Comments
HYEOKJUNE LEE
HYEOKJUNE LEE on 23 Apr 2018
Thanks for answer,
You mentioned, K- matrix. And that is my fault.. so I modified my script, and plus I also modified C- matrix.
And, through your answer, my code be working correctly. Thank you very much.
HYEOKJUNE LEE
HYEOKJUNE LEE on 23 Apr 2018
Torsten, I have more question. When the force is related with the time, then how can I make a code?
You mentioned, "time" in a function scription must be scalar. But the force term could be dependent with the time.
How could I handle this problem?
Best regardes, Lee.

Sign in to comment.

More Answers (0)

Community Treasure Hunt

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

Start Hunting!