N-body problem: Difference between revisions

Line 1,049:
<p>Although in general the n-body problem doesn’t have an explicit solution, certain configurations do, and these 'nice' configurations can be used to validate and perform unbiased comparisons of numeric solvers. A fairly simple configuration with an exact explicit solution is a radially symmetric system of three bodies orbiting at constant velocity around a central axis:
<ul>Given:
<li>n = 3 <em>-- number of bodies</em></li>
<li>g = 1 <em>-- gravitational constant</em></li>
<li>m = 1 <em>-- mass per body</em></li>
<li>r = 1 <em>-- radius</em></li>
</ul><ul>Then:
<li>f = (1/3)^(1/2) <em>-- centripetal force</em></li>
<li>a = (1/3)^(1/2) <em>-- centripetal acceleration</em></li>
<li>v = (1/3)^(1/4) <em>-- velocity</em></li>
</ul></p>
<p>This is a variation on the [Klemperer rosette[https://en.wikipedia.org/wiki/Klemperer_rosette]]. </p>
Line 1,181:
<h5>n>3</h5>
<p>We can generate and execute tests for radially symmetrical configurations of any number of bodies, run the simulations, and compare the results with an explicitly determined result:</p>
<lang J>r =. 4 : '%:+/*:(1-cos 2*pi*y%x),sin 2*pi*y%x' NB. distance to nth body
<lang J>GENM =: 3 : 0"0
f0 =. 4 : '1%((x r y)^2)' NB. force due to nth body
N =: y
f =. 4 : '(x f0 y)*(1-cos 2*pi*y%x)%x r y'"0 NB. centripetal component of force
VS =: (a N)^0.5
a =. 3 : '+/(y f i. y)' NB. centripetal acceleration
TS =: 2*pi%VS
GENM =: 3 : 0"0 NB. accumulated discrepancy over one rotation for n bodies
dt =: 0.001
N =: y NB. number of bodies
IS =: >.2*pi%VS*dt
VS =: (a N)^0.5 NB. velocity for static equilibrium
IS ITER VS GEN N
TS =: 2*pi%VS NB. rotational time
E =: >./>./({.out)-{:out
dt =: 0.001 NB. time increment for simulation
N;VS;TS;E
IS =: >.2*pi%VS*dt NB. number of steps for simulation
IS ITER VS GEN N NB. run the simulation
E =: >./>./({.out)-{:out NB. get max difference between initial and final state
N;VS;TS;E NB. return the results
)
smoutput ' n';'velocity';' period';' error'