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

Commit

Permalink
Browse files Browse the repository at this point in the history 
Simplified symba_frag_pos energy diagnostic to only show the final result, not the intermediate ones.
  • Loading branch information
@daminton
daminton committed May 20, 2021
1 parent eeb2aa5 commit 314a5ff
Showing 1 changed file with 41 additions and 39 deletions.
80 changes: 41 additions & 39 deletions src/symba/symba_frag_pos.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -176,48 +176,48 @@ subroutine symba_frag_pos(param, symba_plA, family, x, v, L_spin, Ip, mass, radi
call define_pre_collisional_family()

do try = 1, ntry
write(*, "(' -------------------------------------------------------------------------------------')")
write(*,*) "Try number: ",try, ' of ',ntry
write(*, "(' -------------------------------------------------------------------------------------')")
write(*, "(' First pass to get angular momentum ')")
!write(*, "(' -------------------------------------------------------------------------------------')")
!write(*,*) "Try number: ",try, ' of ',ntry
!write(*, "(' -------------------------------------------------------------------------------------')")
!write(*, "(' First pass to get angular momentum ')")

call set_fragment_position_vectors()
call set_fragment_tangential_velocities()
call calculate_system_energy(linclude_fragments=.true.)

write(*, "(' -------------------------------------------------------------------------------------')")
write(*, "(' | T_orb T_spin T pe Etot Ltot')")
write(*, "(' -------------------------------------------------------------------------------------')")
write(*,fmtlabel) ' change |',(ke_orb_after - ke_orb_before) / abs(Etot_before), &
(ke_spin_after - ke_spin_before)/ abs(Etot_before), &
(ke_orb_after + ke_spin_after - ke_orb_before - ke_spin_before)/ abs(Etot_before), &
(pe_after - pe_before) / abs(Etot_before), &
dEtot, dLmag

write(*, "(' -------------------------------------------------------------------------------------')")
write(*, "(' If T_offset > 0, this will be a merge')")
write(*,fmtlabel) ' T_offset |',ke_offset / abs(Etot_before)
!write(*, "(' -------------------------------------------------------------------------------------')")
!write(*, "(' | T_orb T_spin T pe Etot Ltot')")
!write(*, "(' -------------------------------------------------------------------------------------')")
!write(*,fmtlabel) ' change |',(ke_orb_after - ke_orb_before) / abs(Etot_before), &
! (ke_spin_after - ke_spin_before)/ abs(Etot_before), &
! (ke_orb_after + ke_spin_after - ke_orb_before - ke_spin_before)/ abs(Etot_before), &
! (pe_after - pe_before) / abs(Etot_before), &
! dEtot, dLmag
!
! write(*, "(' -------------------------------------------------------------------------------------')")
! write(*, "(' If T_offset > 0, this will be a merge')")
! write(*,fmtlabel) ' T_offset |',ke_offset / abs(Etot_before)

call set_fragment_radial_velocities(lmerge)

call calculate_system_energy(linclude_fragments=.true.)

if (.not.lmerge) then
write(*, "(' -------------------------------------------------------------------------------------')")
write(*, "(' Second pass to get energy ')")
write(*, "(' -------------------------------------------------------------------------------------')")
write(*,fmtlabel) ' T_frag targ |',ke_target / abs(Etot_before)
write(*,fmtlabel) ' T_frag calc |',ke_frag / abs(Etot_before)
write(*, "(' T_offset should now be small')")
write(*,fmtlabel) ' T_offset |',ke_offset / abs(Etot_before)
write(*, "(' -------------------------------------------------------------------------------------')")
write(*, "(' | T_orb T_spin T pe Etot Ltot')")
write(*, "(' -------------------------------------------------------------------------------------')")
write(*,fmtlabel) ' change |',(ke_orb_after - ke_orb_before) / abs(Etot_before), &
(ke_spin_after - ke_spin_before)/ abs(Etot_before), &
(ke_orb_after + ke_spin_after - ke_orb_before - ke_spin_before)/ abs(Etot_before), &
(pe_after - pe_before) / abs(Etot_before), &
dEtot, dLmag
!write(*, "(' -------------------------------------------------------------------------------------')")
!write(*, "(' Second pass to get energy ')")
!write(*, "(' -------------------------------------------------------------------------------------')")
!write(*,fmtlabel) ' T_frag targ |',ke_target / abs(Etot_before)
!write(*,fmtlabel) ' T_frag calc |',ke_frag / abs(Etot_before)
!write(*, "(' T_offset should now be small')")
!write(*,fmtlabel) ' T_offset |',ke_offset / abs(Etot_before)
!write(*, "(' -------------------------------------------------------------------------------------')")
!write(*, "(' | T_orb T_spin T pe Etot Ltot')")
!write(*, "(' -------------------------------------------------------------------------------------')")
!write(*,fmtlabel) ' change |',(ke_orb_after - ke_orb_before) / abs(Etot_before), &
! (ke_spin_after - ke_spin_before)/ abs(Etot_before), &
! (ke_orb_after + ke_spin_after - ke_orb_before - ke_spin_before)/ abs(Etot_before), &
! (pe_after - pe_before) / abs(Etot_before), &
! dEtot, dLmag
end if

lmerge = lmerge .or. ((Etot_after - Etot_before) / abs(Etot_before) > 0._DP)
Expand All @@ -227,14 +227,16 @@ subroutine symba_frag_pos(param, symba_plA, family, x, v, L_spin, Ip, mass, radi
write(*, "(' -------------------------------------------------------------------------------------')")
write(*, "(' Final diagnostic')")
write(*, "(' -------------------------------------------------------------------------------------')")
write(*, "(' dL_tot: Should be very small' )")
write(*,fmtlabel) ' dL_tot |', dLmag
write(*, "(' dE_tot: If not a merge, should be negative and equal to Qloss' )")
write(*,fmtlabel) ' dE_tot |', dEtot
write(*,fmtlabel) ' Qloss |', -Qloss / abs(Etot_before)
write(*,fmtlabel) ' dE - Qloss |', (Etot_after - Etot_before + Qloss) / abs(Etot_before)
write(*, "(' -------------------------------------------------------------------------------------')")
write(*,*) "Flag as merge?", lmerge
if (lmerge) then
write(*,*) "symba_frag_pos can't find a solution, so this collision is flagged as a merge"
else
write(*, "(' dL_tot should be very small' )")
write(*,fmtlabel) ' dL_tot |', dLmag
write(*, "(' dE_tot should be negative and equal to Qloss' )")
write(*,fmtlabel) ' dE_tot |', dEtot
write(*,fmtlabel) ' Qloss |', -Qloss / abs(Etot_before)
write(*,fmtlabel) ' dE - Qloss |', (Etot_after - Etot_before + Qloss) / abs(Etot_before)
end if
write(*, "(' -------------------------------------------------------------------------------------')")

call restore_scale_factors()
Expand Down

0 comments on commit 314a5ff

Please sign in to comment.