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

Commit

Permalink
Merge branch 'encounter_storage' into debug
Browse files Browse the repository at this point in the history
  • Loading branch information
daminton committed Dec 14, 2022
2 parents 9150e0f + d556d8a commit 80c9b92
Show file tree
Hide file tree
Showing 7 changed files with 98 additions and 143 deletions.
10 changes: 5 additions & 5 deletions examples/Fragmentation/Fragmentation_Movie.py
Original file line number Diff line number Diff line change
Expand Up @@ -54,7 +54,7 @@
"supercatastrophic_off_axis": [np.array([0.0, 6.28, 0.0]),
np.array([0.5, -6.28, 0.0])],
"hitandrun" : [np.array([0.0, 6.28, 0.0]),
np.array([-0.1, -6.28, 0.0])]
np.array([-0.9, -6.28, 0.0])]
}

rot_vectors = {"disruption_headon" : [np.array([0.0, 0.0, 0.0]),
Expand All @@ -67,7 +67,7 @@

body_Gmass = {"disruption_headon" : [1e-7, 1e-10],
"supercatastrophic_off_axis": [1e-7, 1e-8],
"hitandrun" : [1e-7, 7e-10]
"hitandrun" : [1e-7, 1e-10]
}

density = 3000 * swiftest.AU2M**3 / swiftest.MSun
Expand Down Expand Up @@ -139,7 +139,7 @@ def setup_plot(self):

scale_frame = abs(rhy1) + abs(rhy2)
ax = plt.Axes(fig, [0.1, 0.1, 0.8, 0.8])
self.ax_pt_size = self.figsize[0] * 0.8 * 72 / (np.sqrt(2)*scale_frame)
self.ax_pt_size = self.figsize[0] * 0.7 * 72 / scale_frame
ax.set_xlim(-scale_frame, scale_frame)
ax.set_ylim(-scale_frame, scale_frame)
ax.set_xticks([])
Expand Down Expand Up @@ -203,8 +203,8 @@ def data_stream(self, frame=0):
# Set fragmentation parameters
minimum_fragment_gmass = 0.2 * body_Gmass[style][1] # Make the minimum fragment mass a fraction of the smallest body
gmtiny = 0.99 * body_Gmass[style][1] # Make GMTINY just smaller than the smallest original body. This will prevent runaway collisional cascades
sim.set_parameter(fragmentation=True, encounter_save="trajectory", gmtiny=gmtiny, minimum_fragment_gmass=minimum_fragment_gmass, verbose=False)
sim.run(dt=1e-3, tstop=1.0e-3, istep_out=1, dump_cadence=1)
sim.set_parameter(fragmentation=True, encounter_save="both", gmtiny=gmtiny, minimum_fragment_gmass=minimum_fragment_gmass, verbose=False)
sim.run(dt=1e-3, tstop=1.0e-3, istep_out=1, dump_cadence=0)

print("Generating animation")
anim = AnimatedScatter(sim,movie_filename,movie_titles[style],style,nskip=1)
5 changes: 4 additions & 1 deletion python/swiftest/swiftest/simulation_class.py
Original file line number Diff line number Diff line change
Expand Up @@ -2713,7 +2713,8 @@ def read_output_file(self,read_init_cond : bool = True):
param_tmp['BIN_OUT'] = os.path.join(self.simdir, self.param['BIN_OUT'])
if self.codename == "Swiftest":
self.data = io.swiftest2xr(param_tmp, verbose=self.verbose)
if self.verbose: print('Swiftest simulation data stored as xarray DataSet .data')
if self.verbose:
print('Swiftest simulation data stored as xarray DataSet .data')
if read_init_cond:
if self.verbose:
print("Reading initial conditions file as .init_cond")
Expand All @@ -2725,6 +2726,8 @@ def read_output_file(self,read_init_cond : bool = True):

self.read_encounters()
self.read_collisions()
if self.verbose:
print("Finished reading Swiftest dataset files.")

elif self.codename == "Swifter":
self.data = io.swifter2xr(param_tmp, verbose=self.verbose)
Expand Down
63 changes: 25 additions & 38 deletions src/fraggle/fraggle_io.f90
Original file line number Diff line number Diff line change
Expand Up @@ -98,23 +98,24 @@ module subroutine fraggle_io_initialize_output(self, param)

call check( nf90_def_var(nc%id, nc%rot_varname, nc%out_type,&
[ nc%space_dimid, nc%name_dimid, nc%stage_dimid, nc%event_dimid], nc%rot_varid), "fraggle_io_initialize nf90_def_var rot_varid")

if (param%lenergy) then

call check( nf90_def_var(nc%id, nc%ke_orb_varname, nc%out_type,&
[ nc%stage_dimid, nc%event_dimid], nc%KE_orb_varid), "fraggle_io_initialize_output nf90_def_var KE_orb_varid")

call check( nf90_def_var(nc%id, nc%ke_spin_varname, nc%out_type,&
[ nc%stage_dimid, nc%event_dimid], nc%KE_spin_varid), "fraggle_io_initialize_output nf90_def_var KE_spin_varid" )

call check( nf90_def_var(nc%id, nc%pe_varname, nc%out_type,&
[ nc%stage_dimid, nc%event_dimid], nc%PE_varid), "fraggle_io_initialize_output nf90_def_var PE_varid" )
call check( nf90_def_var(nc%id, nc%ke_orb_varname, nc%out_type,&
[ nc%stage_dimid, nc%event_dimid], nc%KE_orb_varid), "fraggle_io_initialize_output nf90_def_var KE_orb_varid")

call check( nf90_def_var(nc%id, nc%L_orb_varname, nc%out_type, &
[ nc%space_dimid, nc%stage_dimid, nc%event_dimid], nc%L_orb_varid), "fraggle_io_initialize_output nf90_def_var L_orb_varid" )
call check( nf90_def_var(nc%id, nc%ke_spin_varname, nc%out_type,&
[ nc%stage_dimid, nc%event_dimid], nc%KE_spin_varid), "fraggle_io_initialize_output nf90_def_var KE_spin_varid" )

call check( nf90_def_var(nc%id, nc%L_spin_varname, nc%out_type,&
[ nc%space_dimid, nc%stage_dimid, nc%event_dimid], nc%L_spin_varid), "fraggle_io_initialize_output nf90_def_var L_spin_varid" )
call check( nf90_def_var(nc%id, nc%pe_varname, nc%out_type,&
[ nc%stage_dimid, nc%event_dimid], nc%PE_varid), "fraggle_io_initialize_output nf90_def_var PE_varid" )

call check( nf90_def_var(nc%id, nc%L_orb_varname, nc%out_type, &
[ nc%space_dimid, nc%stage_dimid, nc%event_dimid], nc%L_orb_varid), "fraggle_io_initialize_output nf90_def_var L_orb_varid" )

call check( nf90_def_var(nc%id, nc%L_spin_varname, nc%out_type,&
[ nc%space_dimid, nc%stage_dimid, nc%event_dimid], nc%L_spin_varid), "fraggle_io_initialize_output nf90_def_var L_spin_varid" )
end if

call check( nf90_inquire(nc%id, nVariables=nvar), "fraggle_io_initialize nf90_inquire nVariables" )
do varid = 1, nvar
Expand Down Expand Up @@ -202,14 +203,18 @@ module subroutine fraggle_io_write_frame(self, nc, param)
call check( nf90_put_var(nc%id, nc%rot_varid, pl%rot(:,i), start=[1, idslot, stage, eslot], count=[NDIM,1,1,1]), "fraggle_io_write_frame nf90_put_var rotx_varid" )
end do
end do
call check( nf90_put_var(nc%id, nc%ke_orb_varid, fragments%ke_orbit_before, start=[ 1, eslot]), "fraggle_io_write_frame nf90_put_var ke_orb_varid before" )
call check( nf90_put_var(nc%id, nc%ke_orb_varid, fragments%ke_orbit_after, start=[ 2, eslot]), "fraggle_io_write_frame nf90_put_var ke_orb_varid after" )
call check( nf90_put_var(nc%id, nc%pe_varid, fragments%pe_before, start=[ 1, eslot]), "fraggle_io_write_frame nf90_put_var pe_varid before" )
call check( nf90_put_var(nc%id, nc%pe_varid, fragments%pe_after, start=[ 2, eslot]), "fraggle_io_write_frame nf90_put_var pe_varid after" )
call check( nf90_put_var(nc%id, nc%L_orb_varid, fragments%Lorbit_before(:), start=[1, 1, eslot], count=[NDIM, 1, 1]), "fraggle_io_write_frame nf90_put_var L_orb_varid before" )
call check( nf90_put_var(nc%id, nc%L_orb_varid, fragments%Lorbit_after(:), start=[1, 2, eslot], count=[NDIM, 1, 1]), "fraggle_io_write_frame nf90_put_var L_orb_varid after" )
call check( nf90_put_var(nc%id, nc%L_spin_varid, fragments%Lspin_before(:), start=[1, 1, eslot], count=[NDIM, 1, 1]), "fraggle_io_write_frame nf90_put_var L_spin_varid before" )
call check( nf90_put_var(nc%id, nc%L_spin_varid, fragments%Lspin_after(:), start=[1, 2, eslot], count=[NDIM, 1, 1]), "fraggle_io_write_frame nf90_put_var L_spin_varid after" )
if (param%lenergy) then
call check( nf90_put_var(nc%id, nc%ke_orb_varid, fragments%ke_orbit_before, start=[ 1, eslot]), "fraggle_io_write_frame nf90_put_var ke_orb_varid before" )
call check( nf90_put_var(nc%id, nc%ke_orb_varid, fragments%ke_orbit_after, start=[ 2, eslot]), "fraggle_io_write_frame nf90_put_var ke_orb_varid after" )
call check( nf90_put_var(nc%id, nc%ke_spin_varid, fragments%ke_spin_before, start=[ 1, eslot]), "fraggle_io_write_frame nf90_put_var ke_spin_varid before" )
call check( nf90_put_var(nc%id, nc%ke_spin_varid, fragments%ke_spin_after, start=[ 2, eslot]), "fraggle_io_write_frame nf90_put_var ke_spin_varid after" )
call check( nf90_put_var(nc%id, nc%pe_varid, fragments%pe_before, start=[ 1, eslot]), "fraggle_io_write_frame nf90_put_var pe_varid before" )
call check( nf90_put_var(nc%id, nc%pe_varid, fragments%pe_after, start=[ 2, eslot]), "fraggle_io_write_frame nf90_put_var pe_varid after" )
call check( nf90_put_var(nc%id, nc%L_orb_varid, fragments%Lorbit_before(:), start=[1, 1, eslot], count=[NDIM, 1, 1]), "fraggle_io_write_frame nf90_put_var L_orb_varid before" )
call check( nf90_put_var(nc%id, nc%L_orb_varid, fragments%Lorbit_after(:), start=[1, 2, eslot], count=[NDIM, 1, 1]), "fraggle_io_write_frame nf90_put_var L_orb_varid after" )
call check( nf90_put_var(nc%id, nc%L_spin_varid, fragments%Lspin_before(:), start=[1, 1, eslot], count=[NDIM, 1, 1]), "fraggle_io_write_frame nf90_put_var L_spin_varid before" )
call check( nf90_put_var(nc%id, nc%L_spin_varid, fragments%Lspin_after(:), start=[1, 2, eslot], count=[NDIM, 1, 1]), "fraggle_io_write_frame nf90_put_var L_spin_varid after" )
end if

call check( nf90_set_fill(nc%id, old_mode, old_mode) )
end associate
Expand All @@ -219,24 +224,6 @@ module subroutine fraggle_io_write_frame(self, nc, param)
end subroutine fraggle_io_write_frame


module subroutine fraggle_io_log_pl(pl, param)
!! author: David A. Minton
!!
!! Writes a single message to the fraggle log file
implicit none
! Arguments
class(swiftest_pl), intent(in) :: pl !! Swiftest massive body object (only the new bodies generated in a collision)
class(swiftest_parameters), intent(in) :: param !! Current swiftest run configuration parameters
! Internals
integer(I4B) :: i
character(STRMAX) :: errmsg

return
667 continue
write(*,*) "Error writing Fraggle message to log file: " // trim(adjustl(errmsg))
end subroutine fraggle_io_log_pl


module subroutine fraggle_io_log_regime(colliders, frag)
!! author: David A. Minton
!!
Expand Down
19 changes: 18 additions & 1 deletion src/fraggle/fraggle_set.f90
Original file line number Diff line number Diff line change
Expand Up @@ -287,7 +287,6 @@ module subroutine fraggle_set_original_scale_factors(self, colliders)
do i = 1, 2
colliders%rot(:,i) = colliders%L_spin(:,i) * (colliders%mass(i) * colliders%radius(i)**2 * colliders%Ip(3, i))
end do
frag%Qloss = frag%Qloss * frag%Escale

frag%mtot = frag%mtot * frag%mscale
frag%mass = frag%mass * frag%mscale
Expand All @@ -300,6 +299,24 @@ module subroutine fraggle_set_original_scale_factors(self, colliders)
frag%rb(:, i) = frag%x_coll(:, i) + frag%rbcom(:)
frag%vb(:, i) = frag%v_coll(:, i) + frag%vbcom(:)
end do

frag%Qloss = frag%Qloss * frag%Escale

frag%Lorbit_before(:) = frag%Lorbit_before * frag%Lscale
frag%Lspin_before(:) = frag%Lspin_before * frag%Lscale
frag%Ltot_before(:) = frag%Ltot_before * frag%Lscale
frag%ke_orbit_before = frag%ke_orbit_before * frag%Escale
frag%ke_spin_before = frag%ke_spin_before * frag%Escale
frag%pe_before = frag%pe_before * frag%Escale
frag%Etot_before = frag%Etot_before * frag%Escale

frag%Lorbit_after(:) = frag%Lorbit_after * frag%Lscale
frag%Lspin_after(:) = frag%Lspin_after * frag%Lscale
frag%Ltot_after(:) = frag%Ltot_after * frag%Lscale
frag%ke_orbit_after = frag%ke_orbit_after * frag%Escale
frag%ke_spin_after = frag%ke_spin_after * frag%Escale
frag%pe_after = frag%pe_after * frag%Escale
frag%Etot_after = frag%Etot_after * frag%Escale

frag%mscale = 1.0_DP
frag%dscale = 1.0_DP
Expand Down
18 changes: 6 additions & 12 deletions src/modules/fraggle_classes.f90
Original file line number Diff line number Diff line change
Expand Up @@ -86,12 +86,12 @@ module fraggle_classes
real(DP) :: Etot_before, Etot_after !! Before/after total system energy

! Scale factors used to scale dimensioned quantities to a more "natural" system where important quantities (like kinetic energy, momentum) are of order ~1
real(DP) :: dscale !! Distance dimension scale factor
real(DP) :: mscale !! Mass scale factor
real(DP) :: tscale !! Time scale factor
real(DP) :: vscale !! Velocity scale factor (a convenience unit that is derived from dscale and tscale)
real(DP) :: Escale !! Energy scale factor (a convenience unit that is derived from dscale, tscale, and mscale)
real(DP) :: Lscale !! Angular momentum scale factor (a convenience unit that is derived from dscale, tscale, and mscale)
real(DP) :: dscale = 1.0_DP !! Distance dimension scale factor
real(DP) :: mscale = 1.0_DP !! Mass scale factor
real(DP) :: tscale = 1.0_DP !! Time scale factor
real(DP) :: vscale = 1.0_DP !! Velocity scale factor (a convenience unit that is derived from dscale and tscale)
real(DP) :: Escale = 1.0_DP !! Energy scale factor (a convenience unit that is derived from dscale, tscale, and mscale)
real(DP) :: Lscale = 1.0_DP !! Angular momentum scale factor (a convenience unit that is derived from dscale, tscale, and mscale)
contains
procedure :: generate_fragments => fraggle_generate_fragments !! Generates a system of fragments in barycentric coordinates that conserves energy and momentum.
procedure :: accel => fraggle_placeholder_accel !! Placeholder subroutine to fulfill requirement for an accel method
Expand Down Expand Up @@ -166,12 +166,6 @@ module subroutine fraggle_io_write_frame(self, nc, param)
class(swiftest_parameters), intent(inout) :: param !! Current run configuration parameters
end subroutine fraggle_io_write_frame

module subroutine fraggle_io_log_pl(pl, param)
implicit none
class(swiftest_pl), intent(in) :: pl !! Swiftest massive body object (only the new bodies generated in a collision)
class(swiftest_parameters), intent(in) :: param !! Current swiftest run configuration parameters
end subroutine fraggle_io_log_pl

module subroutine fraggle_io_log_regime(colliders, frag)
implicit none
class(fraggle_colliders), intent(in) :: colliders
Expand Down
18 changes: 0 additions & 18 deletions src/modules/symba_classes.f90
Original file line number Diff line number Diff line change
Expand Up @@ -176,8 +176,6 @@ module symba_classes
type, extends(symba_encounter) :: symba_plplenc
contains
procedure :: extract_collisions => symba_collision_extract_collisions_from_encounters !! Processes the pl-pl encounter list remove only those encounters that led to a collision
procedure :: resolve_fragmentations => symba_resolve_collision_fragmentations !! Process list of collisions, determine the collisional regime, and then create fragments
procedure :: resolve_mergers => symba_resolve_collision_mergers !! Process list of collisions and merge colliding bodies together
procedure :: resolve_collision => symba_resolve_collision_plplenc !! Process the pl-pl collision list, then modifiy the massive bodies based on the outcome of the c
end type symba_plplenc

Expand Down Expand Up @@ -232,22 +230,6 @@ module subroutine symba_collision_make_colliders_pl(self,idx)
integer(I4B), dimension(2), intent(in) :: idx !! Array holding the indices of the two bodies involved in the collision
end subroutine symba_collision_make_colliders_pl

module subroutine symba_resolve_collision_fragmentations(self, system, param, t)
implicit none
class(symba_plplenc), intent(inout) :: self !! SyMBA pl-pl encounter list
class(symba_nbody_system), intent(inout) :: system !! SyMBA nbody system object
class(symba_parameters), intent(inout) :: param !! Current run configuration parameters with SyMBA additions
real(DP), intent(in) :: t !! Time of collision
end subroutine symba_resolve_collision_fragmentations

module subroutine symba_resolve_collision_mergers(self, system, param, t)
implicit none
class(symba_plplenc), intent(inout) :: self !! SyMBA pl-pl encounter list
class(symba_nbody_system), intent(inout) :: system !! SyMBA nbody system object
class(symba_parameters), intent(inout) :: param !! Current run configuration parameters with SyMBA additions
real(DP), intent(in) :: t !! Time of collision
end subroutine symba_resolve_collision_mergers

module subroutine symba_resolve_collision_plplenc(self, system, param, t, dt, irec)
implicit none
class(symba_plplenc), intent(inout) :: self !! SyMBA pl-pl encounter list
Expand Down
Loading

0 comments on commit 80c9b92

Please sign in to comment.