diff --git a/cmake/Modules/SetFortranFlags.cmake b/cmake/Modules/SetFortranFlags.cmake index d869e89b6..103ff7e45 100644 --- a/cmake/Modules/SetFortranFlags.cmake +++ b/cmake/Modules/SetFortranFlags.cmake @@ -319,3 +319,7 @@ SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_PROFILE "${CMAKE_Fortran_FLAGS_RELEASE}" "/Qopt-report:5 /traceback -g3" # Windows Intel "-pg -fbacktrace" ) + +SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_PROFILE "${CMAKE_Fortran_FLAGS_PROFILE}" + Fortran "-check bounds,pointers,uninit" # Intel + ) diff --git a/src/fraggle/fraggle_generate.f90 b/src/fraggle/fraggle_generate.f90 index 2d9e1b237..51c4817e9 100644 --- a/src/fraggle/fraggle_generate.f90 +++ b/src/fraggle/fraggle_generate.f90 @@ -420,13 +420,13 @@ module subroutine fraggle_generate_rot_vec(collider, nbody_system, param) fragments%rot(:,1) = L_spin(:) / (fragments%mass(1) * fragments%radius(1)**2 * fragments%Ip(3,1)) fragments%rotmag(1) = .mag.fragments%rot(:,1) ! Add in some random spin noise. The magnitude will be scaled by the before-after amount and the direction will be random - call random_number(fragments%rotmag(2:nfrag)) do i = 2,nfrag call random_number(rotdir) + call random_number(rotmag) rotdir = rotdir - 0.5_DP rotdir = .unit. rotdir - rotmag = random_scale_factor * .mag. (Lbefore(:) - Lafter(:)) / ((nfrag - 1) * fragments%mass(i) * fragments%radius(i)**2 * fragments%Ip(3,i)) - fragments%rot(:,i) = rotmag * fragments%rotmag(i) * rotdir + fragments%rotmag(i) = random_scale_factor * rotmag * .mag.L_spin(:) / ((nfrag - 1) * fragments%mass(i) * fragments%radius(i)**2 * fragments%Ip(3,i)) + fragments%rot(:,i) = fragments%rotmag(i) * rotdir end do end associate @@ -453,8 +453,10 @@ module subroutine fraggle_generate_vel_vec(collider, nbody_system, param, lfailu real(DP) :: vmag, vesc, dE, E_residual, ke_avail, ke_remove, dE_best, E_residual_best, fscale, f_spin, f_orbit, dE_metric integer(I4B), dimension(collider%fragments%nbody) :: vsign real(DP), dimension(collider%fragments%nbody) :: vscale, ke_rot_remove - real(DP), parameter :: vmin_initial_factor = 2.0_DP ! For the initial "guess" of fragment velocities, this is the maximum velocity relative to escape velocity that the fragments will have - real(DP) :: vmax_initial_factor = 8.0_DP + ! For the initial "guess" of fragment velocities, this is the minimum and maximum velocity relative to escape velocity that the fragments will have + real(DP) :: vmin_guess = 1.5_DP + real(DP) :: vmax_guess = 10.0_DP + real(DP) :: delta_v integer(I4B), parameter :: MAXLOOP = 20 integer(I4B), parameter :: MAXTRY = 100 real(DP), parameter :: SUCCESS_METRIC = 1.0e-3_DP @@ -491,23 +493,24 @@ module subroutine fraggle_generate_vel_vec(collider, nbody_system, param, lfailu vesc = sqrt(2 * sum(impactors%Gmass(:)) / sum(impactors%radius(:))) end if - ! Scale the magnitude of the velocity by the distance from the impact point - ! This will reduce the chances of fragments colliding with each other immediately, and is more physically correct - do concurrent(i = 1:nfrag) - rimp(:) = fragments%rc(:,i) - impactors%rbimp(:) - vscale(i) = .mag. rimp(:) / (.mag. (impactors%rb(:,2) - impactors%rb(:,1))) - end do - - ! Set the velocity scale factor to span from vmin/vesc to vmax/vesc - vscale(:) = vscale(:)/minval(vscale(:)) - fscale = log(vmax_initial_factor - vmin_initial_factor + 1.0_DP) / log(maxval(vscale(:))) - vscale(:) = vscale(:)**fscale + vmin_initial_factor - 1.0_DP - E_residual_best = huge(1.0_DP) lfailure = .false. dE_metric = huge(1.0_DP) outer: do try = 1, MAXTRY + ! Scale the magnitude of the velocity by the distance from the impact point + ! This will reduce the chances of fragments colliding with each other immediately, and is more physically correct + do concurrent(i = 1:nfrag) + rimp(:) = fragments%rc(:,i) - impactors%rbimp(:) + vscale(i) = .mag. rimp(:) / (.mag. (impactors%rb(:,2) - impactors%rb(:,1))) + end do + + ! Set the velocity scale factor to span from vmin/vesc to vmax/vesc + vscale(:) = vscale(:)/minval(vscale(:)) + fscale = log(vmax_guess - vmin_guess + 1.0_DP) / log(maxval(vscale(:))) + vscale(:) = vscale(:)**fscale + vmin_guess - 1.0_DP + + ! Set the velocities of all fragments using all of the scale factors determined above do concurrent(i = 1:nfrag) j = fragments%origin_body(i) @@ -588,6 +591,11 @@ module subroutine fraggle_generate_vel_vec(collider, nbody_system, param, lfailu call fraggle_generate_pos_vec(collider_local) call fraggle_generate_rot_vec(collider_local, nbody_system, param) collider_local%fail_scale = collider_local%fail_scale*1.01_DP + + ! Bring the minimum and maximum velocities closer together for the next round + delta_v = 0.125_DP * (vmax_guess - vmin_guess) + vmin_guess = vmin_guess + delta_v + vmax_guess = vmax_guess - delta_v end do outer lfailure = dE_best > 0.0_DP