# Error using ode45

1 view (last 30 days)
Ibro Tutic on 13 Dec 2015
Commented: Walter Roberson on 13 Dec 2015
I am trying to model the position of the Earth relative to the sun. I am using a vector form of the acceleration on a mass due to another mass. I have used this same formula in an Euler-Cromer scheme without error. However, when I try to implement it using the ode45 solver, I get this error:
In an assignment A(:) = B, the number of elements in A and B
must be the same.
Error in orbit (line 10)
dp(2) = -G * Msun * R / (sum(R.*R)^(3/2));
This is my function:
function dp = orbit(t, P)
%where p = [xi yi vxi vyi]
G = 6.673E-11/(1000^3);
Psun = [0 0];
Msun=1.988e30;
R = Psun - P(2);
M = 1.98892E30;
dp = zeros(4,1);
dp(1)= P(2);
dp(2) = -G * Msun * R / (sum(R.*R)^(3/2));
end
And my script to call it:
[T,Y] = ode45(@orbit,[0 12],[1.5e8 0 0 29.75]);

Walter Roberson on 13 Dec 2015
You define
Psun = [0 0];
so Psun is a vector. You then define
R = Psun - P(2);
so R is a vector.
Then when you have
dp(2) = -G * Msun * R / (sum(R.*R)^(3/2));
because you have that vector R on the right hand side, before the division, your result is going to be a vector. You are using /, the mrdivide operator, but (sum(R.*R)^(3/2)) is a scalar so R / (sum(R.*R)^(3/2)) is the same size as the vector R. So your overall expression is a vector, and you are trying to store that vector in the single location dp(2)
##### 2 CommentsShow 1 older commentHide 1 older comment
Walter Roberson on 13 Dec 2015
dp(2:3) = -G * Msun * R / (sum(R.*R)^(3/2));
But then your ode would need an input vector of length 3 (because the output vector length must equal the input vector length)
You could put the two of them together by encoding as a single complex number, but then it would do complex integration on the term.