State-Space for Inerter based suspension system

Hello, I'm doing my project on the comparison of Conventional suspension system with an inerter based suspension system. I'm struck at deriving the state-space equations for the parallel based inerter. Here the acceleration of other mass acts as a loop. I tried but i wasn't able to solve. I'll be grateful if someone could give me a hint to crack this. The Equations are given below.

Where Z0 is an input, Z1 is the Unsprung Mass & Z2 is the Sprung Mass.

These are the Matrices i got.

Amat = [0,1,0,0       % A_mat of the Unsprung and Sprung Mass
    -((Ks/Mu)+(Kt/Mu)), -Cs/Mu, Ks/Mu, Cs/Mu
    0,0,0,1
    Ks/Ms, Cs/Ms, -Ks/Ms, -Cs/Ms];
Bmat = [0;Kt/Mu;0;0]; % B_Matrix of the Unsprung mass and Sprung Mass
Cmat = [1,0,0,0       % C_Matrix of the Unsprung Mass and Sprung Mass
    0,1,0,0
    0,0,1,0
    0,0,0,1];
Dmat = [0];           % D_Matrix of the Unsprung Mass and Sprung Mass

[darova]

What about ode45?

[Aneesur Rahman Shafeer]

No. I mean my task is to model the suspension system using State-Space Equations.

Answer 1

[Ameer Hamza]

See odetovectorfield :https://www.mathworks.com/help/symbolic/odetovectorfield.html. It also have example to convert the symbolic ODE to the system of first order ODE which can be used to write a ss-model.

If you want to do it by hand, then you first need to write the two ODEs in this form

and then convert it into the state-space form

where .

Following shows the code using ode2vectorfield:

syms ms mu ks cs binrt kt z0(t) z1(t) z2(t)

Dz1 = diff(z1);
DDz1 = diff(Dz1);
Dz2 = diff(z2);
DDz2 = diff(Dz2);
eq1 = ms*DDz2 == ks*(z1-z2) + cs*(Dz1-Dz2) + binrt*(DDz1-DDz2);
eq2 = mu*DDz1 == ks*(z1-z2) + cs*(Dz1-Dz2) + binrt*(DDz1-DDz2) + kt*(z0-z1);

[V,S] = odeToVectorField(eq1, eq2)

Result

>> V
V =
                                                                                                                                          Y[2]
 (binrt*kt*z0(t) + kt*ms*z0(t) - binrt*kt*Y[1] + cs*ms*Y[2] - cs*ms*Y[4] + ks*ms*Y[1] - ks*ms*Y[3] - kt*ms*Y[1])/(binrt*mu - binrt*ms + ms*mu)
                                                                                                                                          Y[4]
                            (binrt*kt*z0(t) - binrt*kt*Y[1] + cs*mu*Y[2] - cs*mu*Y[4] + ks*mu*Y[1] - ks*mu*Y[3])/(binrt*mu - binrt*ms + ms*mu)
>> S
S =
  z1
 Dz1
  z2
 Dz2

You can also use this outout to write a state-space model.

[Aneesur Rahman Shafeer]

Thanks for the Answer. What is Y[1], Y[2],Y[3] & Y[4] here ?

[Ameer Hamza]

Y[1], Y[2],Y[3] & Y[4] are same as those shown by vector S.

Y[1] -> z1
Y[2] -> dz1/dt
Y[3] -> z2
Y[4] -> dz2/dt

[Aneesur Rahman Shafeer]

>> V

V =

Y[2]

(binrt*kt*z0(t) + kt*ms*z0(t) - binrt*kt*Y[1] + cs*ms*Y[2] - cs*ms*Y[4] + ks*ms*Y[1] - ks*ms*Y[3] - kt*ms*Y[1])/(binrt*mu - binrt*ms + ms*mu)

Y[4]

(binrt*kt*z0(t) - binrt*kt*Y[1] + cs*mu*Y[2] - cs*mu*Y[4] + ks*mu*Y[1] - ks*mu*Y[3])/(binrt*mu - binrt*ms + ms*mu)

>> S

S =

z1

Dz1

z2

Dz2

Using this ouput as differential equations and then convert it into state-space, right ?

[Ameer Hamza]

Yes, you just need to seperate the cofficients of z1, z2, z1dot, z2dot and z0 to get the matrices A and B.

[Aneesur Rahman Shafeer]

Alright. Wil try as you said.

I tried doing hand calculations and then i got this. I'm not sure if this is right or not

Amat_p = [0,1,0,0

(-(Kt+Ks)/(Mu+b))+((b*Ks)/(Ms*Mu+Ms*b+Mu*b))-((b*b)*(Kt+Ks))/((Ms*Mu+Ms*b+Mu*b)*(Mu+b)), ((-Cs/(Mu+b))+((b*Cs)/(Ms*Mu+Ms*b+Mu*b))-((b*b)*(Cs)/(Ms*Mu+Ms*b+Mu*b)*(Mu+b))), ((Ks/(Mu+b))-((b*Ks)/(Ms*Mu+Ms*b+Mu*b))+((b*b)*(Ks)/(Ms*Mu+Ms*b+Mu*b)*(Mu+b))), ((Cs/(Mu+b))-(b*Cs/(Ms*Mu+Ms*b+Mu*b))+((b*b)*Cs)/((Ms*Mu+Ms*b+Mu*b)*(Mu+b)))

0,0,0,1

(Ks*(Mu+b)/(Ms*Mu+Ms*b+Mu*b))-(b*(Kt+Ks)/(Ms*Mu+Ms*b+Mu*b)), ((Cs*(Mu+b)/(Ms*Mu+Ms*b+Mu*b))-(b*Cs/(Ms*Mu+Ms*b+Mu*b))), ((-Ks*(Mu+b)/(Ms*Mu+Ms*b+Mu*b))+(b*Ks/(Ms*Mu+Ms*b+Mu*b))), ((-Cs*(Mu+b)/(Ms*Mu+Ms*b+Mu*b))+(b*Cs/(Ms*Mu+Ms*b+Mu*b)))];

Bmat_p = [0

(Kt/(Mu+b))+((b*b)*Kt/(Ms*Mu+Ms*b+Mu*b)*(Mu+b))

0

(b*Kt/(Ms*Mu+Ms*b+Mu*b))];

Cmat_p = eye(4,4);

Dmat_p = 0;

[Ameer Hamza]

This looks quite complicate as compared to the answers given by odeToVectorField. I suggest to compare it with that solution to see if there is some issues.

[Aneesur Rahman Shafeer]

I did the state-space and ran the simuation for 10 secs but the output of the state-space turns to be weird.

These are the matrices i got following your method

Amat_p = [0 1 0 0

(-b*Kt+Ks*Mu-Kt*Ms)/(b*Mu-b*Ms+Ms*Mu), (Cs*Ms)/(b*Mu-b*Ms+Ms*Mu),

(-Ks*Ms)/(b*Mu-b*Ms+Ms*Mu), (-Cs*Ms)/(b*Mu-b*Ms+Ms*Mu)

0 0 0 1

(-b*Kt+Ks*Mu)/(b*Mu-b*Ms+Ms*Mu), (Cs*Mu)/(b*Mu-b*Ms+Ms*Mu),

(-Ks*Mu)/(b*Mu-b*Ms+Ms*Mu), (-Cs*Mu)/(b*Mu-b*Ms+Ms*Mu)];

Bmat_p = [0

(b*Kt+Kt*Ms)/(b*Mu-b*Ms+Ms*Mu)

0

b*Kt/(b*Mu-b*Ms+Ms*Mu)];

Cmat_p = eye(4,4);

Dmat_p = [0;0;0;0];

This is the displacement of sprung mass for differential equations.

and this is for state-space

[Ameer Hamza]

I guess there may be a mistake in the derivation of SS equation. Did you derive it by hand?

[Aneesur Rahman Shafeer]

yes

[Ameer Hamza]

I suggest you to write matrix A and B using the output of odeToVectorField as I showed in my answer.

[Aneesur Rahman Shafeer]

Yes. I got the answer now. Though there's a slight variations, eventually will solve this. Thank You so much for your help.

[Ameer Hamza]

I am glad to be of help. I guess you can carefully recheck the parameters again. There might be some minor mistake which causes the difference.