How to solve a system of coupled, differential equations of motions?

Hi,

This follows a previous post that I posted, where I was wrong with formulation and assumptions. I learned the mistake there (I was not forming the equations correctly), and now I'm trying to solve the system of equations following this post and this, but I'm not getting anywhere :(

Here's the equation system, where a1 ~ a21 are just known coefficients and m and I are mass and inertia.

then, following the example on the page I've linked earlier (Solve Equations of Motion for Baton Thrown into Air), I formed the equations as follows. Here, I know some terms are wrong, and that was due to the typing.

I have the function of the mass matrix, although I realize that this is redundant, because I don't necessary solve anything through this. A, B, data are the inputs with the coefficients, mass and inertia values.

function M = learnMass(A, B, data) % Mass matrix function
m = data.m;
I = data.I;
M = zeros(12,12);
M(1,1)  = 1; M(3,3) = 1; M(5,5) = 1; M(7,7) = 1; M(9,9) = 1; M(11,11) = 1;
M(2,2)  = m + A(1,1,1);
M(2,6)  = A(3,1,1);
M(2,10) = A(5,1,1);
M(4,4)  = m + A(2,2,1);
M(4,8)  = A(4,2,1);
M(4,12) = A(6,2,1);
M(6,2)  = A(1,3,1);
M(6,6)  = m + A(3,3,1);
M(6,10) = I(5,5) + A(5,3,1); % add I
M(8,8)  = I(4,4) + A(4,4,1); % add I
M(10,2) = A(1,5,1);
M(10,6) = A(3,5,1);
M(10,10)= A(5,5,1);
M(12,4) = A(2,6,1);
M(12,8) = A(4,6,1);
M(12,12)= I(6,6) + A(6,6,1); % add I
end

Then, I have formed the RHS of the equation system as follows.

function sdot = learnF(A, B, s, data) % Equations to solve

sdot = [s(2);
        data.m*s(4)*s(12) + Ctn*A(2,2,1)*s(4)*s(12) + A(2,4,1)*s(8)*s(12) + A(2,6,1)*s(12)^2;
        s(4);
        -data.m*s(2)*s(12) - A(1,1,1)*s(2)*s(12) - A(2,4,1)*s(8)
        s(6);
        -data.m*data.g;
        s(8);
        data.I(4,4)*s(10)*s(12);
        s(10);
        -data.I(4,4)*s(8)*s(12);
        s(12);
        (data.I(4,4)-data.I(6,6))*s(8)*s(10) + A(1,1,1)*s(2)*s(4);
        -A(2,2,1)*s(2)*s(4) - A(2,4,1)*s(2)*s(8) - A(2,6,1)*s(2)*s(8) - A(6,4,1)*s(8)];

end

I have my tspan=[0,100]; and initial conditions (6 velocities) as y0=[1;0;0;0;0;0];

I don't know how to pass the data and solve the set of equations, and I appreciate any guidence on this.

[Torsten]

You have 12 equations - thus you need 12 initial conditions: for positions and velocities.

[Jake]

Hi, thank you for the response!

Then, suppose I choose o displacement too, as the other initial conditions.

Also, do you think I'm on the right track for solving the aforementioned equation system (I mean, with the matrix approach)?

[Torsten]

With the mass matrix approach, you write the equations as

M*y' = f(t,y)

and you can use several ODE solvers (ODE15S, ODE45,...) for the solution.

With the approach I suggested, you write the equations as

0 = M*y' - f(t,y)

and you can only use ODE15I.

So I think you will be successful.

Mass = @(t,y) learnMass(A, B, data);
Fun = @(t,y) learnF(A, B, y, data)
options = odeset('Mass',Mass);
S0 = a vector of size 12x1 of initial conditions
tspan = a vector of output times
[T,S] = ode15s(Fun,tspan,y0,options)

[Jake]

Thank you so much @Torsten! I'm able to solve the set of EoM now without errors. I will check the accuracy with comparison, but your above responses helped a lot! I'd gladly accept these as the answer, if I could.