How to solve simultaneous second order differential equations

I'm trying to solve 3 simultaneous differential equations equations but I can't get it to work. This is what I have tried so far and it's not working :

kb= 1.888*10^9 ; db= 5.123*10^7 ; mb= 5452000 ; Ib= 1.76667*10^9 ; Rb= 0.281 ;

kt= 1.2519*10 ; dt= 2.869*10^7 ; mt= 6797460 ; It= 3.3428*10^9 ; Rt= 64.2 ;

kTMD= 28805 ; dTMD= 10183 ; mTMD= 40000 ; RTMD= 90.6 ;

syms thb(t) x(t) tht(t)

ode1= Ib*diff(thb, t, 2)== -db*diff(thb,t)-kb*thb-(mb*9.81*Rb*thb)+dt*(diff(tht,t)-diff(thb,t)) ;

ode2= It*diff(tht, t, 2) == mt*9.81*Rt*tht-kt*(tht-thb)-dt*(diff(tht,t)-diff(thb,t))-kTMD*RTMD*(RTMD*tht-x)-dTMD*RTMD*(RTMD*diff(tht,t)-diff(x,t))-mTMD*9.81*(RTMD*tht-x)-RTMD;

ode3= mTMD*diff(x, t, 2) == kTMD*(RTMD*tht- x)+ dTMD*(RTMD*diff(tht,t)-diff(x,t))+mTMD*9.81*tht ;

odes= [ode1; ode2 ;ode3] ;

cond1 = tht(0)== 0;
cond2 = diff(tht, 0)==0;
cond3 = thb(0)==0;
cond4= diff(thb, 0)==0;
cond5 = x==0;
cond6 = diff(x,0)==0;
conds= [cond1; cond2; cond3; cond4; cond5; cond6] ;

[thbSol(t), thtSol(t), xSol(t)] = dsolve(odes, conds)

Any help would be appreciated

[Stephan]

Do you have values for:

 db kb mb Rb dt dTMD mTMD It Ib mt Rt kt kTMD RTMD

If you want to solve this system with only symbolic variables i suspect you will need a huge amount of time...

[Fouad Elias]

Hi,

Yes I have values for all of these. I didn't include them in this post, but I have edited it now. I'm not sure the approach I'm using to solve these 3 simultaneous equations is the correct one.

[Stephan]

If you provide the values maybe we can help

Answer 1

[Stephan]

Hi,

i think an analytical solution wil be hard to find (perhaps impossible) - for numerical solution you can do this:

syms thb(t) tht(t) x(t) Dxt(t) Dthtt(t) Dthbt(t)

kb= 1.888*10^9 ;
db= 5.123*10^7 ;
mb= 5452000 ;
Ib= 1.76667*10^9 ;
Rb= 0.281 ;
kt= 1.2519*10;
dt= 2.869*10^7;
mt= 6797460;
It= 3.3428*10^9 ;
Rt= 64.2;
kTMD= 28805;
dTMD= 10183;
mTMD= 40000;
RTMD= 90.6 ;

ode1= -db*diff(thb,t)-kb*thb-(mb*9.81*Rb*thb)+dt*(diff(tht,t)-diff(thb,t))...
    - Ib*diff(thb, t, 2);
[ode1, var1] = reduceDifferentialOrder(ode1,thb);

ode2 = mt*9.81*Rt*tht-kt*(tht-thb)-dt*(diff(tht,t)-diff(thb,t))...
    -kTMD*RTMD*(RTMD*tht-x)-dTMD*RTMD*(RTMD*diff(tht,t)-diff(x,t))...
    -mTMD*9.81*(RTMD*tht-x)-RTMD-It*diff(tht, t, 2);
[ode2, var2] = reduceDifferentialOrder(ode2,tht);

ode3= kTMD*(RTMD*tht- x)+ dTMD*(RTMD*diff(tht,t)-diff(x,t))+...
    mTMD*9.81*tht-mTMD*diff(x, t, 2);
[ode3, var3] = reduceDifferentialOrder(ode3,x);

ode = [ode1; ode2; ode3];
var = [var1; var2; var3];

[odes, vars] = odeToVectorField(ode);
ode_fun = matlabFunction(odes, 'Vars',{'t','Y'});

y0 = [0 0 0 0 0 0];
tspan = [0 15];

[t,y] = ode45(ode_fun,tspan,y0);

Dthbt = y(:,1);
tht = y(:,2);
thb = y(:,3);
x = y(:,4);
Dthtt = y(:,5);
Dxt = y(:,6);

plot(t,thb,'bo',t,tht,'ro',t,x,'mo',t,Dthbt,'bx',t,Dthtt,'rx',t,Dxt,'mx')
legend({'thb','tht','x','Dthbt','Dthtt','Dxt'}, 'Location','southwest')

This code integrates the system from t=0...15 - you can change the value if needed by changing tspan.

Best regards

Stephan