Skip to content
This repository was archived by the owner on Aug 28, 2024. It is now read-only.

Commit

Permalink
Fixed problem caused by trying to compute energy errors when there ar…
Browse files Browse the repository at this point in the history
…e no massive bodies. Use total energy for the denominator of normalizations rather than orbital energy, which is 0 in that case.
  • Loading branch information
daminton committed Feb 27, 2024
1 parent 7f87b1d commit e673113
Showing 1 changed file with 8 additions and 5 deletions.
13 changes: 8 additions & 5 deletions src/swiftest/swiftest_io.f90
Original file line number Diff line number Diff line change
Expand Up @@ -168,16 +168,19 @@ module subroutine swiftest_io_conservation_report(self, param, lterminal)
end if

if (.not.param%lfirstenergy) then

nbody_system%ke_orbit_error = (ke_orbit_now - nbody_system%ke_orbit_orig) / abs(nbody_system%E_orbit_orig)
nbody_system%ke_spin_error = (ke_spin_now - nbody_system%ke_spin_orig) / abs(nbody_system%E_orbit_orig)
nbody_system%pe_error = (pe_now - nbody_system%pe_orig) / abs(nbody_system%E_orbit_orig)
nbody_system%ke_orbit_error = (ke_orbit_now - nbody_system%ke_orbit_orig) / abs(nbody_system%te_orig)
nbody_system%ke_spin_error = (ke_spin_now - nbody_system%ke_spin_orig) / abs(nbody_system%te_orig)
nbody_system%pe_error = (pe_now - nbody_system%pe_orig) / abs(nbody_system%te_orig)

be_cb_orig = -(3 * cb%GM0**2 / param%GU) / (5 * cb%R0)
nbody_system%be_error = (be_now - nbody_system%be_orig) / abs(nbody_system%te_orig) + (be_cb_now - be_cb_orig) &
/ abs(nbody_system%te_orig)

nbody_system%E_orbit_error = (E_orbit_now - nbody_system%E_orbit_orig) / abs(nbody_system%E_orbit_orig)
if (abs(nbody_system%E_orbit_orig) < 10*tiny(1.0_DP)) then
nbody_system%E_orbit_error = 0.0_DP
else
nbody_system%E_orbit_error = (E_orbit_now - nbody_system%E_orbit_orig) / abs(nbody_system%E_orbit_orig)
end if
nbody_system%Ecoll_error = nbody_system%E_collisions / abs(nbody_system%te_orig)
nbody_system%E_untracked_error = nbody_system%E_untracked / abs(nbody_system%te_orig)
nbody_system%te_error = (nbody_system%te - nbody_system%te_orig - nbody_system%E_collisions - nbody_system%E_untracked)&
Expand Down

0 comments on commit e673113

Please sign in to comment.