Error: ode45 Must return a column vector?

1 view (last 30 days)
I'm trying to solve this problem involving a 2DOF system involving springs and dampers, and i've coded all of this below, including the state space matrix, this is urgent as this is due in 2 days, any help would be appreciated :)
Many thanks in advance
function batch_timeresp
% the total time domain response
% System parameters
k3 = 15; k1 = 10.0; k2 = 10.0; m1 = 2.0; m2 = 6.0; L = 1.0; c = 1; ct = 0.5; g = 9.81;kequi=5;
%ball
zeta1 = (c/(2*sqrt(m2*k3)));
%block
zeta2 = (ct/(2*sqrt(m1*kequi)));
%damped natural freq of block
om02=sqrt(kequi/m2);
omd2=om02*(sqrt(1-(zeta2)^2));
%damped natural freq of ball
om01=sqrt(k3/m1);
omd1=om01*(sqrt(1-(zeta1)^2));
T1 = 0.01*((2*pi)/om02);
T2 = 1.50*((2*pi)/om01);
Tpulse = T2-T1;
% Excitation parameters
Ft = 30;
Im = ((Ft*Tpulse)/2);
% State space system matrix and coupling
A=[0,0,1,0;0,0,0,1;-(k1+k2)/(m1+m2),0,-c/(m1+m2),0;0,-k3/(m1*(L^2)),0,-ct/(m1*(L^2))]; %state space matrix
S=[1,0,0,0;0,1,0,0;0,0,1,(m1*L)/(m1+m2);0,0,(m1*L)/(m1*(L^2)),1]; %coupling
Sinv=inv(S);
gmat = [0;0;-g;g/L]; %gravity terms
Fmat=[0;0;Im/(m1+m2);(2*Im*L)/(m1*(L^2))]; %impulse force
% Numerical integration
tspan=linspace(0,5,1e3);
X0=[0;0;0;0];
[T,X]=ode45(@ssmodel,tspan,X0);
% Results
figure
subplot(2,1,1), plot(T,X(:,1)), ylabel('x [m]')
subplot(2,1,2), plot(T,X(:,2)), ylabel('dx/dt [m/s]'), xlabel('time [s]')
function dx=ssmodel(t,x)
% State-space model of 1 DOF system
dx=(A*x+gmat+Fmat).*Sinv;
end
end

Accepted Answer

Alan Stevens
Alan Stevens on 16 May 2021
Like this?
% the total time domain response
% System parameters
k3 = 15; k1 = 10.0; k2 = 10.0; m1 = 2.0; m2 = 6.0; L = 1.0; c = 1; ct = 0.5; g = 9.81;kequi=5;
%ball
zeta1 = (c/(2*sqrt(m2*k3)));
%block
zeta2 = (ct/(2*sqrt(m1*kequi)));
%damped natural freq of block
om02=sqrt(kequi/m2);
omd2=om02*(sqrt(1-(zeta2)^2));
%damped natural freq of ball
om01=sqrt(k3/m1);
omd1=om01*(sqrt(1-(zeta1)^2));
T1 = 0.01*((2*pi)/om02);
T2 = 1.50*((2*pi)/om01);
Tpulse = T2-T1;
% Excitation parameters
Ft = 30;
Im = ((Ft*Tpulse)/2);
% State space system matrix and coupling
A=[0,0,1,0;0,0,0,1;-(k1+k2)/(m1+m2),0,-c/(m1+m2),0;0,-k3/(m1*(L^2)),0,-ct/(m1*(L^2))]; %state space matrix
S=[1,0,0,0;0,1,0,0;0,0,1,(m1*L)/(m1+m2);0,0,(m1*L)/(m1*(L^2)),1]; %coupling
Sinv=inv(S);
gmat = [0;0;-g;g/L]; %gravity terms
Fmat=[0;0;Im/(m1+m2);(2*Im*L)/(m1*(L^2))]; %impulse force
% Numerical integration
tspan=linspace(0,5,1e3);
X0=[0;0;0;0];
[T,X]=ode45(@(t,x) ssmodel(t,x,A,gmat,Fmat,Sinv),tspan,X0);
% Results
figure
subplot(2,1,1), plot(T,X(:,1)), ylabel('x [m]')
subplot(2,1,2), plot(T,X(:,2)), ylabel('dx/dt [m/s]'), xlabel('time [s]')
function dx=ssmodel(~,x,A,gmat,Fmat,Sinv)
% State-space model of 1 DOF system
dx=Sinv*(A*x+gmat+Fmat); %%%%%%%%%%%%%%%%%%%%%%
end
  5 Comments
Alan Stevens
Alan Stevens on 16 May 2021
You are solving for four parameters, which means the function must return a 4x1 vector. If you have
dx=(A*x+gmat+Fmat).*Sinv;
this returns a 4x4 matrix.
Your state-space equation is presumably
S*dx = A*x+gmat+Fmat;
which means you need to multiply from the left to get x' on its own.
Sinv*S*dx = Sinv*(A*x+gmat+Fmat);
or
dx = Sinv*(A*x+gmat+Fmat);
King Hin Wong
King Hin Wong on 16 May 2021
Solved! Thank you so much, have a nice day :)

Sign in to comment.

More Answers (0)

Categories

Find more on Chassis Systems in Help Center and File Exchange

Community Treasure Hunt

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

Start Hunting!