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

Commit

Permalink
Started implementing the encounter saving with user inputs
Browse files Browse the repository at this point in the history
  • Loading branch information
daminton committed Dec 6, 2022
1 parent f8c4c22 commit 155dc89
Show file tree
Hide file tree
Showing 4 changed files with 81 additions and 66 deletions.
2 changes: 1 addition & 1 deletion src/modules/symba_classes.f90
Original file line number Diff line number Diff line change
Expand Up @@ -222,7 +222,7 @@ module symba_classes
procedure :: dealloc => symba_util_dealloc_system !! Deallocates all allocatable arrays
procedure :: resize_storage => symba_util_resize_storage !! Resizes the encounter history storage object so that it contains enough spaces for the number of snapshots needed
procedure :: snapshot => symba_util_take_encounter_snapshot !! Take a minimal snapshot of the system through an encounter
procedure :: start_encounter => symba_io_start_encounter !! Initializes the new encounter and/or fragmentation save file(s)
procedure :: start_encounter => symba_io_start_encounter !! Initializes the new encounter history
procedure :: stop_encounter => symba_io_stop_encounter !! Saves the encounter and/or fragmentation data to file(s)
final :: symba_util_final_system !! Finalizes the SyMBA nbody system object - deallocates all allocatables
end type symba_nbody_system
Expand Down
1 change: 0 additions & 1 deletion src/setup/setup.f90
Original file line number Diff line number Diff line change
Expand Up @@ -68,7 +68,6 @@ module subroutine setup_construct_system(system, param)
allocate(symba_pltpenc :: system%pltpenc_list)
allocate(symba_plplenc :: system%plplenc_list)
allocate(symba_plplenc :: system%plplcollision_list)
allocate(symba_encounter_storage :: system%encounter_history)
end select
case (RINGMOONS)
write(*,*) 'RINGMOONS-SyMBA integrator not yet enabled'
Expand Down
15 changes: 14 additions & 1 deletion src/symba/symba_io.f90
Original file line number Diff line number Diff line change
Expand Up @@ -356,19 +356,27 @@ module subroutine symba_io_param_writer(self, unit, iotype, v_list, iostat, ioms
write(*,*) "Error writing parameter file for SyMBA: " // trim(adjustl(iomsg))
end subroutine symba_io_param_writer


module subroutine symba_io_start_encounter(self, param, t)
!! author: David A. Minton
!!
!! Initializes the new encounter and/or fragmentation save file(s)
!! Initializes the new encounter and/or fragmentation history
implicit none
! Arguments
class(symba_nbody_system), intent(inout) :: self !! SyMBA nbody system object
class(symba_parameters), intent(inout) :: param !! Current run configuration parameters
real(DP), intent(in) :: t !! Current simulation time

if (allocated(self%encounter_history)) deallocate(self%encounter_history)
allocate(symba_encounter_storage :: self%encounter_history)

! Take the snapshot at the start of the encounter
call self%snapshot(param, t)

return
end subroutine symba_io_start_encounter


module subroutine symba_io_stop_encounter(self, param, t)
!! author: David A. Minton
!!
Expand All @@ -379,6 +387,11 @@ module subroutine symba_io_stop_encounter(self, param, t)
class(symba_parameters), intent(inout) :: param !! Current run configuration parameters
real(DP), intent(in) :: t !! Current simulation time

! Take the final snapshot
call self%snapshot(param, t)
call self%encounter_history%dump(param)
deallocate(self%encounter_history)

return
end subroutine symba_io_stop_encounter

Expand Down
129 changes: 66 additions & 63 deletions src/symba/symba_step.f90
Original file line number Diff line number Diff line change
Expand Up @@ -181,71 +181,74 @@ recursive module subroutine symba_step_recur_system(self, param, t, ireci)
logical :: lencounter, lplpl_collision, lpltp_collision

associate(system => self, plplenc_list => self%plplenc_list, pltpenc_list => self%pltpenc_list)
select type(pl => self%pl)
class is (symba_pl)
select type(tp => self%tp)
class is (symba_tp)
system%irec = ireci
dtl = param%dt / (NTENC**ireci)
dth = 0.5_DP * dtl
IF (dtl / param%dt < VSMALL) THEN
write(*, *) "SWIFTEST Warning:"
write(*, *) " In symba_step_recur_system, local time step is too small"
write(*, *) " Roundoff error will be important!"
call util_exit(FAILURE)
END IF
irecp = ireci + 1
if (ireci == 0) then
nloops = 1
else
nloops = NTENC
end if
do j = 1, nloops
lencounter = plplenc_list%encounter_check(param, system, dtl, irecp) &
.or. pltpenc_list%encounter_check(param, system, dtl, irecp)

call plplenc_list%kick(system, dth, irecp, 1)
call pltpenc_list%kick(system, dth, irecp, 1)
if (ireci /= 0) then
call plplenc_list%kick(system, dth, irecp, -1)
call pltpenc_list%kick(system, dth, irecp, -1)
end if

if (param%lgr) then
call pl%gr_pos_kick(system, param, dth)
call tp%gr_pos_kick(system, param, dth)
end if

call pl%drift(system, param, dtl)
call tp%drift(system, param, dtl)

if (lencounter) call system%recursive_step(param, t+dth,irecp)
select type(param)
class is (symba_parameters)
select type(pl => self%pl)
class is (symba_pl)
select type(tp => self%tp)
class is (symba_tp)
system%irec = ireci

if (param%lgr) then
call pl%gr_pos_kick(system, param, dth)
call tp%gr_pos_kick(system, param, dth)
dtl = param%dt / (NTENC**ireci)
dth = 0.5_DP * dtl
IF (dtl / param%dt < VSMALL) THEN
write(*, *) "SWIFTEST Warning:"
write(*, *) " In symba_step_recur_system, local time step is too small"
write(*, *) " Roundoff error will be important!"
call util_exit(FAILURE)
END IF
irecp = ireci + 1
if (ireci == 0) then
nloops = 1
else
nloops = NTENC
end if

call plplenc_list%kick(system, dth, irecp, 1)
call pltpenc_list%kick(system, dth, irecp, 1)
if (ireci /= 0) then
call plplenc_list%kick(system, dth, irecp, -1)
call pltpenc_list%kick(system, dth, irecp, -1)
end if

if (param%lclose) then
lplpl_collision = plplenc_list%collision_check(system, param, t+dtl, dtl, ireci)
lpltp_collision = pltpenc_list%collision_check(system, param, t+dtl, dtl, ireci)

if (lplpl_collision) call plplenc_list%resolve_collision(system, param, t+dtl, dtl, ireci)
if (lpltp_collision) call pltpenc_list%resolve_collision(system, param, t+dtl, dtl, ireci)
end if
call system%snapshot(param, t+dtl)

call self%set_recur_levels(ireci)

end do
do j = 1, nloops
lencounter = plplenc_list%encounter_check(param, system, dtl, irecp) &
.or. pltpenc_list%encounter_check(param, system, dtl, irecp)

call plplenc_list%kick(system, dth, irecp, 1)
call pltpenc_list%kick(system, dth, irecp, 1)
if (ireci /= 0) then
call plplenc_list%kick(system, dth, irecp, -1)
call pltpenc_list%kick(system, dth, irecp, -1)
end if

if (param%lgr) then
call pl%gr_pos_kick(system, param, dth)
call tp%gr_pos_kick(system, param, dth)
end if

call pl%drift(system, param, dtl)
call tp%drift(system, param, dtl)

if (lencounter) call system%recursive_step(param, t+dth,irecp)
system%irec = ireci

if (param%lgr) then
call pl%gr_pos_kick(system, param, dth)
call tp%gr_pos_kick(system, param, dth)
end if

call plplenc_list%kick(system, dth, irecp, 1)
call pltpenc_list%kick(system, dth, irecp, 1)
if (ireci /= 0) then
call plplenc_list%kick(system, dth, irecp, -1)
call pltpenc_list%kick(system, dth, irecp, -1)
end if

if (param%lclose) then
lplpl_collision = plplenc_list%collision_check(system, param, t+dtl, dtl, ireci)
lpltp_collision = pltpenc_list%collision_check(system, param, t+dtl, dtl, ireci)

if (lplpl_collision) call plplenc_list%resolve_collision(system, param, t+dtl, dtl, ireci)
if (lpltp_collision) call pltpenc_list%resolve_collision(system, param, t+dtl, dtl, ireci)
end if
if (param%lencounter_save) call system%snapshot(param, t+dtl)

call self%set_recur_levels(ireci)

end do
end select
end select
end select
end associate
Expand Down

0 comments on commit 155dc89

Please sign in to comment.