Issues solving an ODE system for the n-body problem (Edited)
Show older comments
I need to solve the ODE system of the n-body problem (in 2D) in order to get future positions and velocities from given initial values of the system and afterwards make a simulation of their behavior.
The ODE system is

where r are the positions, v the velocities, m the masses and Fg(r) the gravitational force.
I use ode45(odefun,tspan,y0). I need to define a function, I call it ode_nbody, which holds the mentioned system and is used in this fashion ode45(ode_nbody,tspan,y0). ode_nbody receive as a parameter the masses. The function I've written is shown below
function [ddt] = ode_nbody(t, y, m)
% y = [r1x, r2x, ..., rnx, r1y, r2y, ..., rny, v1x, v2x, ..., vnx, v1y, v2y, ..., vny]
G = 6.67408E-11;
n = length(m);
r = [y(1:n), y(n+1:2*n)]'; % Vector of positions that makes easier the acceleration computation
for i = 1:n
acc(:,i) = zeros(2,1);
for j = 1:n
if i ~= j
acc(:,i) = acc(:,i) + ...
G * (m(j) .* (r(:,j) - r(:,i)))./(norm(r(:,j) - r(:,i))^3);
end
end
end
% ddt = [v1x, v2x, ..., vnx, v1y, v2y, ..., vny, a1x, a2x, ..., anx]
ddt(1:2*n) = y(2*n+1:4*n);
ddt(2*n+1:4*n) = [acc(1,1:n) acc(2,1:n)];
ddt = ddt';
end
The results is

As you may have noticed, all the planets converge to the sun (that yellow "o" in the center) which, obviusly was not expected.
4 Comments
I recovered the removed content from the Google cache (something which anyone can do). Editing away your question is very rude. Someone spent time reading your question, understanding your issue, figuring out the solution, and writing an answer. Now you repay that kindness by ensuring that the next person with a similar question can't benefit from this answer.
M. Ezequiel
on 20 Jul 2021
Rik
on 20 Jul 2021
Once an answer has been posted, you are no longer able to delete a question. I have understood Mathworks has (or will soon) added a remark explaining this when you post a question. If you did not intend to be rude, then there is no longer an issue, simply leave the question as-is.
If you have already solved it, feel free to post the solution as an answer and accept it. That way you might help someone in the future with a similar issue.
M. Ezequiel
on 20 Jul 2021
Answers (1)
Walter Roberson
on 18 Jul 2021
0 votes
You did not say what the difficulty is ?
Note that you are only tracking the orbits for 10 seconds
1 Comment
Rik
on 20 Jul 2021
Oh yep, sorry! The difficulty is that if you put a bigger interval, for example [0 1000000], the trajectories converge to the sun (the yellow "o" in the image) but they're suppose to depict a circular orbit.

Categories
Find more on Gravitation, Cosmology & Astrophysics in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!