Class | hmc_blk_quark_wilson_class |
In: |
HmcBlkQuarkWilsonClass/hmc_blk_quark_wilson_class.F90
|
used for even-number degenrate-mass flavor simulations
$Id: hmc_blk_quark_wilson_class.F90,v 1.7 2011/11/14 09:26:18 ishikawa Exp ishikawa $
Subroutine : | |||
this : | type(quark_clover), intent(in), target
| ||
this : | type(quark_wilson), intent(in), target
| ||
iout : | integer, intent(in)
| ||
err : | real(DP), intent(in)
| ||
iiter : | integer, intent(out)
| ||
xe(:) : | type(field_quark_eo_wg), intent(inout)
| ||
be(:) : | type(field_quark_eo_wg), intent(in)
| ||
u : | type(vfield_gluon_wg), intent(in), target
| ||
cron(:) : | type(chrolog_alg), intent(inout), optional |
Multiply inverse even/odd-site preconditioned WIlson-Dirac oparator on fermion field. by solving:
D x = b -> x = D^-1 b
ifdef GPU_SOLVER
subroutine assign_inv_mult_eoprec_blk_wd_eo(this,iout,err,iiter,xe,be,u,cron) ! ! Multiply inverse even/odd-site preconditioned WIlson-Dirac oparator on fermion field. ! by solving: ! ! D x = b -> x = D^-1 b ! ! #ifdef _GPU_SOLVER_ use gpu_block_solver_class #endif use blbicggr_mod use chrolog_class use counter_class use timer_class implicit none #ifdef _CLOVER type(quark_clover), intent(in), target :: this ! quark parameters #else type(quark_wilson), intent(in), target :: this ! quark parameters #endif integer, intent(in) :: iout ! solver logfile unit number real(DP), intent(in) :: err ! solver stopping condition integer, intent(out) :: iiter ! total iteration count after convergence type(field_quark_eo_wg), intent(inout) :: xe(:) ! even/odd site solution field type(field_quark_eo_wg), intent(in) :: be(:) ! even/odd site source(rhs) field type(vfield_gluon_wg), intent(in), target :: u ! gauge field type(chrolog_alg), intent(inout), optional :: cron(:) complex(DP), allocatable :: xx(:,:),yy(:,:),vec(:) type(blbicggr) :: solver type(field_quark_eo_wg), allocatable :: we #ifdef _GPU_SOLVER_ type(gpu_block_solver) :: gpu_precnd #endif logical :: is_cron type(counter) :: imult integer :: istat, mode_solver, iter character(len=CHARLEN) :: str type(timer) :: solver_time real(DP) :: etime,rtmp0,rtmp1 real(DP), allocatable :: res_norm(:) integer :: nblk,nblk_x,nblk_cron,iblk, enc_size,ieo nblk = SIZE(be) nblk_x = SIZE(xe) is_cron = .false. if (present(cron)) then do iblk = 1,nblk if (get_ncron(cron(iblk)) > 0) then is_cron = .true. endif enddo endif if (is_cron) then nblk_cron = SIZE(cron) if ( nblk /= nblk_x .or. nblk /= nblk_cron) then write(str,'(A," at ",I5," : assign_inv_mult_eoprec_blk_wd_eo : block size is inconsisntent.")') __FILE__,__LINE__ call error_message(str) write(str,'(" size(be)=",I3," size(xe)=",I3," size(cron)=",I3)')nblk,nblk_x,nblk_cron call error_stop(str) endif else if ( nblk /= nblk_x ) then write(str,'(A," at ",I5," : assign_inv_mult_eoprec_blk_wd_eo : block size is inconsisntent.")') __FILE__,__LINE__ call error_message(str) write(str,'(" size(be)=",I3," size(xe)=",I3)')nblk,nblk_x call error_stop(str) endif endif call new(solver_time) call tic(solver_time) allocate(we) call new(we,be(1)%ieo) do iblk=1,nblk call new(xe(iblk),be(iblk)%ieo) enddo #ifdef _DEBUG do iblk=1,nblk rtmp0 = abs2(be(iblk)) if (0 == nodeid) then write(*,'("%%INV_MULT_EOPREC_WD IEO=",I2)')be(iblk)%ieo write(*,'("%% SOURCE |Be|^2 = ",E24.15,I3)')rtmp0,iblk endif enddo #endif imult = mult_iter !======================================================================= ! Chronological guess solver with Reverse Communication Interface (RCI) !======================================================================= if (is_cron) then do iblk = 1,nblk if (get_ncron(cron(iblk)) > 0) then if (nodeid == 0) write(iout,'("# CRON USED")') endif call pack(be(iblk),cron(iblk)%src_vec) do call solve(cron(iblk)) select case (get_status(cron(iblk))) case (CHRON_OP_NOP) cycle case (CHRON_OP_DO_MATVEC) call unpack(cron(iblk)%src_vec,xe(iblk)) call assign_mult_eoprec_wd(this,we,xe(iblk),u) call pack(we,cron(iblk)%dst_vec) case (CHRON_OP_CONVERGED) exit end select enddo if (get_ncron(cron(iblk)) > 0) then if (nodeid == 0) then write(iout,'("#",I5," ERR_=",E24.16," MULT_=",I10," NPF_=",I4)') 0, get_residual_norm(cron(iblk)),get_count(mult_iter)-get_count(imult),iblk endif endif enddo #ifdef _GPU_SOLVER_ call new(solver,NSIZE=NSITES,NBLK=nblk,tol=err,guess=GUESS_USE,mode=MODE_PRECOND) #else call new(solver,NSIZE=NSITES,NBLK=nblk,tol=err,guess=GUESS_USE) #endif do iblk=1,nblk !$OMP PARALLEL WORKSHARE solver%dst_vec(iblk,:) = cron(iblk)%dst_vec(:) !$OMP END PARALLEL WORKSHARE enddo else #ifdef _GPU_SOLVER_ call new(solver,NSIZE=NSITES,NBLK=nblk,tol=err,guess=GUESS_NO,mode=MODE_PRECOND) #else call new(solver,NSIZE=NSITES,NBLK=nblk,tol=err,guess=GUESS_NO) #endif endif #ifdef _GPU_SOLVER_ gpu_precnd%num_pf = nblk gpu_precnd%tol = err #ifdef _CLOVER call new(gpu_precnd,u,be(1)%ieo,get_ptr_inverse_clover_term(this)) #else call new(gpu_precnd,u,be(1)%ieo) #endif #endif !===================================== ! Bl-BiCGGR solver !===================================== !----------------------------- ! block source set up !----------------------------- allocate(res_norm(nblk)) allocate(vec(NSITES)) do iblk=1,nblk call pack(be(iblk),vec) !$OMP PARALLEL WORKSHARE solver%src_vec(iblk,:) = vec(:) !$OMP END PARALLEL WORKSHARE enddo do call solve(solver) select case(get_status(solver)) case(OP_NOP) cycle case(OP_DO_MATVEC) do iblk=1,nblk !$OMP PARALLEL WORKSHARE vec(:) = solver%src_vec(iblk,:) !$OMP END PARALLEL WORKSHARE call unpack(vec,xe(iblk)) call assign_mult_eoprec_wd(this,we,xe(iblk),u) call pack(we,vec) !$OMP PARALLEL WORKSHARE solver%dst_vec(iblk,:) = vec(:) !$OMP END PARALLEL WORKSHARE enddo case(OP_PRECOND) #ifdef _GPU_SOLVER_ gpu_precnd%tol = err res_norm = get_residual_norm(solver) gpu_precnd%err_norm = MAXVAL(res_norm) call assign_inv_mult_eoprec_wd(gpu_precnd,get_kappa(this),solver%dst_vec,solver%src_vec) #endif case(OP_PRINT_STATUS) #ifdef _DEBUG iiter = get_current_iteration(solver) #ifdef _GPU_SOLVER_ iiter = iiter + gpu_precnd%iter #endif if (nodeid == 0) then res_norm(:) = get_residual_norm(solver) do iblk=1,nblk write(iout,'("#",I5," ERR_=",E24.16," MULT_=",I10," NPF_=",I4)') iiter/nblk,res_norm(iblk),get_count(mult_iter)-get_count(imult),iblk write(* ,'("#",I5," ERR_=",E24.16," MULT_=",I10," NPF_=",I4)') iiter/nblk,res_norm(iblk),get_count(mult_iter)-get_count(imult),iblk enddo endif #endif case(OP_CONVERGED) iiter = get_current_iteration(solver) #ifdef _GPU_SOLVER_ iiter = iiter + gpu_precnd%iter #endif if (nodeid == 0) then res_norm(:) = get_residual_norm(solver) write(iout,'("#",A," iteration converged.")')TRIM(get_name(solver)) do iblk=1,nblk write(iout,'("#",I5," ERR_=",E24.16," MULT_=",I10," NPF_=",I4)',advance='no') iiter/nblk,res_norm(iblk),get_count(mult_iter)-get_count(imult),iblk if (iblk < nblk) write(iout,*) enddo #ifdef _DEBUG write(* ,'("#",A," iteration converged.")')TRIM(get_name(solver)) do iblk=1,nblk write(* ,'("#",I5," ERR_=",E24.16," MULT_=",I10," NPF_=",I4)') iiter/nblk,res_norm(iblk),get_count(mult_iter)-get_count(imult),iblk enddo #endif endif exit case(OP_MAXITER_REACHED) iiter = get_current_iteration(solver) #ifdef _GPU_SOLVER_ iiter = iiter + gpu_precnd%iter #endif if (nodeid == 0) then res_norm = get_residual_norm(solver) do iblk=1,nblk write(iout,'("#",A," iteration does not converge.")')TRIM(get_name(solver)) write(iout,'("#",I5," ERR_=",E24.16," MULT_=",I10," NPF_=",I4)') iiter/nblk,res_norm(iblk),get_count(mult_iter)-get_count(imult),iblk enddo endif write(str,'("solver did not converge.: ",A,I5)')__FILE__,__LINE__ call error_stop(str) end select enddo #ifdef _GPU_SOLVER_ call delete(gpu_precnd) #endif !----------------------------- ! get block solution !----------------------------- do iblk=1,nblk !$OMP PARALLEL WORKSHARE vec(:) = solver%dst_vec(iblk,:) !$OMP END PARALLEL WORKSHARE call unpack(vec,xe(iblk)) enddo deallocate(vec) iiter = get_current_iteration(solver)/nblk if (is_cron) then do iblk=1,nblk !$OMP PARALLEL WORKSHARE cron(iblk)%store_vec(:) = solver%dst_vec(iblk,:) !$OMP END PARALLEL WORKSHARE call store(cron(iblk)) enddo endif #ifdef _DEBUG write(*,*) do iblk=1,nblk rtmp0 = abs2(xe(iblk)) if (0 == nodeid) then write(*,'("%% SOLUTION |Dee \ Be|^2 = ",E24.15," NPF_=",I4)')rtmp0,iblk endif enddo do iblk=1,nblk ! ! compute true residual ! v ! we = Dee xe - be ! call assign_mult_eoprec_wd(this,we,xe(iblk),u) call accum_sub(we,be(iblk)) rtmp0 = abs2(be(iblk)) rtmp1 = abs2(we) rtmp0 = sqrt(rtmp1/rtmp0) if (0 == nodeid) then write(*,'("%",I5," ERR_=",E24.16," MULT_=",I10," NPF_=",I4)',advance='no') iiter,rtmp0,get_count(mult_iter)-get_count(imult),iblk if (iblk<nblk) write(*,*) endif enddo #endif call delete(solver) deallocate(we) deallocate(res_norm) call toc(solver_time) etime = get_elapse(solver_time) if (nodeid == 0) then write(iout,'(" ETIME_=",E24.16)')etime #ifdef _DEBUG write( *,'(" ETIME_=",E24.16)')etime #endif endif return end subroutine
Subroutine : | |
this : | type(hmc_blk_quark_wilson), intent(inout) |
delete HMC quark even/odd preconditioned action algorithm
subroutine delete_wilson_quark_blkhmc(this) ! ! delete HMC quark even/odd preconditioned action algorithm ! implicit none type(hmc_blk_quark_wilson), intent(inout) :: this deallocate(this%y) call delete(this%phys_param) call delete(this%solver_log) return end subroutine
Subroutine : | |||
this : | type(hmc_blk_quark_wilson),
intent(inout)
| ||
BB : | type(vfield_gluon_wg), intent(inout)
| ||
idepth : | integer, intent(in)
| ||
u : | type(vfield_gluon_wg), intent(in)
|
Calc force from hopping part of HMC quark even/odd-site preconditioned pseudo fermion
Action : S_Q v S_Q : Pseudo-fermion HMC action, |(Doo)^-1 yo|^2 (odd site only) v v where Doo = 1 - Moo v Moo = kappa^2 Fo^-1 Moe Fe^-1 Meo
subroutine force_wilson_quark_blkhmc(this,BB,idepth,u) ! ! Calc force from hopping part of HMC quark even/odd-site preconditioned pseudo fermion ! ! Action : S_Q ! v ! S_Q : Pseudo-fermion HMC action, |(Doo)^-1 yo|^2 (odd site only) ! ! v v ! where Doo = 1 - Moo ! v ! Moo = kappa^2 Fo^-1 Moe Fe^-1 Meo ! ! implicit none type(hmc_blk_quark_wilson), intent(inout) :: this ! quark action/algorithm parameters type(vfield_gluon_wg), intent(inout) :: BB ! pre-force contribution from hopping part integer, intent(in) :: idepth ! MD depth for recursive MD integrator type(vfield_gluon_wg), intent(in) :: u ! gauge field type(field_quark_wg), allocatable :: fx(:),fy(:) type(field_quark_eo_wg), allocatable :: xo(:) character(len=CHARLEN) :: str real(DP) :: kappa2 real(DP) :: kappa real(DP) :: err integer :: iout,ipf integer :: iiter #ifdef _CLOVER real(DP) :: csw,dummy complex(DP) :: zickd8 #endif #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 std hmc hopping")')idepth #endif write(str,'("# FORCE SOLVER")') call print(this%solver_log,TRIM(str)) iout = get_file_unit(this%solver_log) allocate(fx(this%num_pf),fy(this%num_pf),xo(this%num_pf)) !************************************************************************ !* Prepare gauge force contributions from pseudo ferimon !* !* v !* Qoo == Doo !* !************************************************************************ !************************************** !* fyo <= (Qoo)^-1 yo !************************************** err = this%force_inv%tol call assign_inv_mult_eoprec_wd(this%phys_param,iout,err,iiter,fy(:)%eo(1),this%y(:)%eo(1),u,this%cron1(:)) call inc(this%force_inv,iiter) ! do ipf=1,this%num_pf ! err = this%force_inv%tol ! call assign_inv_mult_eoprec_wd(this%phys_param,iout,err,iiter,fy(ipf)%eo(1),this%y(ipf)%eo(1),u,this%cron1(ipf)) ! call inc(this%force_inv,iiter) ! enddo !************************************** !* v !* fye <= Meo fyo !* !************************************** do ipf=1,this%num_pf call assign_mult_hop(this%phys_param,fy(ipf)%eo(0),fy(ipf)%eo(1),u) enddo !************************************** !* xo <= (Qoo)^-1 gamma_5 Fo^-1 fyo !************************************** do ipf=1,this%num_pf call assign_mult_gamma5(fx(ipf)%eo(1),fy(ipf)%eo(1)) enddo #ifdef _CLOVER do ipf=1,this%num_pf call accum_mult_clover_term(fx(ipf)%eo(1),this%phys_param) enddo #endif err = this%force_inv%tol call assign_inv_mult_eoprec_wd(this%phys_param,iout,err,iiter,xo(:),fx(:)%eo(1),u,this%cron2(:)) call inc(this%force_inv,iiter) ! do ipf=1,this%num_pf ! err = this%force_inv%tol ! call assign_inv_mult_eoprec_wd(this%phys_param,iout,err,iiter,xo(ipf),fx(ipf)%eo(1),u,this%cron2(ipf)) ! call inc(this%force_inv,iiter) ! enddo !************************************** !* !* fxo <= gamma_5 xo !* v !* fxe <= gamma_5 Meo xo !* !************************************** do ipf=1,this%num_pf call assign_mult_gamma5(fx(ipf)%eo(1),xo(ipf)) enddo do ipf=1,this%num_pf call assign_mult_hop(this%phys_param,fx(ipf)%eo(0),xo(ipf),u) enddo do ipf=1,this%num_pf call accum_mult_gamma5(fx(ipf)%eo(0)) enddo do ipf=1,this%num_pf call copy_boundary(fy(ipf)%eo(0)) call copy_boundary(fx(ipf)%eo(0)) call copy_boundary(fx(ipf)%eo(1)) enddo !====================================== ! force from hopping matrix !====================================== kappa = get_kappa(this%phys_param) kappa2 = kappa**2 do ipf=1,this%num_pf call force_hmc_hopping(BB,kappa2,fx(ipf),fy(ipf)) enddo #ifdef _CLOVER !************************************** !* v !* fyo <= Moe fye !* !************************************** do ipf=1,this%num_pf call assign_mult_hop(this%phys_param,fy(ipf)%eo(1),fy(ipf)%eo(0),u) enddo !====================================== ! force from clover term !====================================== csw = get_csw(this%phys_param) zickd8 = zi*csw*kappa**3/8.0_DP do ipf=1,this%num_pf call force_hmc_clover(this%XX,zickd8,fx(ipf),fy(ipf)) enddo #endif deallocate(fx,fy,xo) endif #ifdef _CLOVER if (idepth == this%depth_clvtrlog) then #ifdef _MD_DEPTH_ if (nodeid==0) write(*,'("depth=",I3," : force quark std hmc 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 ! this is used for interfacing !================================== this%XX => NULL() #endif return end subroutine
Subroutine : | |||
this : | type(hmc_blk_quark_wilson),
intent(inout)
| ||
ifi : | integer, intent(in)
| ||
u : | type(vfield_gluon_wg), intent(in)
| ||
action : | real(DP), intent(out)
|
Calculate Hamiltonian for HMC quark even/odd-site preoconditioned action.
Action : v Shmc = Sum_{j=1,Nf/2}|Doo^-1 yo| - Nf log[det[F]] (yo : lives only on odd field)
subroutine hamil_wilson_quark_blkhmc(this,ifi,u,action) ! ! Calculate Hamiltonian for HMC quark even/odd-site preoconditioned action. ! ! Action : ! v ! Shmc = Sum_{j=1,Nf/2}|Doo^-1 yo| - Nf log[det[F]] ! ! (yo : lives only on odd field) ! use comlib implicit none type(hmc_blk_quark_wilson), intent(inout) :: this ! quark action/algorithm parameters integer, intent(in) :: ifi ! target switch : 0 for initial generating PF field, 1 for final type(vfield_gluon_wg), intent(in) :: u ! gauge field real(DP), intent(out) :: action ! output action value real(DP) :: Spf,Sclv type(hmc_status) :: status type(field_quark_eo_wg), allocatable :: wo,vo(:) character(CHARLEN) :: str real(DP) :: err integer :: iiter,iout,ipf 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 select case(ifi) case(0) ! INITIAL allocate(this%cron1(this%num_pf)) allocate(this%cron2(this%num_pf)) do ipf=1,this%num_pf call new(this%cron1(ipf),NSITES,this%ncron) call new(this%cron2(ipf),NSITES,this%ncron) enddo allocate(wo) call new(wo,1) Spf = 0.0_DP do ipf=1,this%num_pf !===================== ! Heatbath !===================== call set_gaussian_noise(wo) !===================== ! calc. action !===================== Spf=Spf + abs2(wo) !===================== ! v ! yo = Doo wo !===================== call assign_mult_eoprec_wd(this%phys_param,this%y(ipf)%eo(1),wo,u) if (this%ncron > 0) then call pack(wo,this%cron1(ipf)%store_vec) call store(this%cron1(ipf)) endif enddo deallocate(wo) case(1:2) ! FINAL:REVCHECK write(str,'("# HAMIL SOLVER")') call print(this%solver_log,TRIM(str)) iout = get_file_unit(this%solver_log) ! allocate(wo) ! call new(wo,1) allocate(vo(this%num_pf)) do ipf=1,this%num_pf call new(vo(ipf),1) enddo Spf = 0.0_DP !============================ ! v ! solve : wo = Doo^-1 yo !============================ err = this%hamil_inv%tol call assign_inv_mult_eoprec_wd(this%phys_param,iout,err,iiter,vo(:),this%y(:)%eo(1),u,this%cron1(:)) call inc(this%hamil_inv,iiter) ! do ipf=1,this%num_pf ! err = this%hamil_inv%tol ! call assign_inv_mult_eoprec_wd(this%phys_param,iout,err,iiter,wo,this%y(ipf)%eo(1),u,this%cron1(ipf)) ! call inc(this%hamil_inv,iiter) ! enddo do ipf=1,this%num_pf !=================================== ! action = |Doo^-1 yo|^2 = |wo|^2 !=================================== ! Spf = Spf + abs2(wo) Spf = Spf + abs2(vo(ipf)) enddo ! deallocate(wo) deallocate(vo) do ipf=1,this%num_pf call delete(this%cron1(ipf)) call delete(this%cron2(ipf)) enddo deallocate(this%cron1) deallocate(this%cron2) #ifdef _REVCHECK allocate(this%cron1(this%num_pf)) allocate(this%cron2(this%num_pf)) do ipf=1,this%num_pf call new(this%cron1(ipf),NSITES,this%ncron) call new(this%cron2(ipf),NSITES,this%ncron) enddo #endif 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_blk_quark_wilson), intent(in) |
inquire whether this algorithm has GMP method
function has_gmp_wilson_quark_blkhmc(this) result(has_gmp) ! ! inquire whether this algorithm has GMP method ! implicit none type(hmc_blk_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_blk_quark_wilson), intent(in) |
inquire whether this algorithm has reweighting
function has_reweighting_wilson_quark_blkhmc(this) result(has_rew) ! ! inquire whether this algorithm has reweighting ! implicit none type(hmc_blk_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 |
ifdef _CLOVER
Subroutine : | |||
this : | type(hmc_blk_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_blkhmc(this,u1,u0,ihit) ! ! Nothing is done ! implicit none type(hmc_blk_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_blk_quark_wilson), intent(inout) |
id : | integer, intent(in) |
initialize/create HMC quark even/odd preconditioned action algorithm
subroutine new_wilson_quark_blkhmc(this,id) ! ! initialize/create HMC quark even/odd preconditioned action algorithm ! implicit none type(hmc_blk_quark_wilson), intent(inout) :: this integer, intent(in) :: id character(len=CHARLEN) :: str if (allocated(this%y)) deallocate(this%y) this%myid = id call new(this%phys_param,id) call new(this%force_inv) call new(this%hamil_inv) this%depth_hopping = 1 this%depth_clvtrlog = 1 write(str,'(I2.2)')id this%solver_log_fname = TRIM(ADJUSTL(this%solver_log_fname))//TRIM(ADJUSTL(str)) call new(this%solver_log,TRIM(ADJUSTL(this%solver_log_fname))) return end subroutine
Subroutine : | |
this : | type(hmc_blk_quark_wilson), intent(inout) |
print out HMC quark even/odd preconditioned algorithm parameters on disyplay.
subroutine print_wilson_quark_blkhmc(this) ! ! print out HMC quark even/odd preconditioned algorithm parameters on disyplay. ! implicit none type(hmc_blk_quark_wilson), intent(inout) :: this if (nodeid==0) then write(*,'("**** Quark # ",I3," HMC Blocked")')this%myid write(*,'("***** Physical parameters *****")') call print(this%phys_param) write(*,'("***** HMC Blocked Algorithm parameters *****")') write(*,'(" Number of Pseudo-Fermions :",I3)')this%num_pf write(*,'(" Force Stopping Condition :",E16.8)')this%force_inv%tol write(*,'(" Hamil Stopping Condition :",E16.8)')this%hamil_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(*,'(" Chronological History :",I4)') this%ncron write(*,'(80("-"))') endif return end subroutine
Subroutine : | |
this : | type(hmc_blk_quark_wilson), intent(inout) |
Print out HMC quark even/odd-site preconditioned algorithm job statistics on display
subroutine print_stat_wilson_quark_blkhmc(this) ! ! Print out HMC quark even/odd-site preconditioned algorithm job statistics on display ! implicit none type(hmc_blk_quark_wilson), intent(inout) :: this real(DP) :: f_ave,h_ave f_ave = get_average(this%force_inv) h_ave = get_average(this%hamil_inv) if (nodeid==0) then write(*,'("**** Quark # ",I3," HMC Blocked")')this%myid call print(this%phys_param) write(*,'(" HMC algorithm ")') write(*,'(" Force BiCGStab iter")') write(*,'(" (averaged) :",F12.4)') f_ave write(*,'(" Hamiltonian BiCGStab iter")') write(*,'(" (averaged) :",F12.4)') h_ave write(*,'(" Chronological History :",I4)') this%ncron endif return end subroutine
Subroutine : | |
this : | type(hmc_blk_quark_wilson), intent(inout) |
iout : | integer, intent(in) |
Read HMC quark even/odd-site preoconditioned algorithm parameters
subroutine read_wilson_quark_blkhmc(this,iout) ! ! Read HMC quark even/odd-site preoconditioned algorithm parameters ! use comlib implicit none type(hmc_blk_quark_wilson), intent(inout) :: this integer, intent(in) :: iout integer :: nf,ipf character(len=CHARLEN) :: str 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 #ifdef _CLOVER read(iout,*)this%depth_clvtrlog #endif read(iout,*)this%ncron 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) #ifdef _CLOVER call comlib_bcast(this%depth_clvtrlog,0) #endif call comlib_bcast(this%ncron,0) #endif nf = get_nflavor(this%phys_param) if (mod(nf,2) == 0) then this%num_pf = nf/2 else write(str,'("HMC Blocked algorithm works only for even-number flavors. Nf=",I3)')nf call error_stop(str) endif allocate(this%y(this%num_pf)) do ipf=1,this%num_pf call new(this%y(ipf)) enddo return end subroutine
Subroutine : | |||
this : | type(hmc_blk_quark_wilson),
intent(inout)
| ||
u : | type(vfield_gluon_wg), intent(in) |
Nothing is done
subroutine reweighting_wilson_quark_blkhmc(this,u) ! ! Nothing is done ! implicit none type(hmc_blk_quark_wilson), intent(inout) :: this ! quark action/algorihtm parameters type(vfield_gluon_wg), intent(in) :: u return end subroutine
Subroutine : | |||
this : | type(hmc_blk_quark_wilson), intent(in) | ||
iout : | integer, intent(in)
|
Save HMC Quark even/odd-site preoconditioned algorihtm parameters on measurement configuration
subroutine save_config_hmc_blk_quark_wilson(this,iout) ! ! Save HMC Quark even/odd-site preoconditioned algorihtm parameters on measurement configuration ! implicit none type(hmc_blk_quark_wilson), intent(in) :: this integer, intent(in) :: iout ! configuration file unit id call save_config(this%phys_param,iout) return end subroutine
Subroutine : | |
this : | type(hmc_blk_quark_wilson), intent(inout) |
id : | integer, intent(in) |
set quark # id
subroutine set_id_hbqw(this,id) ! ! set quark # id ! implicit none type(hmc_blk_quark_wilson), intent(inout) :: this integer, intent(in) :: id this%myid = id call set_id(this%phys_param,id) return end subroutine