I Want to solve system of three differential equations numerically. How can i solve it with ode45, or are there any other ways?

Hello. I wrote three equations for three different variables. I also have conditions.

syms v1(t) y1(t) y2(t) v2(t) a2(t) phi(t) Omega(t) Omegaac(t) omega(t) omegaac(t)
Omega = diff(phi,t);
Omegaac = diff(Omega,t);
omegaac = diff(omega,t);
v1 = diff(y1,t);
v2 = diff(y2,t);
a2 = diff(v2,t);
v1_y = v2 - (omega*r+Omega*L*cos(beta))*sin(phi);
v1_x = (omega*r+Omega*L*cos(beta))*cos(phi);
v1 = sqrt((v1_y)^2+(v1_x)^2);
cond1 = v1(0) == 0;
cond2 = v1_y(0)==0;
cond3 =v1_x(0)==0;
cond4 = y1(0) == 0;
cond5 = y2(0) == 0;
cond6 = v2(0) == 0;
cond7 = a2(0) == 0;
cond8 = phi(0) == 0;
cond9 = Omega(0)==0;
cond10 = Omegaac(0)==0;
cond11 = omega(0) == 0;
cond12 = omegaac(0)==0;
%equations look a bit like this, coefficients which will be denoted by some
%letters are known. only unknowns are variables.
ode1 = a2==(A+B*(Omega)^2*cos(phi)-C-(D*(v1_y))/v1);
ode2 = omegaac==((E*omega+F*Omega-v2*cos(phi))/v1);
ode2 = Omegaac==(G*sin(phi)+H*a2*sin(phi)-(I*omega+J*Omega-v2*cos(phi))/v1);

[Torsten]

Could you specify the functions by their respective names that are known and that are unknown ?

[Nikoloz]

There are two point of contact in (in motion which theese equations exist). y1 and v1 are one of the points coordinate and velocity, y2 v2 a2 is another points coordinate,velocity, acceleration. phi is some angle tilted, itsderivative with respect to time is \Omega and Omega'ac' is angular acceleration. omega is another variable which denotes angular velocity in another spinning motion and its derivative is omega'ac' is another angular acceleration. every Latin letter I used are just coefficients , just numbers which I can just plug in. v1_y and v1_x are just y and x axis components. at a starting point everything is zero. I used that information to write initial conditions.

[Torsten]

You don't need to explain what the variables mean.

Just state which of them are known and which of them are unknown.

For the unknowns, additionally give the initial values.

Answer 1

[Star Strider]

syms v1(t) y1(t) y2(t) v2(t) a2(t) phi(t) L r A B C D E F G H I J beta Omega(t) Omegaac(t) omega(t) omegaac(t) Y t 
sympref('AbbreviateOutput',false);
Omega = diff(phi,t);
Omegaac = diff(Omega,t);
omegaac = diff(omega,t);
v1 = diff(y1,t);
v2 = diff(y2,t);
a2 = diff(v2,t);
v1_y = v2 - (omega*r+Omega*L*cos(beta))*sin(phi);
v1_x = (omega*r+Omega*L*cos(beta))*cos(phi);
v1 = sqrt((v1_y)^2+(v1_x)^2);
cond1 = v1(0) == 0;
cond2 = v1_y(0)==0;
cond3 =v1_x(0)==0;
cond4 = y1(0) == 0;
cond5 = y2(0) == 0;
cond6 = v2(0) == 0;
cond7 = a2(0) == 0;
cond8 = phi(0) == 0;
cond9 = Omega(0)==0;
cond10 = Omegaac(0)==0;
cond11 = omega(0) == 0;
cond12 = omegaac(0)==0;
%equations look a bit like this, coefficients which will be denoted by some
%letters are known. only unknowns are variables.
ode1 = a2==(A+B*(Omega)^2*cos(phi)-C-(D*(v1_y))/v1);
ode2 = omegaac==((E*omega+F*Omega-v2*cos(phi))/v1);
ode3 = Omegaac==(G*sin(phi)+H*a2*sin(phi)-(I*omega+J*Omega-v2*cos(phi))/v1);

[VF,Sbs] = odeToVectorField([ode1; ode2; ode3])                                 % ADDED

VF = 

Sbs = 

odesfcn = matlabFunction(VF, 'Vars',{t,Y,A,B,C,D,E,F,G,H,I,J,L,beta,r})                                      % ADDED

odesfcn = function_handle with value:
    @(t,Y,A,B,C,D,E,F,G,H,I,J,L,beta,r)[Y(2);(cos(Y(1)).*Y(4).*sqrt(Y(4).^2+r.^2.*cos(Y(1)).^2.*Y(5).^2+r.^2.*sin(Y(1)).^2.*Y(5).^2-r.*sin(Y(1)).*Y(4).*Y(5).*2.0+L.^2.*cos(Y(1)).^2.*cos(beta).^2.*Y(2).^2+L.^2.*sin(Y(1)).^2.*cos(beta).^2.*Y(2).^2-L.*sin(Y(1)).*cos(beta).*Y(2).*Y(4).*2.0+L.*r.*cos(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0+L.*r.*sin(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0)-I.*Y(5).*sqrt(Y(4).^2+r.^2.*cos(Y(1)).^2.*Y(5).^2+r.^2.*sin(Y(1)).^2.*Y(5).^2-r.*sin(Y(1)).*Y(4).*Y(5).*2.0+L.^2.*cos(Y(1)).^2.*cos(beta).^2.*Y(2).^2+L.^2.*sin(Y(1)).^2.*cos(beta).^2.*Y(2).^2-L.*sin(Y(1)).*cos(beta).*Y(2).*Y(4).*2.0+L.*r.*cos(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0+L.*r.*sin(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0)-J.*Y(2).*sqrt(Y(4).^2+r.^2.*cos(Y(1)).^2.*Y(5).^2+r.^2.*sin(Y(1)).^2.*Y(5).^2-r.*sin(Y(1)).*Y(4).*Y(5).*2.0+L.^2.*cos(Y(1)).^2.*cos(beta).^2.*Y(2).^2+L.^2.*sin(Y(1)).^2.*cos(beta).^2.*Y(2).^2-L.*sin(Y(1)).*cos(beta).*Y(2).*Y(4).*2.0+L.*r.*cos(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0+L.*r.*sin(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0)+G.*sin(Y(1)).*(Y(4).^2+r.^2.*cos(Y(1)).^2.*Y(5).^2+r.^2.*sin(Y(1)).^2.*Y(5).^2-r.*sin(Y(1)).*Y(4).*Y(5).*2.0+L.^2.*cos(Y(1)).^2.*cos(beta).^2.*Y(2).^2+L.^2.*sin(Y(1)).^2.*cos(beta).^2.*Y(2).^2-L.*sin(Y(1)).*cos(beta).*Y(2).*Y(4).*2.0+L.*r.*cos(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0+L.*r.*sin(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0)+A.*H.*sin(Y(1)).*(Y(4).^2+r.^2.*cos(Y(1)).^2.*Y(5).^2+r.^2.*sin(Y(1)).^2.*Y(5).^2-r.*sin(Y(1)).*Y(4).*Y(5).*2.0+L.^2.*cos(Y(1)).^2.*cos(beta).^2.*Y(2).^2+L.^2.*sin(Y(1)).^2.*cos(beta).^2.*Y(2).^2-L.*sin(Y(1)).*cos(beta).*Y(2).*Y(4).*2.0+L.*r.*cos(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0+L.*r.*sin(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0)-C.*H.*sin(Y(1)).*(Y(4).^2+r.^2.*cos(Y(1)).^2.*Y(5).^2+r.^2.*sin(Y(1)).^2.*Y(5).^2-r.*sin(Y(1)).*Y(4).*Y(5).*2.0+L.^2.*cos(Y(1)).^2.*cos(beta).^2.*Y(2).^2+L.^2.*sin(Y(1)).^2.*cos(beta).^2.*Y(2).^2-L.*sin(Y(1)).*cos(beta).*Y(2).*Y(4).*2.0+L.*r.*cos(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0+L.*r.*sin(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0)-D.*H.*sin(Y(1)).*Y(4).*sqrt(Y(4).^2+r.^2.*cos(Y(1)).^2.*Y(5).^2+r.^2.*sin(Y(1)).^2.*Y(5).^2-r.*sin(Y(1)).*Y(4).*Y(5).*2.0+L.^2.*cos(Y(1)).^2.*cos(beta).^2.*Y(2).^2+L.^2.*sin(Y(1)).^2.*cos(beta).^2.*Y(2).^2-L.*sin(Y(1)).*cos(beta).*Y(2).*Y(4).*2.0+L.*r.*cos(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0+L.*r.*sin(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0)+B.*H.*cos(Y(1)).*sin(Y(1)).*Y(2).^2.*(Y(4).^2+r.^2.*cos(Y(1)).^2.*Y(5).^2+r.^2.*sin(Y(1)).^2.*Y(5).^2-r.*sin(Y(1)).*Y(4).*Y(5).*2.0+L.^2.*cos(Y(1)).^2.*cos(beta).^2.*Y(2).^2+L.^2.*sin(Y(1)).^2.*cos(beta).^2.*Y(2).^2-L.*sin(Y(1)).*cos(beta).*Y(2).*Y(4).*2.0+L.*r.*cos(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0+L.*r.*sin(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0)+D.*H.*r.*sin(Y(1)).^2.*Y(5).*sqrt(Y(4).^2+r.^2.*cos(Y(1)).^2.*Y(5).^2+r.^2.*sin(Y(1)).^2.*Y(5).^2-r.*sin(Y(1)).*Y(4).*Y(5).*2.0+L.^2.*cos(Y(1)).^2.*cos(beta).^2.*Y(2).^2+L.^2.*sin(Y(1)).^2.*cos(beta).^2.*Y(2).^2-L.*sin(Y(1)).*cos(beta).*Y(2).*Y(4).*2.0+L.*r.*cos(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0+L.*r.*sin(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0)+D.*H.*L.*sin(Y(1)).^2.*cos(beta).*Y(2).*sqrt(Y(4).^2+r.^2.*cos(Y(1)).^2.*Y(5).^2+r.^2.*sin(Y(1)).^2.*Y(5).^2-r.*sin(Y(1)).*Y(4).*Y(5).*2.0+L.^2.*cos(Y(1)).^2.*cos(beta).^2.*Y(2).^2+L.^2.*sin(Y(1)).^2.*cos(beta).^2.*Y(2).^2-L.*sin(Y(1)).*cos(beta).*Y(2).*Y(4).*2.0+L.*r.*cos(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0+L.*r.*sin(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0))./(Y(4).^2+r.^2.*cos(Y(1)).^2.*Y(5).^2+r.^2.*sin(Y(1)).^2.*Y(5).^2-r.*sin(Y(1)).*Y(4).*Y(5).*2.0+L.^2.*cos(Y(1)).^2.*cos(beta).^2.*Y(2).^2+L.^2.*sin(Y(1)).^2.*cos(beta).^2.*Y(2).^2-L.*sin(Y(1)).*cos(beta).*Y(2).*Y(4).*2.0+L.*r.*cos(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0+L.*r.*sin(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0);Y(4);(-D.*Y(4)+A.*sqrt(Y(4).^2+r.^2.*cos(Y(1)).^2.*Y(5).^2+r.^2.*sin(Y(1)).^2.*Y(5).^2-r.*sin(Y(1)).*Y(4).*Y(5).*2.0+L.^2.*cos(Y(1)).^2.*cos(beta).^2.*Y(2).^2+L.^2.*sin(Y(1)).^2.*cos(beta).^2.*Y(2).^2-L.*sin(Y(1)).*cos(beta).*Y(2).*Y(4).*2.0+L.*r.*cos(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0+L.*r.*sin(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0)-C.*sqrt(Y(4).^2+r.^2.*cos(Y(1)).^2.*Y(5).^2+r.^2.*sin(Y(1)).^2.*Y(5).^2-r.*sin(Y(1)).*Y(4).*Y(5).*2.0+L.^2.*cos(Y(1)).^2.*cos(beta).^2.*Y(2).^2+L.^2.*sin(Y(1)).^2.*cos(beta).^2.*Y(2).^2-L.*sin(Y(1)).*cos(beta).*Y(2).*Y(4).*2.0+L.*r.*cos(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0+L.*r.*sin(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0)+B.*cos(Y(1)).*Y(2).^2.*sqrt(Y(4).^2+r.^2.*cos(Y(1)).^2.*Y(5).^2+r.^2.*sin(Y(1)).^2.*Y(5).^2-r.*sin(Y(1)).*Y(4).*Y(5).*2.0+L.^2.*cos(Y(1)).^2.*cos(beta).^2.*Y(2).^2+L.^2.*sin(Y(1)).^2.*cos(beta).^2.*Y(2).^2-L.*sin(Y(1)).*cos(beta).*Y(2).*Y(4).*2.0+L.*r.*cos(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0+L.*r.*sin(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0)+D.*r.*sin(Y(1)).*Y(5)+D.*L.*sin(Y(1)).*cos(beta).*Y(2)).*1.0./sqrt(Y(4).^2+r.^2.*cos(Y(1)).^2.*Y(5).^2+r.^2.*sin(Y(1)).^2.*Y(5).^2-r.*sin(Y(1)).*Y(4).*Y(5).*2.0+L.^2.*cos(Y(1)).^2.*cos(beta).^2.*Y(2).^2+L.^2.*sin(Y(1)).^2.*cos(beta).^2.*Y(2).^2-L.*sin(Y(1)).*cos(beta).*Y(2).*Y(4).*2.0+L.*r.*cos(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0+L.*r.*sin(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0);(-cos(Y(1)).*Y(4)+E.*Y(5)+F.*Y(2)).*1.0./sqrt(Y(4).^2+r.^2.*cos(Y(1)).^2.*Y(5).^2+r.^2.*sin(Y(1)).^2.*Y(5).^2-r.*sin(Y(1)).*Y(4).*Y(5).*2.0+L.^2.*cos(Y(1)).^2.*cos(beta).^2.*Y(2).^2+L.^2.*sin(Y(1)).^2.*cos(beta).^2.*Y(2).^2-L.*sin(Y(1)).*cos(beta).*Y(2).*Y(4).*2.0+L.*r.*cos(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0+L.*r.*sin(Y(1)).^2.*cos(beta).*Y(2).*Y(5).*2.0)]

Refer to it in the ode45 call (or whatever function you want to integrate it with) as:

@(t,Y)odesfcn(t,Y,A,B,C,D,E,F,G,H,I,J,L,beta,r)

Be sure that all the additional arguments (A,B,C,D,E,F,G,H,I,J,L,beta,r) are defined in your workspace first.

.

[Nikoloz]

Hello. I'm a bit new. I just want to ask you something. In which part of your code are initial conditions used, if you didn't need to use them, why? wouldn't they make solving easier? Thank you for your assistance! edit: maybe index (0) is identified whilst solving and it automatically assumes that it is condition for variable?

[Star Strider]

No worries!

I do not use them with odeToVectorField or matlabFunction. They are supplied as arguments when you use ode45 (or whatever solver you choose) with your system. See that documentation for details.