diff --git a/src/collision/collision_generate.f90 b/src/collision/collision_generate.f90 index 333e5d059..27d02359c 100644 --- a/src/collision/collision_generate.f90 +++ b/src/collision/collision_generate.f90 @@ -226,14 +226,6 @@ module subroutine collision_generate_simple(self, nbody_system, param, t) write(*,*) "Error in swiftest_collision, unrecognized collision regime" call util_exit(FAILURE) end select - call collision_io_collider_message(nbody_system%pl, impactors%id, message) - call swiftest_io_log_one_message(COLLISION_LOG_OUT, message) - - ! Collisional fragments will be uniformly distributed around the pre-impact barycenter - call self%set_mass_dist(param) - - ! Generate the position and velocity distributions of the fragments - call self%disrupt(nbody_system, param, t) dpe = self%pe(2) - self%pe(1) nbody_system%Ecollisions = nbody_system%Ecollisions - dpe @@ -379,10 +371,12 @@ module subroutine collision_generate_disrupt(self, nbody_system, param, t, lfail ! Internals real(DP) :: r_max_start + call self%get_energy_and_momentum(nbody_system, param, lbefore=.true.) r_max_start = 1.1_DP * .mag.(self%impactors%rb(:,2) - self%impactors%rb(:,1)) call collision_generate_simple_pos_vec(self, r_max_start) call self%set_coordinate_system() call collision_generate_simple_vel_vec(self) + call self%get_energy_and_momentum(nbody_system, param, lbefore=.false.) return end subroutine collision_generate_disrupt @@ -508,6 +502,16 @@ module subroutine collision_generate_simple_vel_vec(collider) vimp_unit(:) = .unit. (fragments%rc(:,i) + impactors%rbcom(:) - impactors%rbimp(:)) fragments%vc(:,i) = vmag * vimp_unit(:) + vnoise(:,i) end do + do concurrent(i = 1:nfrag) + fragments%vb(:,i) = fragments%vc(:,i) + impactors%vbcom(:) + end do + + impactors%vbcom(:) = 0.0_DP + do concurrent(i = 1:nfrag) + impactors%vbcom(:) = impactors%vbcom(:) + fragments%mass(i) * fragments%vb(:,i) + end do + impactors%vbcom(:) = impactors%vbcom(:) / fragments%mtot + end associate return end subroutine collision_generate_simple_vel_vec diff --git a/src/fraggle/fraggle_generate.f90 b/src/fraggle/fraggle_generate.f90 index 3c48b71bf..0d8331970 100644 --- a/src/fraggle/fraggle_generate.f90 +++ b/src/fraggle/fraggle_generate.f90 @@ -35,7 +35,7 @@ module subroutine fraggle_generate_disrupt(self, nbody_system, param, t, lfailur integer(I4B) :: try real(DP) :: r_max_start, f_spin, dEtot, dLmag integer(I4B), parameter :: MAXTRY = 100 - logical :: lk_plpl + logical :: lk_plpl, lfailure_local logical, dimension(size(IEEE_ALL)) :: fpe_halting_modes, fpe_quiet_modes logical, dimension(size(IEEE_USUAL)) :: fpe_flag character(len=STRMAX) :: message @@ -59,7 +59,7 @@ module subroutine fraggle_generate_disrupt(self, nbody_system, param, t, lfailur if (nfrag < NFRAG_MIN) then write(message,*) "Fraggle needs at least ",NFRAG_MIN," fragments, but only ",nfrag," were given." call swiftest_io_log_one_message(COLLISION_LOG_OUT, message) - lfailure = .true. + lfailure_local = .true. return end if @@ -79,19 +79,19 @@ module subroutine fraggle_generate_disrupt(self, nbody_system, param, t, lfailur ! Start out the fragments close to the initial separation distance. This will be increased if there is any overlap or we fail to find a solution r_max_start = 1.1_DP * .mag.(impactors%rb(:,2) - impactors%rb(:,1)) - lfailure = .false. + lfailure_local = .false. try = 1 f_spin = F_SPIN_FIRST do while (try < MAXTRY) write(message,*) try call swiftest_io_log_one_message(COLLISION_LOG_OUT, "Fraggle try " // trim(adjustl(message))) - if (lfailure) then + if (lfailure_local) then call fragments%restructure(impactors, try, f_spin, r_max_start) call fragments%reset() try = try + 1 end if - lfailure = .false. + lfailure_local = .false. call ieee_set_flag(ieee_all, .false.) ! Set all fpe flags to quiet call collision_generate_simple_pos_vec(self, r_max_start) @@ -107,20 +107,20 @@ module subroutine fraggle_generate_disrupt(self, nbody_system, param, t, lfailur call collision_generate_simple_vel_vec(self) - call fraggle_generate_tan_vel(self, lfailure) - if (lfailure) then + call fraggle_generate_tan_vel(self, lfailure_local) + if (lfailure_local) then call swiftest_io_log_one_message(COLLISION_LOG_OUT, "Fraggle failed to find tangential velocities") cycle end if - call fraggle_generate_rad_vel(self, lfailure) - if (lfailure) then + call fraggle_generate_rad_vel(self, lfailure_local) + if (lfailure_local) then call swiftest_io_log_one_message(COLLISION_LOG_OUT, "Fraggle failed to find radial velocities") cycle end if - call fraggle_generate_spins(self, f_spin, lfailure) - if (lfailure) then + call fraggle_generate_spins(self, f_spin, lfailure_local) + if (lfailure_local) then call swiftest_io_log_one_message(COLLISION_LOG_OUT, "Fraggle failed to find spins") cycle end if @@ -130,16 +130,16 @@ module subroutine fraggle_generate_disrupt(self, nbody_system, param, t, lfailur dLmag = .mag. (self%Ltot(:,2) - self%Ltot(:,1)) exit - lfailure = ((abs(dEtot + impactors%Qloss) > FRAGGLE_ETOL) .or. (dEtot > 0.0_DP)) - if (lfailure) then + lfailure_local = ((abs(dEtot + impactors%Qloss) > FRAGGLE_ETOL) .or. (dEtot > 0.0_DP)) + if (lfailure_local) then write(message, *) dEtot, abs(dEtot + impactors%Qloss) / FRAGGLE_ETOL call swiftest_io_log_one_message(COLLISION_LOG_OUT, "Fraggle failed due to high energy error: " // & trim(adjustl(message))) cycle end if - lfailure = ((abs(dLmag) / (.mag.self%Ltot(:,1))) > FRAGGLE_LTOL) - if (lfailure) then + lfailure_local = ((abs(dLmag) / (.mag.self%Ltot(:,1))) > FRAGGLE_LTOL) + if (lfailure_local) then write(message,*) dLmag / (.mag.self%Ltot(:,1)) call swiftest_io_log_one_message(COLLISION_LOG_OUT, "Fraggle failed due to high angular momentum error: " // & trim(adjustl(message))) @@ -148,14 +148,14 @@ module subroutine fraggle_generate_disrupt(self, nbody_system, param, t, lfailur ! Check if any of the usual floating point exceptions happened, and fail the try if so call ieee_get_flag(ieee_usual, fpe_flag) - lfailure = any(fpe_flag) - if (.not.lfailure) exit + lfailure_local = any(fpe_flag) + if (.not.lfailure_local) exit write(message,*) "Fraggle failed due to a floating point exception: ", fpe_flag call swiftest_io_log_one_message(COLLISION_LOG_OUT, message) end do write(message,*) try - if (lfailure) then + if (lfailure_local) then call swiftest_io_log_one_message(COLLISION_LOG_OUT, "Fraggle fragment generation failed after " // & trim(adjustl(message)) // " tries") else @@ -167,6 +167,7 @@ module subroutine fraggle_generate_disrupt(self, nbody_system, param, t, lfailur ! Restore the big array if (lk_plpl) call pl%flatten(param) + if (present(lfailure)) lfailure = lfailure_local end associate end select end select