From 002dfc36b7cabfc8e17d00c3505365b1b27d7aa1 Mon Sep 17 00:00:00 2001 From: David A Minton Date: Sun, 12 Feb 2023 12:22:14 -0500 Subject: [PATCH] Let the upper limit on the number of fragments scale with the fragment reduction factor --- src/fraggle/fraggle_util.f90 | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/src/fraggle/fraggle_util.f90 b/src/fraggle/fraggle_util.f90 index 1bf4fa961..98ebe0654 100644 --- a/src/fraggle/fraggle_util.f90 +++ b/src/fraggle/fraggle_util.f90 @@ -86,14 +86,14 @@ module subroutine fraggle_util_set_mass_dist(self, param) class(collision_fraggle), intent(inout) :: self !! Fraggle collision system object class(swiftest_parameters), intent(in) :: param !! Current Swiftest run configuration parameters ! Internals - integer(I4B) :: i, j, jproj, jtarg, nfrag, istart + integer(I4B) :: i, j, jproj, jtarg, nfrag, istart, nfragmax real(DP), dimension(2) :: volume real(DP), dimension(NDIM) :: Ip_avg real(DP) :: mfrag, mremaining, mtot, mcumul, G, mass_noise real(DP), dimension(:), allocatable :: mass real(DP), parameter :: BETA = 2.85_DP integer(I4B), parameter :: MASS_NOISE_FACTOR = 5 !! The number of digits of random noise that get added to the minimum mass value to prevent identical masses from being generated in a single run - integer(I4B), parameter :: NFRAGMAX = 1000 !! Maximum number of fragments that can be generated + integer(I4B), parameter :: NFRAGMAX_UNSCALED = 3000 !! Maximum number of fragments that can be generated integer(I4B), parameter :: NFRAGMIN = 1 !! Minimum number of fragments that can be generated (set by the fraggle_generate algorithm for constraining momentum and energy) integer(I4B), parameter :: iMlr = 1 integer(I4B), parameter :: iMslr = 2 @@ -107,6 +107,7 @@ module subroutine fraggle_util_set_mass_dist(self, param) mtot = sum(impactors%mass(:)) G = impactors%Gmass(1) / impactors%mass(1) Ip_avg(:) = (impactors%mass(1) * impactors%Ip(:,1) + impactors%mass(2) * impactors%Ip(:,2)) / mtot + nfragmax = ceiling(NFRAGMAX_UNSCALED / param%nfrag_reduction) if (impactors%mass(1) >= impactors%mass(2)) then jtarg = 1 @@ -135,7 +136,7 @@ module subroutine fraggle_util_set_mass_dist(self, param) mfrag = max((nfrag - 1)**(-3._DP / BETA) * impactors%mass_dist(iMslr), min_mfrag) mremaining = mremaining - mfrag end do - nfrag = min(ceiling(nfrag / param%nfrag_reduction), NFRAGMAX) + nfrag = min(ceiling(nfrag / param%nfrag_reduction), nfragmax) call self%setup_fragments(nfrag) case (COLLRESOLVE_REGIME_MERGE, COLLRESOLVE_REGIME_GRAZE_AND_MERGE)