diff --git a/src/swiftest/swiftest_drift.f90 b/src/swiftest/swiftest_drift.f90 index fe315fc74..0fdb0f63f 100644 --- a/src/swiftest/swiftest_drift.f90 +++ b/src/swiftest/swiftest_drift.f90 @@ -207,6 +207,12 @@ pure subroutine swiftest_drift_dan(mu, rx0, ry0, rz0, vx0, vy0, vz0, dt0, iflag) end if end if + ! For debugging overflow error + + if(r0 .ge. 11970204779.00) then + f = 1.0_DP + end if + call swiftest_drift_kepu(dt, r0, mu, alpha, u, fp, c1, c2, c3, iflag) if (iflag == 0) then f = 1.0_DP - mu / r0 * c2 @@ -462,10 +468,10 @@ pure subroutine swiftest_drift_kepu_new(s, dt, r0, mu, alpha, u, fp, c1, c2, c3, x = s*s*alpha ! for debugging overflow error - if (abs(x) .ge. HUGE(0.0_DP)) then - write(*,*) "big x" + if (abs(x) .ge. HUGE(0.0_DP) .or. x < 0) then + f = 1.0_DP ! "big x" end if - + call swiftest_drift_kepu_stumpff(x, c0, c1, c2, c3) c1 = c1*s c2 = c2*s*s @@ -559,28 +565,12 @@ pure subroutine swiftest_drift_kepu_stumpff(x, c0, c1, c2, c3) integer(I4B) :: i, n real(DP) :: xm - - ! for debugging Floating overflow error - if (abs(x) .ge. HUGE(0.0_DP)) then - write(*,*) "big x" - end if - - if (abs(c0) .ge. HUGE(0.0_DP)) then - write(*,*) "big c0" - end if - - if (abs(c1) .ge. HUGE(0.0_DP)) then - write(*,*) "big c1" - end if - - if (abs(c2) .ge. HUGE(0.0_DP)) then - write(*,*) "big c2" - end if - - if (abs(c3) .ge. HUGE(0.0_DP)) then - write(*,*) "big c3" + ! for debugging overflow error + if (abs(x) .ge. HUGE(0.0_DP) .or. x < 0) then + x = x ! "big x" + n = 0 end if - + n = 0 xm = 0.1_DP do while (abs(x) >= xm) @@ -605,6 +595,25 @@ pure subroutine swiftest_drift_kepu_stumpff(x, c0, c1, c2, c3) end do end if + ! for debugging Floating overflow error + + if (abs(c0) .ge. HUGE(0.0_DP)) then + xm = 0.1_DP ! "big c0" + end if + + if (abs(c1) .ge. HUGE(0.0_DP)) then + xm = 0.1_DP ! "big c1" + end if + + if (abs(c2) .ge. HUGE(0.0_DP)) then + xm = 0.1_DP ! "big c2" + end if + + if (abs(c3) .ge. HUGE(0.0_DP)) then + xm = 0.1_DP ! "big c3" + end if + + return end subroutine swiftest_drift_kepu_stumpff