diff --git a/src/base/base_module.f90 b/src/base/base_module.f90 index 0610ece9e..b55fe3d98 100644 --- a/src/base/base_module.f90 +++ b/src/base/base_module.f90 @@ -29,6 +29,7 @@ module base real(DP) :: tstop = -1.0_DP !! Integration stop time real(DP) :: dt = -1.0_DP !! Time step integer(I8B) :: iloop = 0_I8B !! Main loop counter + integer(I8B) :: nloops = 0_I8B !! Total number of loops to execute character(STRMAX) :: incbfile = CB_INFILE !! Name of input file for the central body character(STRMAX) :: inplfile = PL_INFILE !! Name of input file for massive bodies character(STRMAX) :: intpfile = TP_INFILE !! Name of input file for test particles diff --git a/src/swiftest/swiftest_driver.f90 b/src/swiftest/swiftest_driver.f90 index d4ed0cc9a..0e9f369df 100644 --- a/src/swiftest/swiftest_driver.f90 +++ b/src/swiftest/swiftest_driver.f90 @@ -16,30 +16,17 @@ program swiftest_driver !! Adapted from Swifter by David E. Kaufmann's Swifter driver programs swifter_[bs,helio,ra15,rmvs,symba,tu4,whm].f90 !! Adapted from Hal Levison and Martin Duncan's Swift driver programs use swiftest - use symba implicit none - class(swiftest_nbody_system), allocatable :: nbody_system !! Polymorphic object containing the nbody system to be integrated + class(swiftest_nbody_system), allocatable :: nbody_system !! Polymorphic object containing the nbody system to be integrated class(swiftest_parameters), allocatable :: param !! Run configuration parameters character(len=:), allocatable :: integrator !! Integrator type code (see globals for symbolic names) character(len=:), allocatable :: param_file_name !! Name of the file containing user-defined parameters character(len=:), allocatable :: display_style !! Style of the output display {"STANDARD", "COMPACT", "PROGRESS"}). Default is "STANDARD" integer(I8B) :: istart !! Starting index for loop counter - integer(I8B) :: nloops !! Number of steps to take in the simulation integer(I4B) :: iout !! Output cadence counter integer(I4B) :: idump !! Dump cadence counter type(walltimer) :: integration_timer !! Object used for computing elapsed wall time - real(DP) :: tfrac !! Fraction of total simulation time completed - type(progress_bar) :: pbar !! Object used to print out a progress bar - character(*), parameter :: statusfmt = '("Time = ", ES12.5, "; fraction done = ", F6.3, ' // & - '"; Number of active pl, tp = ", I6, ", ", I6)' - character(*), parameter :: symbastatfmt = '("Time = ", ES12.5, "; fraction done = ", F6.3, ' // & - '"; Number of active plm, pl, tp = ", I6, ", ", I6, ", ", I6)' - character(*), parameter :: pbarfmt = '("Time = ", ES12.5," of ",ES12.5)' - character(len=64) :: pbarmessage - - character(*), parameter :: symbacompactfmt = '(";NPLM",ES22.15,$)' - !type(base_storage(nframes=:)), allocatable :: system_history call swiftest_io_get_args(integrator, param_file_name, display_style) @@ -54,9 +41,9 @@ program swiftest_driver dt => param%dt, & tstop => param%tstop, & iloop => param%iloop, & + nloops => param%nloops, & istep_out => param%istep_out, & dump_cadence => param%dump_cadence, & - display_style => param%display_style, & display_unit => param%display_unit) ! Set up loop and output cadence variables @@ -84,29 +71,18 @@ program swiftest_driver associate (system_history => param%system_history) ! If this is a new run, compute energy initial conditions (if energy tracking is turned on) and write the initial conditions to file. + call nbody_system%display_run_information(param, integration_timer, phase="FIRST") if (param%lenergy) then if (param%lrestart) then call nbody_system%get_t0_values(param) - call nbody_system%conservation_report(param, lterminal=.true.) else call nbody_system%conservation_report(param, lterminal=.false.) ! This will save the initial values of energy and momentum end if end if + call nbody_system%conservation_report(param, lterminal=.true.) call system_history%take_snapshot(param,nbody_system) call nbody_system%dump(param) - write(display_unit, *) " *************** Main Loop *************** " - - if (display_style == "PROGRESS") then - call pbar%reset(nloops) - write(pbarmessage,fmt=pbarfmt) t0, tstop - call pbar%update(1,message=pbarmessage) - else if (display_style == "COMPACT") then - write(*,*) "SWIFTEST START " // param%integrator - call nbody_system%compact_output(param,integration_timer) - end if - - do iloop = istart, nloops !> Step the nbody_system forward in time call integration_timer%start() @@ -132,25 +108,10 @@ program swiftest_driver end if - tfrac = (nbody_system%t - t0) / (tstop - t0) - - select type(pl => nbody_system%pl) - class is (symba_pl) - write(display_unit, symbastatfmt) nbody_system%t, tfrac, pl%nplm, pl%nbody, nbody_system%tp%nbody - class default - write(display_unit, statusfmt) nbody_system%t, tfrac, pl%nbody, nbody_system%tp%nbody - end select - if (param%lenergy) call nbody_system%conservation_report(param, lterminal=.true.) call integration_timer%report(message="Integration steps:", unit=display_unit, nsubsteps=istep_out) - - if (display_style == "PROGRESS") then - write(pbarmessage,fmt=pbarfmt) nbody_system%t, tstop - call pbar%update(1,message=pbarmessage) - else if (display_style == "COMPACT") then - call nbody_system%compact_output(param,integration_timer) - end if - + call nbody_system%display_run_information(param, integration_timer) call integration_timer%reset() + if (param%lenergy) call nbody_system%conservation_report(param, lterminal=.true.) end if end if @@ -159,7 +120,8 @@ program swiftest_driver ! Dump any remaining history if it exists call nbody_system%dump(param) call system_history%dump(param) - if (display_style == "COMPACT") write(*,*) "SWIFTEST STOP" // param%integrator + call nbody_system%display_run_information(param, integration_timer, phase="LAST") + end associate end associate diff --git a/src/swiftest/swiftest_io.f90 b/src/swiftest/swiftest_io.f90 index a9a9a873d..e49bdd3f0 100644 --- a/src/swiftest/swiftest_io.f90 +++ b/src/swiftest/swiftest_io.f90 @@ -194,6 +194,74 @@ module subroutine swiftest_io_conservation_report(self, param, lterminal) end subroutine swiftest_io_conservation_report + module subroutine swiftest_io_display_run_information(self, param, integration_timer, phase) + !! author: David A. Minton + !! + !! Displays helpful information while a run is executing. The format of the output depends on user parameters + implicit none + ! Arguments + class(swiftest_nbody_system), intent(inout) :: self !! Swiftest nbody system object + class(swiftest_parameters), intent(inout) :: param !! Current run configuration parameters + type(walltimer), intent(inout) :: integration_timer !! Object used for computing elapsed wall time + character(len=*), optional, intent(in) :: phase !! One of "FIRST" or "LAST" + ! Internals + integer(I4B) :: phase_val + real(DP) :: tfrac !! Fraction of total simulation time completed + type(progress_bar), save :: pbar !! Object used to print out a progress bar + character(len=64) :: pbarmessage + character(*), parameter :: symbacompactfmt = '(";NPLM",ES22.15,$)' + character(*), parameter :: statusfmt = '("Time = ", ES12.5, "; fraction done = ", F6.3, ' // & + '"; Number of active pl, tp = ", I6, ", ", I6)' + character(*), parameter :: symbastatfmt = '("Time = ", ES12.5, "; fraction done = ", F6.3, ' // & + '"; Number of active plm, pl, tp = ", I6, ", ", I6, ", ", I6)' + character(*), parameter :: pbarfmt = '("Time = ", ES12.5," of ",ES12.5)' + + phase_val = 1 + if (present(phase)) then + if (phase == "FIRST") then + phase_val = 0 + else if (phase == "LAST") then + phase_val = -1 + end if + end if + + tfrac = (self%t - param%t0) / (param%tstop - param%t0) + + if (phase_val == 0) then + if (param%lrestart) then + write(param%display_unit, *) " *************** Swiftest restart " // param%integrator // " *************** " + else + write(param%display_unit, *) " *************** Swiftest start " // param%integrator // " *************** " + end if + if (param%display_style == "PROGRESS") then + call pbar%reset(param%nloops) + else if (param%display_style == "COMPACT") then + write(*,*) "SWIFTEST START " // param%integrator + end if + end if + + if (param%display_style == "PROGRESS") then + write(pbarmessage,fmt=pbarfmt) self%t, param%tstop + call pbar%update(1,message=pbarmessage) + else if (param%display_style == "COMPACT") then + call self%compact_output(param,integration_timer) + end if + + if (self%pl%nplm > 0) then + write(param%display_unit, symbastatfmt) self%t, tfrac, self%pl%nplm, self%pl%nbody, self%tp%nbody + else + write(param%display_unit, statusfmt) self%t, tfrac, self%pl%nbody, self%tp%nbody + end if + + if (phase_val == -1) then + write(param%display_unit, *)" *************** Swiftest stop " // param%integrator // " *************** " + if (param%display_style == "COMPACT") write(*,*) "SWIFTEST STOP" // param%integrator + end if + + return + end subroutine swiftest_io_display_run_information + + module subroutine swiftest_io_dump_param(self, param_file_name) !! author: David A. Minton !! @@ -504,56 +572,54 @@ module subroutine swiftest_io_netcdf_get_t0_values_system(self, param) real(DP), dimension(NDIM) :: rot0, Ip0, Lnow real(DP) :: KE_orb_orig, KE_spin_orig, PE_orig, BE_orig - associate (nc => param%system_history%nc, cb => self%cb) + associate (nc => param%system_history%nc, cb => self%cb, tslot => param%system_history%nc%tslot) call nc%open(param, readonly=.true.) + call nc%find_tslot(param%t0) call netcdf_io_check( nf90_inquire_dimension(nc%id, nc%time_dimid, len=itmax), "netcdf_io_get_t0_values_system time_dimid" ) call netcdf_io_check( nf90_inquire_dimension(nc%id, nc%name_dimid, len=idmax), "netcdf_io_get_t0_values_system name_dimid" ) allocate(vals(idmax)) - call netcdf_io_check( nf90_get_var(nc%id, nc%time_varid, rtemp, start=[1], count=[1]), "netcdf_io_get_t0_values_system time_varid" ) + call netcdf_io_check( nf90_get_var(nc%id, nc%time_varid, rtemp, start=[tslot], count=[1]), "netcdf_io_get_t0_values_system time_varid" ) if (param%lenergy) then - call netcdf_io_check( nf90_get_var(nc%id, nc%KE_orb_varid, rtemp, start=[1], count=[1]), "netcdf_io_get_t0_values_system KE_orb_varid" ) + call netcdf_io_check( nf90_get_var(nc%id, nc%KE_orb_varid, rtemp, start=[tslot], count=[1]), "netcdf_io_get_t0_values_system KE_orb_varid" ) KE_orb_orig = rtemp(1) - call netcdf_io_check( nf90_get_var(nc%id, nc%KE_spin_varid, rtemp, start=[1], count=[1]), "netcdf_io_get_t0_values_system KE_spin_varid" ) + call netcdf_io_check( nf90_get_var(nc%id, nc%KE_spin_varid, rtemp, start=[tslot], count=[1]), "netcdf_io_get_t0_values_system KE_spin_varid" ) KE_spin_orig = rtemp(1) - call netcdf_io_check( nf90_get_var(nc%id, nc%PE_varid, rtemp, start=[1], count=[1]), "netcdf_io_get_t0_values_system PE_varid" ) + call netcdf_io_check( nf90_get_var(nc%id, nc%PE_varid, rtemp, start=[tslot], count=[1]), "netcdf_io_get_t0_values_system PE_varid" ) PE_orig = rtemp(1) - call netcdf_io_check( nf90_get_var(nc%id, nc%BE_varid, rtemp, start=[1], count=[1]), "netcdf_io_get_t0_values_system BE_varid" ) + call netcdf_io_check( nf90_get_var(nc%id, nc%BE_varid, rtemp, start=[tslot], count=[1]), "netcdf_io_get_t0_values_system BE_varid" ) BE_orig = rtemp(1) - call netcdf_io_check( nf90_get_var(nc%id, nc%E_collisions_varid, self%E_collisions, start=[1]), "netcdf_io_get_t0_values_system E_collisions_varid" ) - call netcdf_io_check( nf90_get_var(nc%id, nc%E_untracked_varid, self%E_untracked, start=[1]), "netcdf_io_get_t0_values_system E_untracked_varid" ) + call netcdf_io_check( nf90_get_var(nc%id, nc%E_collisions_varid, self%E_collisions, start=[tslot]), "netcdf_io_get_t0_values_system E_collisions_varid" ) + call netcdf_io_check( nf90_get_var(nc%id, nc%E_untracked_varid, self%E_untracked, start=[tslot]), "netcdf_io_get_t0_values_system E_untracked_varid" ) self%E_orbit_orig = KE_orb_orig + KE_spin_orig + PE_orig + BE_orig + self%E_collisions + self%E_untracked - call netcdf_io_check( nf90_get_var(nc%id, nc%L_orbit_varid, self%L_orbit_orig(:), start=[1,1], count=[NDIM,1]), "netcdf_io_get_t0_values_system L_orbit_varid" ) - call netcdf_io_check( nf90_get_var(nc%id, nc%L_spin_varid, self%L_spin_orig(:), start=[1,1], count=[NDIM,1]), "netcdf_io_get_t0_values_system L_spin_varid" ) - call netcdf_io_check( nf90_get_var(nc%id, nc%L_escape_varid, self%L_escape(:), start=[1,1], count=[NDIM,1]), "netcdf_io_get_t0_values_system L_escape_varid" ) + call netcdf_io_check( nf90_get_var(nc%id, nc%L_orbit_varid, self%L_orbit_orig(:), start=[1,tslot], count=[NDIM,1]), "netcdf_io_get_t0_values_system L_orbit_varid" ) + call netcdf_io_check( nf90_get_var(nc%id, nc%L_spin_varid, self%L_spin_orig(:), start=[1,tslot], count=[NDIM,1]), "netcdf_io_get_t0_values_system L_spin_varid" ) + call netcdf_io_check( nf90_get_var(nc%id, nc%L_escape_varid, self%L_escape(:), start=[1,tslot], count=[NDIM,1]), "netcdf_io_get_t0_values_system L_escape_varid" ) self%L_total_orig(:) = self%L_orbit_orig(:) + self%L_spin_orig(:) + self%L_escape(:) - call netcdf_io_check( nf90_get_var(nc%id, nc%Gmass_varid, vals, start=[1,1], count=[idmax,1]), "netcdf_io_get_t0_values_system Gmass_varid" ) - call netcdf_io_check( nf90_get_var(nc%id, nc%GMescape_varid, self%GMescape, start=[1]), "netcdf_io_get_t0_values_system GMescape_varid" ) + call netcdf_io_check( nf90_get_var(nc%id, nc%Gmass_varid, vals, start=[1,tslot], count=[idmax,1]), "netcdf_io_get_t0_values_system Gmass_varid" ) + call netcdf_io_check( nf90_get_var(nc%id, nc%GMescape_varid, self%GMescape, start=[tslot]), "netcdf_io_get_t0_values_system GMescape_varid" ) self%GMtot_orig = vals(1) + sum(vals(2:idmax), vals(2:idmax) == vals(2:idmax)) + self%GMescape cb%GM0 = vals(1) cb%dGM = cb%Gmass - cb%GM0 - call netcdf_io_check( nf90_get_var(nc%id, nc%radius_varid, rtemp, start=[1,1], count=[1,1]), "netcdf_io_get_t0_values_system radius_varid" ) + call netcdf_io_check( nf90_get_var(nc%id, nc%radius_varid, rtemp, start=[1,tslot], count=[1,1]), "netcdf_io_get_t0_values_system radius_varid" ) cb%R0 = rtemp(1) if (param%lrotation) then - call netcdf_io_check( nf90_get_var(nc%id, nc%rot_varid, rot0, start=[1,1,1], count=[NDIM,1,1]), "netcdf_io_get_t0_values_system rot_varid" ) - call netcdf_io_check( nf90_get_var(nc%id, nc%Ip_varid, Ip0, start=[1,1,1], count=[NDIM,1,1]), "netcdf_io_get_t0_values_system Ip_varid" ) + call netcdf_io_check( nf90_get_var(nc%id, nc%rot_varid, rot0, start=[1,1,tslot], count=[NDIM,1,1]), "netcdf_io_get_t0_values_system rot_varid" ) + call netcdf_io_check( nf90_get_var(nc%id, nc%Ip_varid, Ip0, start=[1,1,tslot], count=[NDIM,1,1]), "netcdf_io_get_t0_values_system Ip_varid" ) cb%L0(:) = Ip0(3) * cb%GM0 * cb%R0**2 * rot0(:) - - Lnow(:) = cb%Ip(3) * cb%Gmass * cb%radius**2 * cb%rot(:) - cb%dL(:) = Lnow(:) - cb%L0(:) end if end if @@ -2198,7 +2264,7 @@ module subroutine swiftest_io_param_reader(self, unit, iotype, v_list, iostat, i call param%set_display(param%display_style) ! Print the contents of the parameter file to standard output - call param%writer(unit = param%display_unit, iotype = "none", v_list = [0], iostat = iostat, iomsg = iomsg) + if (.not.param%lrestart) call param%writer(unit = param%display_unit, iotype = "none", v_list = [0], iostat = iostat, iomsg = iomsg) end associate diff --git a/src/swiftest/swiftest_module.f90 b/src/swiftest/swiftest_module.f90 index c38470b69..ce74907b5 100644 --- a/src/swiftest/swiftest_module.f90 +++ b/src/swiftest/swiftest_module.f90 @@ -374,19 +374,20 @@ module swiftest procedure(abstract_step_system), deferred :: step ! Concrete classes that are common to the basic integrator (only test particles considered for discard) - procedure :: discard => swiftest_discard_system !! Perform a discard step on the nbody_system - procedure :: compact_output => swiftest_io_compact_output !! Prints out out terminal output when display_style is set to COMPACT - procedure :: conservation_report => swiftest_io_conservation_report !! Compute energy and momentum and print out the change with time - procedure :: dump => swiftest_io_dump_system !! Dump the state of the nbody_system to a file - procedure :: get_t0_values => swiftest_io_netcdf_get_t0_values_system !! Validates the dump file to check whether the dump file initial conditions duplicate the last frame of the netcdf output. - procedure :: read_frame => swiftest_io_netcdf_read_frame_system !! Read in a frame of input data from file - procedure :: write_frame_netcdf => swiftest_io_netcdf_write_frame_system !! Write a frame of input data from file - procedure :: write_frame_system => swiftest_io_write_frame_system !! Write a frame of input data from file - procedure :: read_hdr => swiftest_io_netcdf_read_hdr_system !! Read a header for an output frame in NetCDF format - procedure :: write_hdr => swiftest_io_netcdf_write_hdr_system !! Write a header for an output frame in NetCDF format - procedure :: read_in => swiftest_io_read_in_system !! Reads the initial conditions for an nbody system - procedure :: read_particle_info => swiftest_io_netcdf_read_particle_info_system !! Read in particle metadata from file - procedure :: obl_pot => swiftest_obl_pot_system !! Compute the contribution to the total gravitational potential due solely to the oblateness of the central body + procedure :: discard => swiftest_discard_system !! Perform a discard step on the nbody_system + procedure :: compact_output => swiftest_io_compact_output !! Prints out out terminal output when display_style is set to COMPACT + procedure :: conservation_report => swiftest_io_conservation_report !! Compute energy and momentum and print out the change with time + procedure :: display_run_information => swiftest_io_display_run_information !! Displays helpful information about the run + procedure :: dump => swiftest_io_dump_system !! Dump the state of the nbody_system to a file + procedure :: get_t0_values => swiftest_io_netcdf_get_t0_values_system !! Validates the dump file to check whether the dump file initial conditions duplicate the last frame of the netcdf output. + procedure :: read_frame => swiftest_io_netcdf_read_frame_system !! Read in a frame of input data from file + procedure :: write_frame_netcdf => swiftest_io_netcdf_write_frame_system !! Write a frame of input data from file + procedure :: write_frame_system => swiftest_io_write_frame_system !! Write a frame of input data from file + procedure :: read_hdr => swiftest_io_netcdf_read_hdr_system !! Read a header for an output frame in NetCDF format + procedure :: write_hdr => swiftest_io_netcdf_write_hdr_system !! Write a header for an output frame in NetCDF format + procedure :: read_in => swiftest_io_read_in_system !! Reads the initial conditions for an nbody system + procedure :: read_particle_info => swiftest_io_netcdf_read_particle_info_system !! Read in particle metadata from file + procedure :: obl_pot => swiftest_obl_pot_system !! Compute the contribution to the total gravitational potential due solely to the oblateness of the central body procedure :: initialize => swiftest_util_setup_initialize_system !! Initialize the nbody_system from input files procedure :: init_particle_info => swiftest_util_setup_initialize_particle_info_system !! Initialize the nbody_system from input files ! procedure :: step_spin => tides_step_spin_system !! Steps the spins of the massive & central bodies due to tides. @@ -574,6 +575,14 @@ module subroutine swiftest_io_conservation_report(self, param, lterminal) logical, intent(in) :: lterminal !! Indicates whether to output information to the terminal screen end subroutine swiftest_io_conservation_report + module subroutine swiftest_io_display_run_information(self, param, integration_timer, phase) + implicit none + class(swiftest_nbody_system), intent(inout) :: self !! Swiftest nbody system object + class(swiftest_parameters), intent(inout) :: param !! Current run configuration parameters + type(walltimer), intent(inout) :: integration_timer !! Object used for computing elapsed wall time + character(len=*), optional, intent(in) :: phase !! One of "FIRST" or "LAST" + end subroutine swiftest_io_display_run_information + module subroutine swiftest_io_dump_param(self, param_file_name) implicit none class(swiftest_parameters),intent(in) :: self !! Output collection of parameters