4.3. Intermezzo
From this point on, Alice and Bob spent quite a bit of time coding up
an individual time step algorithm, using their five new classes. As
one can imagine, they ran into all kinds of problems, and went through
many a debug session. Finally, after several days, it all worked. We
will join them now, while they look at their first finished product,
in file world1.rb. This will be the first of a series of such
implementations, and the detailed layout will change quite a bit.
However, the current code contains all the basic ingredients for a
worldline-based N-body code.
Bob: That was quite a lot of work. The new file world1.rb
is a complete overhaul of what we had in nbody_ind1.rb.
It is twice a long, too.
Alice: But it is a lot more structured, and I'm sure it will be much
easier to extend to more complicated situations. Besides, it has
significantly more functionality.
Bob: Before we congratulate ourselves too much, we should make sure
that we get the same answers as before. Let us apply our standard test
again, using the by now familiar initial conditions for our 4-particle
system:
|gravity> kali mkplummer.rb -n 4 -s 1 | kali nbody_set_id.rb > test.in
==> Plummer's Model Builder <==
Number of particles : N = 4
pseudorandom number seed given : 1
Floating point precision : precision = 16
Incremental indentation : add_indent = 2
Screen Output Verbosity Level : verbosity = 1
ACS Output Verbosity Level : acs_verbosity = 1
Floating point precision : precision = 16
Incremental indentation : add_indent = 2
==> Takes an N-body system, and gives each body a unique ID <==
value of @body_id for 1st body : n = 1
Floating point precision : precision = 16
Screen Output Verbosity Level : verbosity = 1
ACS Output Verbosity Level : acs_verbosity = 1
Floating point precision : precision = 16
Incremental indentation : add_indent = 2
actual seed used : 1
Here is our comparison run for the older code:
|gravity> kali nbody_ind1.rb -t 1 -c 0.001 < test.in > tmpi.out
==> Individual Time Step Hermite Code <==
Parameter to determine time step size: dt_param = 0.001
Interval between diagnostics output: dt_dia = 1.0
Time interval between snapshot output: dt_out = 1.0
Duration of the integration : t = 1.0
Screen Output Verbosity Level : verbosity = 1
ACS Output Verbosity Level : acs_verbosity = 1
Floating point precision : precision = 16
Incremental indentation : add_indent = 2
at time t = 0, after 0 steps :
E_kin = 0.25 , E_pot = -0.5 , E_tot = -0.25
E_tot - E_init = 0
(E_tot - E_init) / E_init = -0
at time t = 1, after 42565 steps :
E_kin = 0.313 , E_pot = -0.563 , E_tot = -0.25
E_tot - E_init = -1.42e-12
(E_tot - E_init) / E_init = 5.68e-12
And here is what our new code gives:
|gravity> kali world1.rb -t 1 -c 0.001 < test.in > tmpw.out
==> Individual Time Step, Individual Integration Scheme Code <==
Determines the time step size : dt_param = 0.001
Duration of an era : dt_era = 0.01
Maximum time step (units dt_era): dt_max_param = 1.0
Diagnostics output interval : dt_dia = 1.0
Snapshot output interval : dt_out = 1.0
Duration of the integration : t = 1.0
Screen Output Verbosity Level : verbosity = 1
ACS Output Verbosity Level : acs_verbosity = 1
Floating point precision : precision = 16
Incremental indentation : add_indent = 2
at time t = 0 (from interpolation after 0 steps to time 0):
E_kin = 0.25 , E_pot = -0.5 , E_tot = -0.25
E_tot - E_init = 0
(E_tot - E_init) / E_init = -0
at time t = 1 (from interpolation after 42578 steps to time 1):
E_kin = 0.313 , E_pot = -0.563 , E_tot = -0.25
E_tot - E_init = -1.05e-12
(E_tot - E_init) / E_init = 4.21e-12
Aha, at least the number of time steps is very similar, and so is the
total energy error. Here is the distance in phase space
|gravity> cat tmpi.out tmpw.out | kali nbody_diff.rb
==> 6N-dimensional phase space distance between two N-body systems <==
Floating point precision : precision = 2
Screen Output Verbosity Level : verbosity = 1
ACS Output Verbosity Level : acs_verbosity = 1
Floating point precision : precision = 16
Incremental indentation : add_indent = 2
6N-dim. phase space dist. for two 4-body systems: 4.2386924923927053e-12
4.4. Almost Perfect Scaling
Alice: I'm glad to see that: the distance is comparable to the energy error.
In fact, the correspondance is so good, that I wonder whether we've
been a bit lucky. I'd feel more comfortable if we would redo the run
with somewhat lower accuracy, just to check whether the phase space
distance will become larger.
Bob: Why not. When we take time steps that are ten times longer, we get:
|gravity> kali nbody_ind1.rb -t 1 -c 0.01 < test.in > tmpi.out
==> Individual Time Step Hermite Code <==
Parameter to determine time step size: dt_param = 0.01
Interval between diagnostics output: dt_dia = 1.0
Time interval between snapshot output: dt_out = 1.0
Duration of the integration : t = 1.0
Screen Output Verbosity Level : verbosity = 1
ACS Output Verbosity Level : acs_verbosity = 1
Floating point precision : precision = 16
Incremental indentation : add_indent = 2
at time t = 0, after 0 steps :
E_kin = 0.25 , E_pot = -0.5 , E_tot = -0.25
E_tot - E_init = 0
(E_tot - E_init) / E_init = -0
at time t = 1, after 4257 steps :
E_kin = 0.313 , E_pot = -0.563 , E_tot = -0.25
E_tot - E_init = 7.24e-10
(E_tot - E_init) / E_init = -2.9e-09
I bet the new code will show a similar error:
|gravity> kali world1.rb -t 1 -c 0.01 < test.in > tmpw.out
==> Individual Time Step, Individual Integration Scheme Code <==
Determines the time step size : dt_param = 0.01
Duration of an era : dt_era = 0.01
Maximum time step (units dt_era): dt_max_param = 1.0
Diagnostics output interval : dt_dia = 1.0
Snapshot output interval : dt_out = 1.0
Duration of the integration : t = 1.0
Screen Output Verbosity Level : verbosity = 1
ACS Output Verbosity Level : acs_verbosity = 1
Floating point precision : precision = 16
Incremental indentation : add_indent = 2
at time t = 0 (from interpolation after 0 steps to time 0):
E_kin = 0.25 , E_pot = -0.5 , E_tot = -0.25
E_tot - E_init = 0
(E_tot - E_init) / E_init = -0
at time t = 1 (from interpolation after 4291 steps to time 1):
E_kin = 0.313 , E_pot = -0.563 , E_tot = -0.25
E_tot - E_init = -1.98e-10
(E_tot - E_init) / E_init = 7.92e-10
Alice: And so it does. So now let us see whether the phase space difference
is also be comparable:
|gravity> cat tmpi.out tmpw.out | kali nbody_diff.rb
==> 6N-dimensional phase space distance between two N-body systems <==
Floating point precision : precision = 2
Screen Output Verbosity Level : verbosity = 1
ACS Output Verbosity Level : acs_verbosity = 1
Floating point precision : precision = 16
Incremental indentation : add_indent = 2
6N-dim. phase space dist. for two 4-body systems: 1.0927259036948475e-08
Bob: Which it is. Okay, we seem to have done something right.