From 4cfb29a45f4bb9c3938d2651a2dbad80bf8e9757 Mon Sep 17 00:00:00 2001 From: David Minton Date: Wed, 19 May 2021 23:17:12 -0400 Subject: [PATCH] Fixed some minor issues to help with convergence and consistency --- src/symba/symba_frag_pos.f90 | 5 ++--- src/util/util_minimize_bfgs.f90 | 8 ++++++-- 2 files changed, 8 insertions(+), 5 deletions(-) diff --git a/src/symba/symba_frag_pos.f90 b/src/symba/symba_frag_pos.f90 index d1ba8ac8a..a1d8d900b 100644 --- a/src/symba/symba_frag_pos.f90 +++ b/src/symba/symba_frag_pos.f90 @@ -500,7 +500,7 @@ subroutine set_fragment_position_vectors() ! Place the fragments into a region that is big enough that we should usually not have overlapping bodies ! An overlapping bodies will collide in the next time step, so it's not a major problem if they do (it just slows the run down) - r_max = 1.5_DP + r_max = 1.0_DP ! We will treat the first fragment of the list as a special case. It gets initialized the maximum distance along the original impactor distance vector. ! This is done because in a regular disruption, the first body is the largest fragment. @@ -522,9 +522,8 @@ subroutine set_fragment_position_vectors() loverlap(i) = loverlap(i) .or. (dis <= (rad_frag(i) + rad_frag(j))) end do end do + r_max = r_max + 0.2_DP end do - - x_frag(:, :) = x_frag(:, :) * r_max call shift_vector_to_origin(m_frag, x_frag) do i = 1, nfrag diff --git a/src/util/util_minimize_bfgs.f90 b/src/util/util_minimize_bfgs.f90 index d529dc11b..7aae0f89e 100644 --- a/src/util/util_minimize_bfgs.f90 +++ b/src/util/util_minimize_bfgs.f90 @@ -28,7 +28,7 @@ function util_minimize_bfgs(f, N, x0_d, eps_d, lerr) result(x1_d) real(DP), dimension(:), allocatable :: x1_d ! Internals integer(I4B) :: i, j, k, l, conv, num - integer(I4B), parameter :: MAXLOOP = 100 !! Maximum number of loops before method is determined to have failed + integer(I4B), parameter :: MAXLOOP = 500 !! Maximum number of loops before method is determined to have failed real(QP), parameter :: gradeps = 1e-5_QP !! Tolerance for gradient calculations real(QP), dimension(N) :: S !! Direction vectors real(QP), dimension(N) :: Snorm !! normalized direction @@ -110,13 +110,17 @@ function util_minimize_bfgs(f, N, x0_d, eps_d, lerr) result(x1_d) ! update H matrix H(:,:) = H(:,:) + ((1._QP - yHy / Py) * PP(:,:) - PyH(:,:) - HyP(:,:)) / Py if (any(H(:,:) > sqrt(huge(1._QP)) / N)) then + lerr = .true. write(*,*) 'BFGS did not converge after ',i,'iterations: H too big' exit end if ! Stop everything if there are any exceptions to allow the routine to fail gracefully call ieee_get_flag(ieee_usual, fpe_flag) if (any(fpe_flag)) exit - if (i == MAXLOOP) write(*,*) "BFGS ran out of loops!" + if (i == MAXLOOP) then + lerr = .true. + write(*,*) "BFGS ran out of loops!" + end if end do x1_d = x1 call ieee_get_flag(ieee_usual, fpe_flag)