the goal with this code is to simulate numerically the movement of three bodies in three dimensions. All of the bodies have a constant masses. We are going to determinate the evolution of the position of the bodies in tn relation to time. However when we run the code the RHS function malfuctions becasue "not enough inputs have been given" and we cant figure out why.
% initial conditions and parameters
clc;
clear all;
M1=1; %mass of body 1
M2=2; %mass of body 2
M3=0.5; %mass of body 3
% the position of the three bodies in the initial conditions
P1_0=[0,0,0];
P2_0=[1,0,0];
P3_0=[0,0,1];
P1=P1_0;
P2=P2_0;
P3=P3_0;
% the initial speed of the three bodies
a=0.05;
V1_0=[a,a,a];
V2_0=[-a,a,a];
V3_0=[-a,-a,-a];
V1=V1_0;
V2=V2_0;
V3=V3_0;
X0=[P1_0,V1_0,P2_0,V2_0,P3_0,V3_0]'; % the initial positions and speeds of the three bodies in
one column vector.
[t, X]= ode45(@(t,X) RHS(t,X), [0 3], X0);
h=0.1;
K1=(h*RHS(X0));
K2=h*(RHS(X0+0.5*K1));
K3=h*(RHS(X0-K1+(2*K2)));
X1=X0+(1/6)*(K1+4*K2+K3);
function Y=RHS(t,X)
V1=X(4:6);
V2=X(10:12);
V3=X(16:18);
P1=X(1:3);
P2=X(7:9);
P3=X(13:15);
M1=1;
M2=2;
M3=0.5;
Y=[V1;Force(P1,M1,P2,M2,P3,M3)/M1;V2;Force(P2,M2,P1,M1,P3,M3)/M2;V3;Force(P3,M3,P1,M1,P2,M2)/M3];
end
function F=Force(P,M,Pa,Ma,Pb,Mb)
F=(-M*Ma*(P-Pa))/(norm(P-Pa)^3)-(M*Mb*(P-Pb))/(norm(P-Pb)^3);
end
Thomas:
I just ran your code on my end and it works flawlessly up to and including the callout to ode45
. I did, however, get the same error that you claimed in the question. The error is triggered by the line K1=(h*RHS(X0))
, which I believe is the beginning of a block where you're trying to set up a Runge-Kutta method. The reason you get the error is that in that line you only give your RHS
function a single input, yet, you've defined RHS
to take two inputs (t
, and X
):
function Y=RHS(t,X)
V1=X(4:6);
V2=X(10:12);
V3=X(16:18);
P1=X(1:3);
P2=X(7:9);
P3=X(13:15);
M1=1;
M2=2;
M3=0.5;
Y=[V1;Force(P1,M1,P2,M2,P3,M3)/M1;V2;Force(P2,M2,P1,M1,P3,M3)/M2;V3;Force(P3,M3,P1,M1,P2,M2)/M3];
end
Since the RHS
function does not use t
at all (autonomous system of ODE's), the following will circumvent the error:
h=0.1;
K1=(h*RHS([],X0));
K2=h*(RHS([],X0+0.5*K1));
K3=h*(RHS([],X0-K1+(2*K2)));
X1=X0+(1/6)*(K1+4*K2+K3);
MATLAB was just looking for that additional input, that's all.
NOTE: There was nothing special about inputting an empty vector []
, really this "fix" would have worked with any inputs (including a string
or a complex number). The trick was to occupy that placeholder input.