MinimizeMHD3D_descent Subroutine

public subroutine MinimizeMHD3D_descent(sf, abstol, maxIter_in)

Uses

  • proc~~minimizemhd3d_descent~~UsesGraph proc~minimizemhd3d_descent t_minimizer_mhd3d%MinimizeMHD3D_descent MODgvec_MHD3D_visu MODgvec_MHD3D_visu proc~minimizemhd3d_descent->MODgvec_MHD3D_visu module~modgvec_analyze MODgvec_Analyze proc~minimizemhd3d_descent->module~modgvec_analyze module~modgvec_mhd3d_evalfunc MODgvec_MHD3D_evalFunc proc~minimizemhd3d_descent->module~modgvec_mhd3d_evalfunc module~modgvec_restart MODgvec_Restart proc~minimizemhd3d_descent->module~modgvec_restart module~modgvec_sol_var_mhd3d MODgvec_sol_var_MHD3D proc~minimizemhd3d_descent->module~modgvec_sol_var_mhd3d module~modgvec_globals MODgvec_Globals module~modgvec_analyze->module~modgvec_globals module~modgvec_mhd3d_evalfunc->module~modgvec_globals module~sll_m_spline_matrix sll_m_spline_matrix module~modgvec_mhd3d_evalfunc->module~sll_m_spline_matrix module~sll_m_spline_matrix_banded sll_m_spline_matrix_banded module~modgvec_mhd3d_evalfunc->module~sll_m_spline_matrix_banded module~modgvec_restart->module~modgvec_globals module~modgvec_c_sol_var MODgvec_c_sol_var module~modgvec_sol_var_mhd3d->module~modgvec_c_sol_var module~modgvec_sol_var_mhd3d->module~modgvec_globals module~modgvec_c_sol_var->module~modgvec_globals iso_fortran_env iso_fortran_env module~modgvec_globals->iso_fortran_env module~sll_m_spline_matrix->module~sll_m_spline_matrix_banded module~sll_m_errors sll_m_errors module~sll_m_spline_matrix->module~sll_m_errors module~sll_m_spline_matrix_base sll_m_spline_matrix_base module~sll_m_spline_matrix->module~sll_m_spline_matrix_base module~sll_m_spline_matrix_dense sll_m_spline_matrix_dense module~sll_m_spline_matrix->module~sll_m_spline_matrix_dense module~sll_m_working_precision sll_m_working_precision module~sll_m_spline_matrix->module~sll_m_working_precision module~sll_m_spline_matrix_banded->iso_fortran_env module~sll_m_assert sll_m_assert module~sll_m_spline_matrix_banded->module~sll_m_assert module~sll_m_spline_matrix_banded->module~sll_m_errors module~sll_m_spline_matrix_banded->module~sll_m_spline_matrix_base module~sll_m_spline_matrix_banded->module~sll_m_working_precision module~sll_m_errors->iso_fortran_env module~sll_m_spline_matrix_base->module~sll_m_working_precision module~sll_m_spline_matrix_dense->iso_fortran_env module~sll_m_spline_matrix_dense->module~sll_m_assert module~sll_m_spline_matrix_dense->module~sll_m_errors module~sll_m_spline_matrix_dense->module~sll_m_spline_matrix_base module~sll_m_spline_matrix_dense->module~sll_m_working_precision

Core minimization routine

Iterates until the (force-)tolerance or the maximum number of iterations is reached. Also writes output files after set number of iterations

Type Bound

t_minimizer_mhd3d

Arguments

Type IntentOptional Attributes Name
class(t_minimizer_mhd3d), intent(inout) :: sf

minimizer
real(kind=wp), intent(in) :: abstol

tolerance on the forces. If reached terminaters the minimization
integer, intent(in) :: maxIter_in

maximum number of iterations after which the iterations are terminated

Calls

proc~~minimizemhd3d_descent~~CallsGraph proc~minimizemhd3d_descent t_minimizer_mhd3d%MinimizeMHD3D_descent analyze analyze proc~minimizemhd3d_descent->analyze interface~writestate WriteState proc~minimizemhd3d_descent->interface~writestate none~set_to t_sol_var_MHD3D%set_to proc~minimizemhd3d_descent->none~set_to norm_2 norm_2 proc~minimizemhd3d_descent->norm_2 proc~evalenergy EvalEnergy proc~minimizemhd3d_descent->proc~evalenergy proc~evalforce EvalForce proc~minimizemhd3d_descent->proc~evalforce proc~logging_mhd3d t_minimizer_mhd3d%Logging_MHD3D proc~minimizemhd3d_descent->proc~logging_mhd3d proc~minimizemhd3d_resetdescent t_minimizer_mhd3d%MinimizeMHD3d_ResetDescent proc~minimizemhd3d_descent->proc~minimizemhd3d_resetdescent proc~sol_var_mhd3d_axby t_sol_var_MHD3D%sol_var_MHD3D_AXBY proc~minimizemhd3d_descent->proc~sol_var_mhd3d_axby proc~startlogging_mhd3d t_minimizer_mhd3d%StartLogging_MHD3D proc~minimizemhd3d_descent->proc~startlogging_mhd3d tau tau proc~minimizemhd3d_descent->tau velocity velocity proc~minimizemhd3d_descent->velocity writesfloutfile writesfloutfile proc~minimizemhd3d_descent->writesfloutfile proc~writestatetoascii WriteStateToASCII interface~writestate->proc~writestatetoascii proc~sol_var_mhd3d_set_to_scalar t_sol_var_MHD3D%sol_var_MHD3D_set_to_scalar none~set_to->proc~sol_var_mhd3d_set_to_scalar proc~sol_var_mhd3d_set_to_solvar t_sol_var_MHD3D%sol_var_MHD3D_set_to_solvar none~set_to->proc~sol_var_mhd3d_set_to_solvar interface~par_allreduce par_AllReduce proc~evalenergy->interface~par_allreduce proc~evalaux EvalAux proc~evalenergy->proc~evalaux interface~par_bcast par_Bcast proc~evalforce->interface~par_bcast interface~par_ibcast par_IBcast proc~evalforce->interface~par_ibcast interface~par_ireduce par_IReduce proc~evalforce->interface~par_ireduce par_wait par_wait proc~evalforce->par_wait proc~applybc_fstrong ApplyBC_Fstrong proc~evalforce->proc~applybc_fstrong proc~applyprecond ApplyPrecond proc~evalforce->proc~applyprecond proc~buildprecond BuildPrecond proc~evalforce->proc~buildprecond proc~evalforce->proc~evalaux proc~fbase_projectiptodof_tens t_fBase%fBase_projectIPtoDOF_tens proc~evalforce->proc~fbase_projectiptodof_tens proc~sbase_applybctorhs t_sBase%sBase_applyBCtoRHS proc~evalforce->proc~sbase_applybctorhs solve_inplace solve_inplace proc~evalforce->solve_inplace proc~logging_mhd3d->none~set_to checkdistance checkdistance proc~logging_mhd3d->checkdistance checkpos checkpos proc~logging_mhd3d->checkpos vnorm vnorm proc~logging_mhd3d->vnorm proc~minimizemhd3d_resetdescent->none~set_to proc~minimizemhd3d_resetdescent->proc~evalenergy proc~minimizemhd3d_resetdescent->proc~evalforce proc~minimizemhd3d_resetdescent->tau proc~minimizemhd3d_resetdescent->proc~evalaux proc~sol_var_mhd3d_norm_2 t_sol_var_MHD3D%sol_var_MHD3D_norm_2 proc~minimizemhd3d_resetdescent->proc~sol_var_mhd3d_norm_2 proc~startlogging_mhd3d->checkpos interface~getfreeunit GETFREEUNIT proc~startlogging_mhd3d->interface~getfreeunit proc~startlogging_mhd3d->vnorm interface~getfreeunit->interface~getfreeunit proc~par_allreduce_array1d par_AllReduce_array1D interface~par_allreduce->proc~par_allreduce_array1d proc~par_allreduce_array2d par_AllReduce_array2D interface~par_allreduce->proc~par_allreduce_array2d proc~par_allreduce_scalar par_AllReduce_scalar interface~par_allreduce->proc~par_allreduce_scalar proc~par_allreduce_scalar_int par_AllReduce_scalar_int interface~par_allreduce->proc~par_allreduce_scalar_int proc~par_bcast_array1d par_Bcast_array1D interface~par_bcast->proc~par_bcast_array1d proc~par_bcast_array1d_int par_Bcast_array1D_int interface~par_bcast->proc~par_bcast_array1d_int proc~par_bcast_array1d_str par_Bcast_array1D_str interface~par_bcast->proc~par_bcast_array1d_str proc~par_bcast_array2d par_Bcast_array2D interface~par_bcast->proc~par_bcast_array2d proc~par_bcast_scalar par_Bcast_scalar interface~par_bcast->proc~par_bcast_scalar proc~par_bcast_scalar_int par_Bcast_scalar_int interface~par_bcast->proc~par_bcast_scalar_int proc~par_bcast_scalar_str par_Bcast_scalar_str interface~par_bcast->proc~par_bcast_scalar_str proc~par_ibcast_array1d par_IBcast_array1D interface~par_ibcast->proc~par_ibcast_array1d proc~par_ibcast_array2d par_IBcast_array2D interface~par_ibcast->proc~par_ibcast_array2d proc~par_ireduce_array1d par_IReduce_array1D interface~par_ireduce->proc~par_ireduce_array1d proc~par_ireduce_array2d par_IReduce_array2D interface~par_ireduce->proc~par_ireduce_array2d proc~sbase_applybctodof_lgm t_sBase%sBase_applyBCtoDOF_LGM proc~applybc_fstrong->proc~sbase_applybctodof_lgm proc~s_spline_matrix_banded__solve_inplace sll_t_spline_matrix_banded%s_spline_matrix_banded__solve_inplace proc~applyprecond->proc~s_spline_matrix_banded__solve_inplace proc~buildprecond->interface~par_allreduce add_element add_element proc~buildprecond->add_element factorize factorize proc~buildprecond->factorize get_element get_element proc~buildprecond->get_element reset reset proc~buildprecond->reset set_element set_element proc~buildprecond->set_element proc~evalaux->interface~par_allreduce eval_all eval_all proc~evalaux->eval_all proc~base_evaldof t_base%base_evalDOF proc~evalaux->proc~base_evaldof proc~base_evaldof_all t_base%base_evalDOF_all proc~evalaux->proc~base_evaldof_all dgemm dgemm proc~fbase_projectiptodof_tens->dgemm proc~writestatetoascii->interface~getfreeunit proc~evaltotals EvalTotals proc~writestatetoascii->proc~evaltotals proc~rprofile_eval_at_rho c_rProfile%rProfile_eval_at_rho proc~writestatetoascii->proc~rprofile_eval_at_rho proc~fbase_evaldof_ip_tens t_fBase%fBase_evalDOF_IP_tens proc~base_evaldof->proc~fbase_evaldof_ip_tens proc~base_evaldof_all->proc~fbase_evaldof_ip_tens proc~evaltotals->proc~evalaux interface~par_reduce par_Reduce proc~evaltotals->interface~par_reduce eval_at_rho2 eval_at_rho2 proc~rprofile_eval_at_rho->eval_at_rho2 proc~rho2_derivative rho2_derivative proc~rprofile_eval_at_rho->proc~rho2_derivative proc~rprofile_drho2 c_rProfile%rProfile_drho2 proc~rprofile_eval_at_rho->proc~rprofile_drho2 proc~rprofile_drho3 c_rProfile%rProfile_drho3 proc~rprofile_eval_at_rho->proc~rprofile_drho3 proc~rprofile_drho4 c_rProfile%rProfile_drho4 proc~rprofile_eval_at_rho->proc~rprofile_drho4 dgbtrs dgbtrs proc~s_spline_matrix_banded__solve_inplace->dgbtrs proc~sll_s_error_handler sll_s_error_handler proc~s_spline_matrix_banded__solve_inplace->proc~sll_s_error_handler proc~solve SOLVE proc~sbase_applybctodof_lgm->proc~solve proc~par_reduce_array1d par_Reduce_array1D interface~par_reduce->proc~par_reduce_array1d proc~par_reduce_array2d par_Reduce_array2D interface~par_reduce->proc~par_reduce_array2d proc~par_reduce_scalar par_Reduce_scalar interface~par_reduce->proc~par_reduce_scalar proc~par_reduce_scalar_int par_Reduce_scalar_int interface~par_reduce->proc~par_reduce_scalar_int proc~fbase_evaldof_ip_tens->dgemm proc~fbase_evaldof_xn t_fBase%fBase_evalDOF_xn proc~fbase_evaldof_ip_tens->proc~fbase_evaldof_xn proc~poly_derivative_prefactor poly_derivative_prefactor proc~rho2_derivative->proc~poly_derivative_prefactor proc~rprofile_drho2->eval_at_rho2 proc~rprofile_drho2->proc~rho2_derivative proc~rprofile_drho3->eval_at_rho2 proc~rprofile_drho3->proc~rho2_derivative proc~rprofile_drho4->eval_at_rho2 proc~rprofile_drho4->proc~rho2_derivative interface~c_abort~2 c_abort proc~sll_s_error_handler->interface~c_abort~2 proc~errout errout proc~sll_s_error_handler->proc~errout dgetrf dgetrf proc~solve->dgetrf dgetrs dgetrs proc~solve->dgetrs dgemv dgemv proc~fbase_evaldof_xn->dgemv proc~fbase_eval_xn t_fBase%fBase_eval_xn proc~fbase_evaldof_xn->proc~fbase_eval_xn

Called by

proc~~minimizemhd3d_descent~~CalledByGraph proc~minimizemhd3d_descent t_minimizer_mhd3d%MinimizeMHD3D_descent proc~minimizemhd3d t_functional_mhd3d%MinimizeMHD3D proc~minimizemhd3d->proc~minimizemhd3d_descent proc~minimize minimize proc~minimize->proc~minimizemhd3d proc~rungvec rungvec proc~rungvec->proc~minimizemhd3d proc~start_rungvec start_rungvec proc~start_rungvec->proc~rungvec program~gvec GVEC program~gvec->proc~rungvec

Source Code

    SUBROUTINE MinimizeMHD3D_descent(sf, abstol, maxIter_in)
        USE MODgvec_MHD3D_EvalFunc, ONLY: EvalAux, EvalForce, EvalEnergy
        USE MODgvec_Analyze, ONLY:analyze
        USE MODgvec_Restart, ONLY:WriteState
        USE MODgvec_MHD3D_visu, ONLY:WriteSFLoutfile
        USE MODgvec_sol_var_MHD3D, ONLY: t_sol_var_MHD3D
        IMPLICIT NONE

        REAL(wp), INTENT(IN) :: abstol !! tolerance on the forces. If reached terminaters the minimization
        INTEGER, INTENT(IN)  :: maxIter_in !! maximum number of iterations after which the iterations are terminated
        !-----------------------------------------------------------------------------------------------------------------------------------
        ! OUTPUT VARIABLES
        CLASS(t_minimizer_mhd3d), INTENT(INOUT) :: sf !! minimizer
        !-----------------------------------------------------------------------------------------------------------------------------------
        ! LOCAL VARIABLES
        INTEGER :: maxIter !! maximum number of iterations
        !=================================================================================================================================
        ASSOCIATE(vars=>sf%vars)
            maxIter = vars%iter + maxIter_in
            DO WHILE(vars%iter.LT.maxIter)
                IF((vars%restart_iter))THEN
                    CALL sf%reset()
                END IF !before first iteration or after restart Jac<0

                IF(.NOT.vars%logger_is_initialized)THEN
                    CALL sf%StartLogging()
                END IF
                !COMPUTE NEW SOLUTION P(1) as a prediction
                SELECT TYPE(vars)
                TYPE IS(t_gradient_descent_vars) !gradient descent, previously used for minimizerType=0
                    CALL vars%temp_dofs(1)%AXBY(1.0_wp,vars%dofs(0),vars%dt,vars%force(0)) !overwrites P(1), predicts solution U(1)
                TYPE IS(t_accelerated_gradient_descent_vars) !hirshman method
                    !tau is damping parameter
                    vars%tau(1:vars%ndamp-1) = vars%tau(2:vars%ndamp) !save old

                    !ln(|F_n|^2/|F_{n-1}|^2), Fnorm=|F_X1|,|F_X2|,|F_LA|
                    vars%tau(vars%ndamp)  = MIN(0.15_wp,ABS(LOG(SUM(vars%Fnorm**2)/SUM(vars%Fnorm_old**2))))/vars%dt

                    vars%tau_bar = 0.5_wp*vars%dt*SUM(vars%tau)/REAL(vars%ndamp,wp)   !=1/2 * tauavg
                    CALL vars%velocity(1)%AXBY(((1.0_wp-vars%tau_bar)/(1.0_wp+vars%tau_bar)),vars%velocity(0),(vars%dt/(1.0_wp+vars%tau_bar)),vars%force(0)) !velocity V(1)
                    CALL vars%temp_dofs(1)%AXBY(1.0_wp,vars%dofs(0),vars%dt,vars%velocity(1)) !overwrites P(1), predicst solution U(1)
                    vars%Vnorm=SQRT(vars%velocity(1)%norm_2())
                END SELECT !Type

                vars%JacCheck=2 !no abort,if detJ<0, JacCheck<=-1

                vars%temp_dofs(1)%W_MHD3D=EvalEnergy(vars%temp_dofs(1),.TRUE.,vars%JacCheck)
                IF(vars%JacCheck.LE.-1)THEN
                    vars%nstepDecreased=vars%nStepDecreased+1
                    IF(vars%iter.GT.100) THEN
                        CALL WriteState(vars%dofs(0),vars%iter)
                        CALL abort(__STAMP__,&
                                "detJ<0. This is an abort during minimization (iteration count>100), where detJ<0 occurred.",&
                                TypeInfo="RuntimeError")
                    END IF
                    IF(vars%nStepDecreased.GT.130) THEN
                        CALL WriteState(vars%dofs(0),vars%iter)
                        CALL abort(__STAMP__,&
                                "detJ<0. This is an abort, since the timestep has been decreased by a factor of 1e-6 which still produces detJ<0.",&
                                RealInfo=vars%dt, TypeInfo="InitializationError")
                    END IF
                    vars%dt=0.9_wp*vars%dt
                    vars%nSkip_Jac=vars%nSkip_Jac+1
                    vars%restart_iter=.TRUE.
                    CALL vars%dofs(0)%set_to(vars%dofs(-3)) !reset to initial state
                    SWRITE(UNIT_stdOut,'(8X,I8,A,E11.4,A)')vars%iter,'...detJac<0, iter<=100, decrease stepsize to dt=',vars%dt,  ' and RESTART simulation!!!!!!!'
                ELSE
                    !detJ>0
                    vars%deltaW=vars%temp_dofs(1)%W_MHD3D-vars%dofs(0)%W_MHD3D!should be <=0,
                    IF(vars%deltaW.LE.vars%dW_allowed*vars%W_MHD3D_0)THEN !valid step /hirshman method accept W increase!
                        IF(ALL(vars%Fnorm.LE.abstol))THEN
                            CALL sf%Logging(.FALSE.)
                            SWRITE(UNIT_stdOut,'(4x,A)')'==>Iteration finished, |force| in relative tolerance'
                            EXIT !DO LOOP
                        END IF
                        vars%iter=vars%iter+1
                        vars%t_pseudo=vars%t_pseudo+vars%dt
                        ! for simple gradient & hirshman
                        CALL vars%dofs(-1)%set_to(vars%dofs(0))
                        CALL vars%dofs(0)%set_to(vars%temp_dofs(1))
                        ! for hirshman method
                        SELECT TYPE(vars)
                        TYPE IS(t_accelerated_gradient_descent_vars)
                            CALL vars%velocity(-1)%set_to(vars%velocity(0))
                            CALL vars%velocity(0)%set_to(vars%velocity(1))
                        END SELECT
                        CALL EvalForce(vars%temp_dofs(1),.FALSE.,vars%JacCheck,vars%force(0)) !evalAux was already called on P(1)=U(0), so that its set false here.
                        vars%Fnorm_old=vars%Fnorm
                        vars%Fnorm=SQRT(vars%force(0)%norm_2())
                        vars%nstepDecreased=0
                        vars%min_dt_out=MIN(vars%min_dt_out,vars%dt)
                        vars%max_dt_out=MAX(vars%max_dt_out,vars%dt)
                        vars%min_dW_out=MIN(vars%min_dW_out,vars%deltaW)
                        vars%max_dW_out=MAX(vars%max_dW_out,vars%deltaW)
                        vars%sum_dW_out=vars%sum_dW_out+vars%deltaW
                        IF(MOD(vars%iter,vars%logIter_ramp).EQ.0)THEN
                            CALL sf%Logging(.NOT.((vars%logIter_ramp.GE.vars%logIter).AND.(MOD(vars%logscreen,vars%nLogScreen).EQ.0)))
                            IF(.NOT.(vars%logIter_ramp.LT.vars%logIter))THEN !only reset for logIter
                                vars%logscreen=vars%logscreen+1
                                vars%min_dt_out=1.0e+30_wp
                                vars%max_dt_out=0.0_wp
                                vars%min_dW_out=1.0e+30_wp
                                vars%max_dW_out=-1.0e+30_wp
                                vars%sum_dW_out=0.0_wp
                                vars%nSkip_dW =0
                            END IF
                            vars%logIter_ramp=MIN(vars%logIter,vars%logIter_ramp*2)
                        END IF
                    ELSE !not a valid step, decrease timestep and skip P(1)
                        vars%dt=0.9_wp*vars%dt
                        vars%nstepDecreased=vars%nStepDecreased+1
                        vars%nSkip_dW=vars%nSkip_dW+1
                        vars%restart_iter=.TRUE.
                        SWRITE(UNIT_stdOut,'(8X,I8,A,E8.2,A,E8.1,A,E11.4)')vars%iter,'...deltaW=',vars%deltaW,'>',vars%dW_allowed,&
                        '*W_MHD3D_0, skip step and decrease stepsize to dt=',vars%dt
                        IF(vars%nStepDecreased.GT.130) THEN ! 0.9^130 ~10^-6
                            SWRITE(UNIT_stdOut,'(A,E21.11)')'Iteration stopped since no decent direction found with decreased timestep: ', vars%dt
                            SWRITE(UNIT_stdOut,fmt_sep)
                            RETURN
                        END IF
                    END IF !deltaW>dW_allowed*W_MHD_0 check
                END IF !JacCheck

                IF((MOD(vars%iter,vars%outputIter).EQ.0).AND.(vars%lastoutputIter.NE.vars%iter))THEN
                    __PERFON('output')
                    SWRITE(UNIT_stdOut,'(A)')'##########################  OUTPUT ##################################'
                    CALL Analyze(vars%iter, vars%dofs(0), vars%force(0))
                    CALL WriteState(vars%dofs(0),vars%iter)
                    SWRITE(UNIT_stdOut,'(A)')'#####################################################################'
                    vars%lastOutputIter=vars%iter
                    __PERFOFF('output')
                END IF
            END DO !iter

            IF(vars%iter.GE.MaxIter)THEN
                SWRITE(UNIT_stdOut,'(A,E21.11)')"maximum iteration count exceeded"
            END IF
            SWRITE(UNIT_stdOut,'(A)') "... DONE."
            SWRITE(UNIT_stdOut,fmt_sep)
            IF(vars%lastoutputIter.NE.vars%iter)THEN
                CALL Analyze(MIN(vars%iter,MaxIter), vars%dofs(0), vars%force(0))
                CALL WriteState(vars%dofs(0),MIN(vars%iter,MaxIter))
            END IF
            CALL writeSFLoutfile(vars%dofs(0),MIN(vars%iter,MaxIter))
        END ASSOCIATE !vars
    END SUBROUTINE MinimizeMHD3D_descent