!> A module which contains methods for checking whether nonlinear !> runs have reached a saturated steady state module diagnostics_nonlinear_convergence implicit none private public :: check_nonlin_convergence, check_nonlin_convergence_calc public :: init_nonlinear_convergence, finish_nonlinear_convergence integer :: conv_istep = 0 integer :: conv_isteps_converged = 0 real, allocatable, dimension(:) :: conv_heat real :: heat_sum_av = 0, heat_av = 0, heat_av_test = 0 contains !> FIXME : Add documentation subroutine init_nonlinear_convergence(conv_nstep_av, nwrite) implicit none integer, intent(in) :: conv_nstep_av, nwrite if(.not. allocated(conv_heat)) allocate(conv_heat(0:conv_nstep_av/nwrite-1)) heat_sum_av = 0 ; heat_av = 0 ; heat_av_test = 0 end subroutine init_nonlinear_convergence !> FIXME : Add documentation subroutine finish_nonlinear_convergence() if (allocated(conv_heat)) deallocate(conv_heat) end subroutine finish_nonlinear_convergence subroutine check_nonlin_convergence_calc(istep, nwrite, heat_flux, exit, & exit_when_converged, conv_nstep_av, conv_nsteps_converged, conv_test_multiplier, & conv_min_step, conv_max_step) use mp, only: proc0, broadcast, job use gs2_time, only: user_time logical, intent(in out) :: exit logical, intent(in) :: exit_when_converged integer, intent(in) :: istep, nwrite, conv_nstep_av, conv_nsteps_converged, conv_min_step integer, intent(in) :: conv_max_step real, intent(in) :: heat_flux, conv_test_multiplier integer :: place, iwrite, nwrite_av real :: heat_av_new, heat_av_diff logical, parameter :: debug = .false. if (proc0 .and. .not. (conv_istep /= 0 .and. istep == 0)) then ! Total number of steps including restarted trinity runs if (istep > 0) conv_istep = conv_istep + nwrite iwrite = conv_istep / nwrite ! Number of diagnostic write steps written nwrite_av = conv_nstep_av / nwrite ! Number of diagnostic write steps to average over ! Cumulative average of heat flux heat_sum_av = (heat_sum_av * iwrite + heat_flux) / (iwrite + 1) place = mod(conv_istep, conv_nstep_av) / nwrite !Just mod(iwrite, conv_nstep_av)? conv_heat(place) = heat_flux if (debug) write(6,'(A,I5,A,e11.4,A,I6,I6)') 'Job ',job, & ' time = ',user_time, ' step = ',conv_istep if (debug) write(6,'(A,I5,A,e11.4,A,e11.4)') 'Job ',job, & ' heat = ',heat_flux, ' heatsumav = ',heat_av if (conv_istep >= conv_nstep_av) then heat_av_new = sum(conv_heat) / nwrite_av heat_av_diff = heat_av_new - heat_av if(debug) write(6,'(A,I5,A,e11.4,A,e11.4)') 'Job ',job, & ' heat_sum_av_diff = ',heat_sum_av heat_av = heat_av_new ! Convergence test - needs to be met conv_nsteps_converged/nwrite times in succession if (abs(heat_av_diff) < heat_av_test) then conv_isteps_converged = conv_isteps_converged + 1 else conv_isteps_converged = 0 heat_av_test = heat_sum_av * conv_test_multiplier end if if ((conv_isteps_converged >= conv_nsteps_converged / nwrite) .and. & (conv_istep >= conv_min_step)) then if (debug) write(6,'(A,I5,A,I6,I3)')'Job ',job, & ' Reached convergence condition after step ',conv_istep exit = .true. .and. exit_when_converged endif if (conv_istep > conv_max_step) then write(6,'(A,I5,A,I7)') '*** Warning. Job ',job, & ' did not meet the convergence condition after ',conv_istep exit = .true. .and. exit_when_converged endif endif endif call broadcast(exit) end subroutine check_nonlin_convergence_calc !> Nonlinear convergence condition - simple and experimental look !> for the averaged differential of the summed averaged heat flux to !> drop below a threshold subroutine check_nonlin_convergence(gnostics) use diagnostics_config, only: diagnostics_type type (diagnostics_type), intent(inout) :: gnostics call check_nonlin_convergence_calc(gnostics%istep, gnostics%nwrite, & gnostics%current_results%total_heat_flux, gnostics%exit, & gnostics%exit_when_converged, & gnostics%conv_nstep_av, gnostics%conv_nsteps_converged, & gnostics%conv_test_multiplier, gnostics%conv_min_step, gnostics%conv_max_step) end subroutine check_nonlin_convergence end module diagnostics_nonlinear_convergence