2 degrees of freedom mass-spring system

Hi, I tried to code a script and a function in order to resolve a system composed by two masses and two springs of different values, but I had some troubles in writing the function that resolves the equations of motion.

m1=1; %mass 1 [kg]
m2=2; %mass 2 [kg]
k1=100; %spring 1 [N/m]
k2=150; %spring 2 [N/m]


M=[m1 0;0 m2]; %mass matrix
K=[k1+k2 -k2; -k2 k2]; %stiffness matrix

%% Solving system

tRange=[0 20];
x0=[0.3 0.1 1 4] ;
odefun=@(t,y) solve_sys_2nd_order(t,y,M,K);
[tsol,ysol]=ode45(odefun,tRange,x0);

%%Plot of the response
subplot(2,1,1)
plot(tsol,ysol(:,1),'b');
xlabel('tempo [s]');
ylabel('spostamento [m]');
subplot(2,1,2)
plot(tsol,ysol(:,2),'r')
xlabel('tempo [s]');
ylabel('velocità [m/s]');

The function I wrote is given by

function dy=solve_sys_2nd_order(t,y,m1,m2,k1,k2)

dy=zeros(4,1);

dy(1)=y(3);
dy(2)=y(4);

dy=[y(2);y(4);dy(3);-inv(M)*K*y];

But when I run the code, an error appears saying:

Error using *

Incorrect dimensions for matrix multiplication. Check that the number of columns in the first matrix matches the number of

rows in the second matrix. To perform elementwise multiplication, use '.*'.

Error in solve_sys_2nd_order (line 8)

dy=[y(2);y(4);dy(3);-inv(M)*K*y];

Error in two_dof_system>@(t,y)solve_sys_2nd_order(t,y,M,K)

Error in odearguments (line 90)

f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.

Error in ode45 (line 115)

odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);

Error in two_dof_system (line 24)

[tsol,ysol]=ode45(odefun,tRange,x0);

Can anyone help me please?

Answer 1

[James Tursa]

In this line

dy=[y(2);y(4);dy(3);-inv(M)*K*y];

looks like M and K are both 2x2 but y is 4x1, hence the error.

Also, your arguments don't match in these two line:

odefun=@(t,y) solve_sys_2nd_order(t,y,M,K);
function dy=solve_sys_2nd_order(t,y,m1,m2,k1,k2)

In the first case you have M and K, but in the second case you have m1, m2, k1, k2.

[federico gottardi]

Ok thanks, I corrected the last thing you noticed, but have you any idea about how can I fix the dy vector? I thought to write each single differential equation but I would like to have a matrix form to simplify things if I add another dof or spring.

[James Tursa]

The state vector y that gets passed in will be a vector, but your other variables M and K can be passed in as matrices. You can always reshape y inside your derivative function if needed. What specific difficulty are you having with this? What does your latest code look like?

[federico gottardi]

Ok, probably I fixed, could you please say if I did the right thing?

(I modified the function)

function dy=solve_sys_2nd_order(t,y,M,K)

dy=zeros(4,1);

dy(1)=y(3);
dy(2)=y(4);

dy=[y(2);y(4);-inv(M)*K*[y(1);y(2)];

the script is

clc;
clear;
close all;
% Sistema composto da due masse unite da due molle --> si hanno 2 gdl

%% Proprietà del sistema

m1=150; %massa 1 [kg]
m2=200; %massa 2 [kg]
k1=100; %ridgidezza della molla 1 [N/m]
k2=150; %rigidezza della molla 2 [N/m]

%% Scrittura delle matrici di massa e rigidezza ricavate dal PLV

M=[m1 0;0 m2];
K=[k1+k2,-k2; -k2,k2];
C=[0;0];

%% Risoluzione del sistema di equazioni differenziali

tRange=[0 20];
x0=[0 0 0 0.5] ;
odefun=@(t,y) solve_sys_2nd_order(t,y,M,K);
[tsol,ysol]=ode45(odefun,tRange,x0);

%% Grafico della risposta del sistema
subplot(2,1,1)
plot(tsol,ysol(:,1),'b'); hold on
plot(tsol,ysol(:,2),'g'); hold off
xlabel('tempo [s]');
ylabel('spostamento [m]');
legend('x1','x2');
subplot(2,1,2)
plot(tsol,ysol(:,3),'r'); hold on
plot(tsol,ysol(:,4),'k'); hold off
xlabel('tempo [s]');
ylabel('velocità [m/s]');
legend('xdot1','xdot2');

[LEKIREDDY BHARADHWAJ REDDY]

it's not working