Class | hmc_ratio_quark_wilson_class |
In: |
HmcRationalQuarkWilsonClass/hmc_ratio_quark_wilson_class.F90
|
Action = sum_{j=1,Nf} ye_{j}^+ R_n(Qee) ye_{j} - Nf TrLog[F]
where
R_n(x) = (1/x)^(1/2) v v Qee = Dee^+ Dee : Hermitian operator v Dee = 1ee - kapap^2 Fee \ Meo Foo \ Moe, even/odd site preconditioned operator v Dee^dag = 1ee - kapap^2 gamma5 Meo Foo \ Moe Fee \ gamma5 v = gamma5 Fee Dee Fee \ gamma5
Rational Approximation : R_n(x) = P_n(x)/Q_n(x), x in [a,b]
Partial Fraction form is used.
R_n(x) = (1/x)^(1/2) = alp_0 + \sum_{j=1..n} alp_j/(x-bet_j) alp_j : residue bet_j : pole
$Id: hmc_ratio_quark_wilson_class.F90,v 1.26 2011/06/21 03:29:27 ishikawa Exp $
Subroutine : | |
this : | type(hmc_ratio_quark_wilson), intent(inout) |
delete wilson/clover quark RHMC algorithm
subroutine delete_wilson_quark_rhmc(this) ! ! delete wilson/clover quark RHMC algorithm ! implicit none type(hmc_ratio_quark_wilson), intent(inout) :: this if (allocated(this%y)) deallocate(this%y) call delete(this%phys_param) call delete(this%solver_log) call delete(this%eigen_log) call delete(this%eigen_force) call delete(this%eigen_hamil) call delete(this%eigen_hamil_hist(0)) call delete(this%eigen_hamil_hist(1)) call delete(this%eigen_hamil_hist(2)) call delete(this%eigen_union) return end subroutine
Subroutine : | |||
this : | type(hmc_ratio_quark_wilson),
intent(inout)
| ||
BB : | type(vfield_gluon_wg), intent(inout)
| ||
idepth : | integer, intent(in)
| ||
u : | type(vfield_gluon_wg), intent(in)
|
Calc MD force for Rational HMC action.
Action : Srhmc Srhm = \sum_{j=1,Nf} ye(j)^+ R_n(Qee) ye(j) (even sites only) - Nf*TrLog[F]
where
R_n(x) = (1/x)^(Nf/2) = alp_0 + \sum_{k=1,Nratio} alp_k/(x-bet_k) v v Q = Dee^+ Dee
subroutine force_wilson_quark_rhmc(this,BB,idepth,u) ! ! Calc MD force for Rational HMC action. ! ! Action : Srhmc ! ! Srhm = \sum_{j=1,Nf} ye(j)^+ R_n(Qee) ye(j) (even sites only) ! ! - Nf*TrLog[F] ! ! where ! ! R_n(x) = (1/x)^(Nf/2) = alp_0 + \sum_{k=1,Nratio} alp_k/(x-bet_k) ! ! v v ! Q = Dee^+ Dee ! ! use pfcoef_class implicit none type(hmc_ratio_quark_wilson), intent(inout) :: this ! quark action/algorithm parameters integer, intent(in) :: idepth ! MD depth for recuresive MD integrator type(vfield_gluon_wg), intent(in) :: u ! gauge field type(vfield_gluon_wg), intent(inout) :: BB ! pre-force contribution (\dot{u}) real(DP) :: kappa2,fcoef, rtmp type(field_quark_wg), allocatable :: fy,fx type(field_quark_eo_wg), allocatable :: fye(:),we #ifdef _CLOVER real(DP) :: dummy real(DP) :: csw complex(DP) :: zickd8,ctmp #endif type(partial_fraction_coef) :: pfcoef character(len=CHARLEN) :: str integer :: iout,iiter integer :: ir,nr integer :: iflv,nf real(DP) :: kappa,err #ifdef _CLOVER if (idepth == this%depth_hopping .or. idepth == this%depth_clvtrlog) then !======================================================= ! comp. inverce clover term matrix ! fclinv = [1 - csw kappa/2 sigma_{mu,nu}F_{mu,nu}]^-1 !======================================================= call make_inverse_clover_term(this%phys_param,0,dummy) endif #endif if (idepth == this%depth_hopping) then #ifdef _MD_DEPTH_ if (nodeid==0) write(*,'("depth=",I3," : force quark poly hopping")')idepth #endif write(str,'("# FORCE SOLVER")') call print(this%solver_log,TRIM(str)) iout = get_file_unit(this%solver_log) nf = get_nflavor(this%phys_param) kappa = get_kappa(this%phys_param) nr = this%nratio_force if (allocated(this%alp)) deallocate(this%alp) if (allocated(this%bet)) deallocate(this%bet) call new(pfcoef,TYPE_INV_SQUARE_ROOT,nr,get_min(this%eigen_force),get_max(this%eigen_force)) allocate(this%alp(0:nr)) allocate(this%bet(1:nr)) this%alp(:) = pfcoef%alp(:) this%bet(:) = pfcoef%bet(:) call assign(this%eigen_func_eval,this%eigen_force) write(iout,'("# order=",I3," minmaxerr=",E24.15," [a,b]=",2E24.15)')nr,pfcoef%minimaxerror, get_min(this%eigen_force), get_max(this%eigen_force) call delete(pfcoef) allocate(fye(1:nr),fx,fy,we) do ir=1,nr call new(fye(ir),0) enddo call new(fx) call new(fy) call new(we,0) !================ ! flavor loop !================ do iflv = 1,nf !======================================== ! ! fye(j) = (Q-bet_j) \ ye(iflv) ! ! using blocked/multi-shift solver !======================================== err = this%force_inv%tol call assign_inv_mult_Qee_shift(this%phys_param,iout,err,iiter, fye,this%y(iflv)%eo(0),u,this%bet) call inc(this%force_inv,iiter) do ir=1,nr !======================================== ! v ! fyo = Moe fye(j) ! v ! fxe = Dee fye(j) ! ! fxe = Fee \ fxe ! v ! fxo = gamma5 Moe gamma5 fxe !======================================== call assign(fy%eo(0),fye(ir)) call assign_mult_hop( this%phys_param,fy%eo(1),fy%eo(0),u) call assign_mult_eoprec_wd(this%phys_param,fx%eo(0),fy%eo(0),u) #ifdef _CLOVER call accum_mult_clover_term(fx%eo(0),this%phys_param) #endif call assign_mult_gamma5(we,fx%eo(0)) call assign_mult_hop(this%phys_param,fx%eo(1),we,u) call accum_mult_gamma5(fx%eo(1)) !====================================== ! force from hopping matrix !====================================== call copy_boundary(fx) call copy_boundary(fy) fcoef = (kappa**2)*this%alp(ir) call force_hmc_hopping(BB,fcoef,fx,fy) #ifdef _CLOVER call assign_mult_hop(this%phys_param,fy%eo(0),fy%eo(1),u) !====================================== ! force from clover term in hopping !====================================== csw = get_csw(this%phys_param) zickd8 = ((zi*csw*kappa**3)/8.0_DP)*this%alp(ir) call force_hmc_clover(this%XX,zickd8,fx,fy) #endif enddo ! end of do ir enddo ! end of do iflv deallocate(fye,fx,fy,we) endif #ifdef _CLOVER if (idepth == this%depth_clvtrlog) then #ifdef _MD_DEPTH_ if (nodeid==0) write(*,'("depth=",I3," : force quark poly trlog")')idepth #endif !====================================== ! calc force from -Nf*log(det[F]) !====================================== call force_clover_trlog(this%phys_param,this%XX) endif call delete_inverse_clover_term(this%phys_param) !====================================== ! unlink clover term contribution ! XX is used for interfacing. !====================================== this%XX => NULL() #endif return end subroutine
Subroutine : | |||
phys_param : | type(quark_clover), intent(inout) | ||
phys_param : | type(quark_wilson), intent(inout) | ||
iout : | integer, intent(in)
| ||
err : | real(DP), intent(in)
| ||
iiter : | integer, intent(inout)
| ||
u : | type(vfield_gluon_wg), intent(in)
| ||
ieo : | integer, intent(in)
| ||
eval : | real(DP), intent(out)
| ||
shift : | real(DP), optional, intent(in)
|
Computes maximum eigenvalue lambda for
v
(Qee - shift)^(-1) v = lambda v.
and return eval for v Qee v = eval v with eval = 1/lambda + shift
subroutine get_eigen_Qee(phys_param,iout,err,iiter,u,ieo,eval,shift) ! ! Computes maximum eigenvalue lambda for ! v ! (Qee - shift)^(-1) v = lambda v. ! ! and return eval for ! v ! Qee v = eval v with eval = 1/lambda + shift ! use ks_class use timer_class implicit none #ifdef _CLOVER type(quark_clover), intent(inout) :: phys_param #else type(quark_wilson), intent(inout) :: phys_param #endif integer, intent(in) :: iout ! log file unit real(DP), intent(in) :: err ! eigenvalue precison goal integer, intent(inout) :: iiter ! iteration counter type(vfield_gluon_wg), intent(in) :: u ! gauge field integer, intent(in) :: ieo ! even/odd site index real(DP), intent(out) :: eval ! eigenvalue real(DP), optional, intent(in) :: shift ! shift value character(len=CHARLEN) :: str type(field_quark_eo_wg), allocatable :: we(:),ve(:),re type(ks_alg), allocatable :: eigen_solver type(timer) :: solver_time real(DP), allocatable :: sft(:) integer :: jiter,nnev,i logical :: is_shift real(DP) :: tol real(DP) :: tol_slv integer, parameter :: NEV=4, MEV=20 integer, parameter :: NSITE=COL*SPIN*NTH*NZ*NY*NX real(DP) :: eev(NEV),rtmp,rtmp2,etime logical :: eigen_debug call new(solver_time) call tic(solver_time) is_shift = present(shift) if (is_shift) then allocate(sft(1)) sft(1) = shift endif tol = err ! eigen solver tolerance tol_slv = err*1.0e-1_DP ! inversion solver tolerance allocate(we(1),ve(1)) call new(we(1),ieo) call new(ve(1),ieo) allocate(eigen_solver) iiter = 0 write(iout,'("# NEV=",I3)')NEV #ifdef _DEBUG eigen_debug = .TRUE. #else eigen_debug = .FALSE. #endif call new(eigen_solver,n=NSITE,m=MEV,NEV=nev,maxiter=400,tol=tol,debug=eigen_debug,mode=ESLV_MODE_HERMITE) do call solve(eigen_solver) select case(get_status(eigen_solver)) case (ESLV_OP_NOP) cycle case (ESLV_OP_DO_MATVEC) call unpack(eigen_solver%src_vec,we(1)) if (is_shift) then call assign_inv_mult_Qee_shift(phys_param,iout,tol_slv,jiter,ve,we(1),u,sft) else call assign_inv_mult_Qee_shift(phys_param,iout,tol_slv,jiter,ve,we(1),u) endif call pack(ve(1),eigen_solver%dst_vec) iiter = iiter + jiter case (ESLV_OP_CONVERGED) exit case (ESLV_OP_FAIL) write(str,'("Error in eigen solver: get_eigen_Qee, file:",A," line:",I5)')__FILE__,__LINE__ call error_stop(TRIM(str)) end select enddo if (0==nodeid) then do i=1,NEV write(iout,'("%EVS",I3,2E24.15)')i,eigen_solver%eval(i) enddo endif if (is_shift) then do i=1,NEV eev(i) = 1.0_DP/REAL(eigen_solver%eval(i),kind=DP)+sft(1) enddo eval = MINVAL(eev(:)) else do i=1,NEV eev(i) = 1.0_DP/REAL(eigen_solver%eval(i),kind=DP) enddo eval = MINVAL(eev(:)) endif if (0==nodeid) then do i=1,NEV write(iout,'("%EEV",I3,E24.15)')i,eev(i) enddo endif #ifdef _DEBUG !==================================================== ! check residual !==================================================== ! v v ! ve = Dee^+ Dee we(1) = Qee we(1) ! ! re = Qee we(1) - eval we(1) ! !==================================================== allocate(re) call unpack(eigen_solver%V(:,1),we(1)) call assign_mult_eoprec_wd( phys_param,ve(1),we(1),u) call assign_mult_eoprec_wd_dag(phys_param,re ,ve(1),u) rtmp = - eval call accum_add_mult(re,we(1),rtmp) ! re = re + we(1)*(-eval) rtmp2 = abs2(we(1)) rtmp = sqrt(abs2(re)/rtmp2/eval**2) deallocate(re) if (0 == nodeid) then write(*,'("%REL_ERR_=",E24.15)')rtmp write(*,'("%EVEC2_=",E24.15)')rtmp2 endif !==================================================== #endif deallocate(eigen_solver) deallocate(ve,we) if (is_shift) deallocate(sft) call toc(solver_time) etime = get_elapse(solver_time) if (0 == nodeid) then write(iout,'("%GET_EIGEN_ETIME_=",E24.16)')etime endif return end subroutine
Subroutine : | |||
this : | type(hmc_ratio_quark_wilson),
intent(inout)
| ||
ifi : | integer, intent(in)
| ||
u : | type(vfield_gluon_wg), intent(in)
| ||
action : | real(DP), intent(out)
|
Calc Hamil for Rational HMC action.
Action : Srhmc Srhmc = \sum_{j=1,Nf} ye(j)^+ R_n(Qee) ye(j) (even sites only) - Nf*TrLog[F]
where
R_n(x) = (1/x)^(1/2) = alp_0 + \sum_{j=1,Nratio} alp_j/(x-bet_j) v v Qee = Dee^+ Dee
Initial Generation:
ye = G_n(Qee) eta eta : Gaussina noise G_n(x) = (x)^(1/4) : precise rational approximation for quad root.
Fianal computation
Srhmc = |H_n(Qee)ye|^2, H_n(x) = (1/x)^(1/4) : precise rational approximation for inverse quad root.
subroutine hamil_wilson_quark_rhmc(this,ifi,u,action) ! ! Calc Hamil for Rational HMC action. ! ! Action : Srhmc ! ! Srhmc = \sum_{j=1,Nf} ye(j)^+ R_n(Qee) ye(j) (even sites only) ! ! - Nf*TrLog[F] ! ! where ! R_n(x) = (1/x)^(1/2) = alp_0 + \sum_{j=1,Nratio} alp_j/(x-bet_j) ! v v ! Qee = Dee^+ Dee ! ! Initial Generation: ! ! ye = G_n(Qee) eta ! ! eta : Gaussina noise ! ! G_n(x) = (x)^(1/4) : precise rational approximation for quad root. ! ! Fianal computation ! ! Srhmc = |H_n(Qee)ye|^2, ! ! H_n(x) = (1/x)^(1/4) : precise rational approximation for inverse quad root. ! use comlib implicit none type(hmc_ratio_quark_wilson), intent(inout) :: this ! quark action/algorithm paramters integer, intent(in) :: ifi ! target switch, ! * 0 for initial with generating PF field, ! * 1 for final action, ! * 2 for revcheck type(vfield_gluon_wg), intent(in) :: u ! gauge field real(DP), intent(out) :: action ! action value type(field_quark_eo_wg), allocatable :: we type(field_quark_eo_wg), allocatable :: ve integer :: nf,iflv real(DP) :: Spf,Sclv,kappa character(CHARLEN) :: str type(hmc_status) :: status integer :: iout,iiter real(DP) :: eval_min,eval_max,shift_max,shift_min,err,rtmp0,rtmp1,rtmp2 logical :: do_get_eigen complex(DP) :: ctmp0 nf = get_nflavor(this%phys_param) Spf =0.0_DP Sclv=0.0_DP #ifdef _CLOVER !=========================================== ! calc. F^-1 and ! actaion value of clover term trace log !=========================================== call hamil_clover_trlog(this%phys_param,Sclv) if (0 == this%depth_clvtrlog) Sclv = 0.0_DP #endif !======================================================== ! Compute eigenvalue of Q to get ! approximation iterval [min,max] !======================================================== call print(this%eigen_log,"%MIN and MAX EIGEN VALUES for Qee") iout = get_file_unit(this%eigen_log) if ( (.not.is_null(this%eigen_hamil_hist(0))) .and. (.not.is_null(this%eigen_hamil_hist(1))) .and. 0 == ifi) then ! ! if previous traj is accepted, then ! the initial eigen is identical to the previous last eigen, else ! is identical to the previous initial eigen. ! if (get_previous_status(status)) then ! accepted case call assign(this%eigen_hamil,this%eigen_hamil_hist(1)) else ! rejected case call assign(this%eigen_hamil,this%eigen_hamil_hist(0)) endif eval_min = get_min(this%eigen_hamil) eval_max = get_max(this%eigen_hamil) else err = 1.0e-4_DP iiter = 0 if ( is_null(this%eigen_hamil_hist(0)) ) then shift_min = get_min(this%eigen_force) shift_max = get_max(this%eigen_force) else shift_min = get_min(this%eigen_hamil_hist(0)) shift_max = get_max(this%eigen_hamil_hist(0)) endif shift_min = shift_min/1.1_DP shift_max = shift_max*1.1_DP call get_eigen_Qee(this%phys_param,iout,err,iiter,u,ieo=0,eval=eval_min,shift=shift_min) call get_eigen_Qee(this%phys_param,iout,err,iiter,u,ieo=0,eval=eval_max,shift=shift_max) call set(this%eigen_hamil,eval_min,eval_max) endif call assign(this%eigen_hamil_hist(ifi),this%eigen_hamil) ! save for eval history call assign(this%eigen_func_eval,this%eigen_hamil) ! for function evaluation this%eigen_union = union(this%eigen_union,this%eigen_hamil) ! to get largest eval range if (0 == nodeid) then write(iout,'("%EVAL_MIN_MAX:",I9,I2,2E12.4)')get_trajectory_number(status),ifi,eval_min,eval_max endif !!================== ! call check_rational_approx(this,iout,u) ! stop !!================== if ( (eval_min < get_min(this%eigen_force)) .or. (get_max(this%eigen_force) < eval_max) ) then write(str,'(" Force approximation range does not contains the evaluated range.")') call error_message(TRIM(str)) write(str,'(" Evaluated range [a,b]=[",F12.4,",",F12.4,"]")')eval_min,eval_max call error_message(TRIM(str)) write(str,'(" Force range [a,b]=[",F12.4,",",F12.4,"]")')get_min(this%eigen_force), get_max(this%eigen_force) call error_message(TRIM(str)) call error_stop("Adjust eigen_min_force and eigen_max_force") endif write(str,'("# HAMIL SOLVER")') call print(this%solver_log,TRIM(str)) iout = get_file_unit(this%solver_log) select case(ifi) case(0) ! INITIAL allocate(we) call new(we,0) Spf = 0.0_DP do iflv = 1,nf call set_gaussian_noise(we) !============================= ! calc. pseudo fermion action !============================= Spf = Spf + abs2(we) !============================ ! ! ye(iflv) = R_n(Q) we ! ! R_n(x) = (x)^(1/4) ! ! with good approximation ! ! This needs |Q| range ! !============================ call assign_mult_function_Qee(this,TYPE_QUAD_ROOT,iout,this%y(iflv)%eo(0),we,u) enddo deallocate(we) case(1:2) ! FINAL : REVCHECK allocate(we) call new(we,0) Spf = 0.0_DP do iflv = 1,nf !============================ ! ! we = R_n[Q] ye(iflv) ! ! R_n[x] = (1/x)^(1/4) ! ! with good approximation ! ! This needs |Q| range ! !============================ call assign_mult_function_Qee(this,TYPE_INV_QUAD_ROOT,iout,we,this%y(iflv)%eo(0),u) Spf = Spf + abs2(we) enddo deallocate(we) end select if (0 == this%depth_hopping) Spf = 0.0_DP #ifdef _CLOVER call delete_inverse_clover_term(this%phys_param) #endif action=Spf+Sclv write(str,'("@",I8,I2," Quark: ID:",I3," SQPF:",E24.16," SCLV:",E24.16)') get_trajectory_number(status),ifi,this%myid,Spf,Sclv call print_log_action(status,TRIM(str)) return end subroutine
Function : | |
has_gmp : | logical |
this : | type(hmc_ratio_quark_wilson), intent(in) |
inquire whether this algorithm has GMP
function has_gmp_wilson_quark_rhmc(this) result(has_gmp) ! ! inquire whether this algorithm has GMP ! implicit none type(hmc_ratio_quark_wilson), intent(in) :: this logical :: has_gmp has_gmp = this%has_global_metropolis_test return end function
Function : | |
has_rew : | logical |
this : | type(hmc_ratio_quark_wilson), intent(in) |
inquire whether this algorithm has reweighting
function has_reweighting_wilson_quark_rhmc(this) result(has_rew) ! ! inquire whether this algorithm has reweighting ! implicit none type(hmc_ratio_quark_wilson), intent(in) :: this logical :: has_rew has_rew = this%has_reweighting return end function
Derived Type : | |||
phys_param : | type(quark_clover), public | ||
phys_param : | type(quark_wilson), public | ||
XX => NULL() : | type(tfield_gluon_wg), public, pointer | ||
phys_param : | type(quark_clover), public | ||
phys_param : | type(quark_wilson), public | ||
eigen_force : | type(eigen_value_range)
| ||
eigen_hamil : | type(eigen_value_range)
| ||
eigen_hamil_hist(0:2) : | type(eigen_value_range)
| ||
eigen_union : | type(eigen_value_range)
| ||
eigen_func_eval : | type(eigen_value_range)
|
HMC quark even/odd site preconditioned Rational HMC algorithm action
Subroutine : | |||
this : | type(hmc_ratio_quark_wilson),
intent(inout)
| ||
u1 : | type(vfield_gluon_wg), intent(in)
| ||
u0 : | type(vfield_gluon_wg), intent(in)
| ||
ihit : | integer, intent(inout)
|
Nothing is done
subroutine metropolis_wilson_quark_rhmc(this,u1,u0,ihit) ! ! Nothing is done ! implicit none type(hmc_ratio_quark_wilson), intent(inout) :: this ! quark action/algorihtm parameters type(vfield_gluon_wg), intent(in) :: u1 ! MD final gauge field type(vfield_gluon_wg), intent(in) :: u0 ! MD inital gauge field integer, intent(inout) :: ihit ! 1 for accept, 0 for reject return end subroutine
Subroutine : | |||
this : | type(hmc_ratio_quark_wilson), intent(inout) | ||
id : | integer, intent(in)
|
initialize wilson/clover quark RHMC algorithm
subroutine new_wilson_quark_rhmc(this,id) ! ! initialize wilson/clover quark RHMC algorithm ! implicit none type(hmc_ratio_quark_wilson), intent(inout) :: this integer, intent(in) :: id ! quark id character(len=CHARLEN) :: str this%myid = id if (allocated(this%y)) deallocate(this%y) call new(this%phys_param,id) call new(this%hamil_inv) call new(this%force_inv) call new(this%eigen_force) call new(this%eigen_hamil) call new(this%eigen_hamil_hist(0)) call new(this%eigen_hamil_hist(1)) call new(this%eigen_hamil_hist(2)) call new(this%eigen_union) write(str,'(I2.2)')id this%solver_log_fname = TRIM(this%solver_log_fname)//TRIM(str) this%eigen_log_fname = TRIM(this%eigen_log_fname)//TRIM(str) call new(this%solver_log,this%solver_log_fname) call new(this%eigen_log,this%eigen_log_fname) return end subroutine
Subroutine : | |
this : | type(hmc_ratio_quark_wilson), intent(inout) |
print out quark RHMC algorithm parameters on display
subroutine print_wilson_quark_rhmc(this) ! ! print out quark RHMC algorithm parameters on display ! implicit none type(hmc_ratio_quark_wilson), intent(inout) :: this integer :: imb if (nodeid==0) then write(*,'("**** Quark # ",I3," RHMC")')this%myid write(*,'("***** Physical parameters *****")') call print(this%phys_param) write(*,'("***** RHMC Algorithm parameters *****")') write(*,'(" Hamiltonian Stopping Condition :",E16.8," (pseudo-fermion)")')this%hamil_inv%tol write(*,'(" Force Stopping Condition :",E16.8," (pseudo-fermion)")')this%force_inv%tol write(*,'(" MD force depth PSfermion part :",I3)')this%depth_hopping #ifdef _CLOVER write(*,'(" MD force depth CLV trlog part :",I3)')this%depth_clvtrlog #endif write(*,'("Rational Approx Order for Force :",I4)') this%nratio_force write(*,'("Rational Approx Range for Force :")') write(*,'(" Qee in [a,b]=[",E24.15,",",E24.15,"]")') get_min(this%eigen_force), get_max(this%eigen_force) write(*,'(80("-"))') endif return end subroutine
Subroutine : | |
this : | type(hmc_ratio_quark_wilson), intent(inout) |
print out quark RHMC algorithm job statistics on display
subroutine print_stat_wilson_quark_rhmc(this) ! ! print out quark RHMC algorithm job statistics on display ! implicit none type(hmc_ratio_quark_wilson), intent(inout) :: this character(CHARLEN) :: name real(DP) :: hamil_ave,force_ave hamil_ave = get_average(this%hamil_inv) force_ave = get_average(this%force_inv) if (nodeid==0) then write(*,'("**** Quark # ",I3," RHMC")')this%myid call print(this%phys_param) write(*,'(" RHMC algorithm ")') write(*,'(" Hamil CG iter (averaged) :",F12.4)') hamil_ave write(*,'(" Force CG iter (averaged) :",F12.4)') force_ave write(*,'(" Minimum Eigenvalue for Qee :",E24.15)') get_min(this%eigen_union) write(*,'(" Maximum Eigenvalue for Qee :",E24.15)') get_max(this%eigen_union) endif return end subroutine
Subroutine : | |
this : | type(hmc_ratio_quark_wilson), intent(inout) |
iout : | integer, intent(in) |
Read RHMC quark algorithm parameters and Initialize algorithm parameters.
subroutine read_wilson_quark_rhmc(this,iout) ! ! Read RHMC quark algorithm parameters and Initialize algorithm parameters. ! use comlib implicit none type(hmc_ratio_quark_wilson), intent(inout) :: this integer, intent(in) :: iout real(DP) :: rtmp0 real(DP) :: eigen_min_force real(DP) :: eigen_max_force integer :: imb call new(this%hamil_inv) call new(this%force_inv) call read(this%phys_param,iout) if (nodeid==0) then read(iout,*) read(iout,*)this%hamil_inv%tol read(iout,*)this%force_inv%tol read(iout,*)this%depth_hopping read(iout,*)this%depth_clvtrlog read(iout,*)this%approx_range_factor_hamil read(iout,*)this%tol_approx_hamil read(iout,*)eigen_min_force read(iout,*)eigen_max_force read(iout,*)this%nratio_force endif #ifndef _singlePU call comlib_bcast(this%hamil_inv%tol,0) call comlib_bcast(this%force_inv%tol,0) call comlib_bcast(this%depth_hopping,0) call comlib_bcast(this%depth_clvtrlog,0) call comlib_bcast(this%approx_range_factor_hamil) call comlib_bcast(this%tol_approx_hamil) call comlib_bcast(eigen_min_force) call comlib_bcast(eigen_max_force) call comlib_bcast(this%nratio_force) #endif call set(this%eigen_force,eigen_min_force,eigen_max_force) this%nflavor = get_nflavor(this%phys_param) if (allocated(this%y)) deallocate(this%y) allocate(this%y(this%nflavor)) ! allocate multiple flavor return end subroutine
Subroutine : | |||
this : | type(hmc_ratio_quark_wilson),
intent(inout)
| ||
u : | type(vfield_gluon_wg), intent(in) |
Nothing is done
subroutine reweighting_wilson_quark_rhmc(this,u) ! ! Nothing is done ! implicit none type(hmc_ratio_quark_wilson), intent(inout) :: this ! quark action/algorihtm parameters type(vfield_gluon_wg), intent(in) :: u return end subroutine
Subroutine : | |
this : | type(hmc_ratio_quark_wilson), intent(in) |
iout : | integer, intent(in) |
Save Quark parameter on measurement configuration
subroutine save_config_hmc_ratio_quark_wilson(this,iout) ! ! Save Quark parameter on measurement configuration ! implicit none type(hmc_ratio_quark_wilson), intent(in) :: this integer, intent(in) :: iout call save_config(this%phys_param,iout) return end subroutine
Subroutine : | |
this : | type(hmc_ratio_quark_wilson), intent(inout) |
id : | integer, intent(in) |
set quark # id
subroutine set_id_hrqw(this,id) ! ! set quark # id ! implicit none type(hmc_ratio_quark_wilson), intent(inout) :: this integer, intent(in) :: id this%myid = id call set_id(this%phys_param,id) return end subroutine