Class | hmc_std_quark_wilson_class |
In: |
HmcStdQuarkWilsonClass/hmc_std_quark_wilson_class.F90
|
$Id: hmc_std_quark_wilson_class.F90,v 1.32 2011/09/28 07:16:16 ishikawa Exp $
Subroutine : | |||
this : | type(quark_clover), intent(in)
| ||
this : | type(quark_wilson), intent(in)
| ||
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)
| ||
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_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_solver_class ! use single_solver_class #endif use solver_class use chrolog_class use counter_class use timer_class implicit none #ifdef _CLOVER type(quark_clover), intent(in) :: this ! quark parameters #else type(quark_wilson), intent(in) :: 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) :: u ! gauge field type(chrolog_alg), intent(inout), optional :: cron type(_HMC_SOLVER_), pointer :: solver type(field_quark_eo_wg), allocatable :: we type(counter) :: imult integer :: istat, mode_solver, iter character(len=CHARLEN) :: str type(timer) :: solver_time real(DP) :: etime, rtmp0,rtmp1,res_norm #ifdef _GPU_SOLVER_ type(gpu_solver) :: gpu_precnd ! type(single_solver) :: gpu_precnd #endif call new(solver_time) call tic(solver_time) allocate(solver) allocate(we) call new(we,be%ieo) call new(xe,be%ieo) #ifdef _GPU_SOLVER_ mode_solver = MODE_PRECOND #else mode_solver = MODE_NORMAL #endif #ifdef _DEBUG rtmp0 = abs2(be) if (0 == nodeid) then write(*,'("%%INV_MULT_EOPREC_WD IEO=",I2)')be%ieo write(*,'("%% SOURCE |Be|^2 = ",E24.15)')rtmp0 endif #endif imult = mult_iter if (present(cron)) then if (get_ncron(cron) > 0) then if (nodeid == 0) then write(iout,'("# CRON USED")') endif endif call pack(be,cron%src_vec) do call solve(cron) select case (get_status(cron)) case (CHRON_OP_NOP) cycle case (CHRON_OP_DO_MATVEC) call unpack(cron%src_vec,xe) call assign_mult_eoprec_wd(this,we,xe,u) call pack(we,cron%dst_vec) case (CHRON_OP_CONVERGED) exit end select enddo call new(solver,NSIZE=NSITES,max_iter=2500,guess=GUESS_USE,mode=mode_solver,tol=err) ! call new(solver,NSIZE=NSITES,max_iter=2500,guess=GUESS_USE,mode=mode_solver,tol=err,debug=.true.) solver%src_vec(:) = cron%src_vec(:) solver%dst_vec(:) = cron%dst_vec(:) if (get_ncron(cron) > 0) then if (nodeid == 0) then write(iout,'("#",I5," ERR_=",E24.16," MULT_=",I10)') 0, get_residual_norm(cron),get_count(mult_iter)-get_count(imult) endif endif else call new(solver,NSIZE=NSITES,max_iter=2500,guess=GUESS_NO,mode=mode_solver,tol=err) ! call new(solver,NSIZE=NSITES,max_iter=2500,guess=GUESS_NO,mode=mode_solver,tol=err,debug=.true.) call pack(be,solver%src_vec) endif #ifdef _GPU_SOLVER_ #ifdef _CLOVER call new(gpu_precnd,u,get_ptr_inverse_clover_term(this)) #else call new(gpu_precnd,u) #endif #endif do call solve(solver) istat = get_status(solver) select case (istat) case (OP_NOP) cycle case (OP_DO_MATVEC) !==================== ! v ! we = Dee xe !==================== call unpack(solver%src_vec,xe) call assign_mult_eoprec_wd(this,we,xe,u) call pack(we,solver%dst_vec) #ifdef _GPU_SOLVER_ case (OP_PRECOND) !============================= ! Preconditioner ! v ! dst_vec <= Dee \ src_vec !============================= gpu_precnd%tol = err gpu_precnd%err_norm = get_residual_norm(solver) call assign_inv_mult_eoprec_wd(gpu_precnd,get_kappa(this),solver%dst_vec,solver%src_vec,xe%ieo) #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) write(iout,'("#",I5," ERR_=",E24.16," MULT_=",I10)') iiter,res_norm,get_count(mult_iter)-get_count(imult) write(* ,'("#",I5," ERR_=",E24.16," MULT_=",I10)') iiter,res_norm,get_count(mult_iter)-get_count(imult) endif #endif cycle 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)) write(iout,'("#",I5," ERR_=",E24.16," MULT_=",I10)',advance='no') iiter,res_norm,get_count(mult_iter)-get_count(imult) #ifdef _DEBUG write(* ,'("#",A," iteration converged.")')TRIM(get_name(solver)) write(* ,'("#",I5," ERR_=",E24.16," MULT_=",I10)',advance='no') iiter,res_norm,get_count(mult_iter)-get_count(imult) #endif endif exit case (OP_MAXITER_REACHED) write(*,*)"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) write(iout,'("#",A," iteration does not converge.")')TRIM(get_name(solver)) write(iout,'("#",I5," ERR_=",E24.16," MULT_=",I10)') iiter,res_norm,get_count(mult_iter)-get_count(imult) endif write(str,'("solver did not converge.: ",A,I4)')__FILE__,__LINE__ call error_stop(str) case default write(str,'("something is wrong in solver.: ",A,I4)')__FILE__,__LINE__ call error_stop(str) end select enddo #ifdef _GPU_SOLVER_ call delete(gpu_precnd) #endif if (present(cron)) then cron%store_vec(:) = solver%dst_vec(:) call store(cron) endif call unpack(solver%dst_vec,xe) call delete(solver) deallocate(solver) #ifdef _DEBUG rtmp0 = abs2(xe) if (0 == nodeid) then write(*,'("%% SOLUTION |Dee \ Be|^2 = ",E24.15)')rtmp0 endif ! ! compute true residual ! v ! we = Dee xe - be ! call assign_mult_eoprec_wd(this,we,xe,u) call accum_sub(we,be) rtmp0 = abs2(be) rtmp1 = abs2(we) rtmp0 = sqrt(rtmp1/rtmp0) if (0 == nodeid) then write(*,'("%",I5," ERR_=",E24.16," MULT_=",I10)') iiter,rtmp0,get_count(mult_iter)-get_count(imult) endif #endif deallocate(we) call toc(solver_time) etime = get_elapse(solver_time) if (nodeid == 0) then write(iout,'(" ETIME_=",E24.16)')etime endif return end subroutine
Subroutine : | |||
this : | type(quark_clover), intent(in)
| ||
yde : | type(field_quark_eo_wg), intent(inout)
| ||
ye : | type(field_quark_eo_wg), intent(inout)
| ||
u : | type(vfield_gluon_wg), intent(in)
|
Multiply even/odd-site and clover term preconditioned Wilson-Clover-Dirac operator(even/odd sites only):
v yde <= Dee ye, where v Dee = 1ee - kappa^2 Fee \ Meo Foo \ Moe
Original external subprogram is quark_clover_class#assign_mult_eoprec_wd
Subroutine : | |||
this : | type(quark_wilson), intent(in)
| ||
yde : | type(field_quark_eo_wg), intent(inout)
| ||
ye : | type(field_quark_eo_wg), intent(inout)
| ||
u : | type(vfield_gluon_wg), intent(in)
|
v yde = Dee ye, where v Dee = 1ee - kappa^2 Meo Moe
Original external subprogram is quark_wilson_class#assign_mult_eoprec_wd
Subroutine : | |
this : | type(hmc_std_quark_wilson), intent(inout) |
delete HMC quark even/odd preconditioned action algorithm
subroutine delete_wilson_quark_hmc(this) ! ! delete HMC quark even/odd preconditioned action algorithm ! implicit none type(hmc_std_quark_wilson), intent(inout) :: this deallocate(this%y) NULLIFY(this%y) call delete(this%phys_param) call delete(this%solver_log) return end subroutine
Subroutine : | |||
this : | type(hmc_std_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_hmc(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_std_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 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,fy,xo) !************************************************************************ !* 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) !************************************** !* v !* fye <= Meo fyo !* !************************************** call assign_mult_hop(this%phys_param,fy%eo(0),fy%eo(1),u) !************************************** !* xo <= (Qoo)^-1 gamma_5 Fo^-1 fyo !************************************** call assign_mult_gamma5(fx%eo(1),fy%eo(1)) #ifdef _CLOVER call accum_mult_clover_term(fx%eo(1),this%phys_param) #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) !************************************** !* !* fxo <= gamma_5 xo !* v !* fxe <= gamma_5 Meo xo !* !************************************** call assign_mult_gamma5(fx%eo(1),xo) call assign_mult_hop(this%phys_param,fx%eo(0),xo,u) call accum_mult_gamma5(fx%eo(0)) call copy_boundary(fy%eo(0)) call copy_boundary(fx%eo(0)) call copy_boundary(fx%eo(1)) !====================================== ! force from hopping matrix !====================================== kappa = get_kappa(this%phys_param) kappa2 = kappa**2 call force_hmc_hopping(BB,kappa2,fx,fy) #ifdef _CLOVER !************************************** !* v !* fyo <= Moe fye !* !************************************** call assign_mult_hop(this%phys_param,fy%eo(1),fy%eo(0),u) !====================================== ! force from clover term !====================================== csw = get_csw(this%phys_param) zickd8 = zi*csw*kappa**3/8.0_DP call force_hmc_clover(this%XX,zickd8,fx,fy) #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 interafaceing !================================== this%XX => NULL() #endif return end subroutine
Subroutine : | |||
this : | type(hmc_std_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 = |Doo^-1 yo| - 2 log[det[F]] (yo : lives only on odd field)
subroutine hamil_wilson_quark_hmc(this,ifi,u,action) ! ! Calculate Hamiltonian for HMC quark even/odd-site preoconditioned action. ! ! Action : ! v ! Shmc = |Doo^-1 yo| - 2 log[det[F]] ! ! (yo : lives only on odd field) ! use comlib implicit none type(hmc_std_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 character(CHARLEN) :: str real(DP) :: err integer :: iiter,iout 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 if (this%ncron > 0) then allocate(this%cron1) allocate(this%cron2) call new(this%cron1,NSITES,this%ncron) call new(this%cron2,NSITES,this%ncron) endif allocate(wo) call new(wo,1) !===================== ! Heatbath !===================== call set_gaussian_noise(wo) !===================== ! calc. action !===================== Spf=abs2(wo) !===================== ! v ! yo = Doo wo !===================== call assign_mult_eoprec_wd(this%phys_param,this%y%eo(1),wo,u) if (this%ncron > 0) then call pack(wo,this%cron1%store_vec) call store(this%cron1) endif 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) !============================ ! v ! solve : wo = Doo^-1 yo !============================ err = this%hamil_inv%tol call assign_inv_mult_eoprec_wd(this%phys_param,iout,err,iiter,wo,this%y%eo(1),u,this%cron1) call inc(this%hamil_inv,iiter) !=================================== ! action = |Doo^-1 yo|^2 = |wo|^2 !=================================== Spf=abs2(wo) deallocate(wo) if (this%ncron > 0) then call delete(this%cron1) call delete(this%cron2) deallocate(this%cron1) deallocate(this%cron2) endif #ifdef _REVCHECK if (this%ncron > 0) then allocate(this%cron1) allocate(this%cron2) call new(this%cron1,NSITES,this%ncron) call new(this%cron2,NSITES,this%ncron) endif #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_std_quark_wilson), intent(in) |
inquire whether this algorithm has GMP method
function has_gmp_wilson_quark_hmc(this) result(has_gmp) ! ! inquire whether this algorithm has GMP method ! implicit none type(hmc_std_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_std_quark_wilson), intent(in) |
inquire whether this algorithm has reweighting
function has_reweighting_wilson_quark_hmc(this) result(has_rew) ! ! inquire whether this algorithm has reweighting ! implicit none type(hmc_std_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_std_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_hmc(this,u1,u0,ihit) ! ! Nothing is done ! implicit none type(hmc_std_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_std_quark_wilson), intent(inout) |
id : | integer, intent(in) |
initialize/create HMC quark even/odd preconditioned action algorithm
subroutine new_wilson_quark_hmc(this,id) ! ! initialize/create HMC quark even/odd preconditioned action algorithm ! implicit none type(hmc_std_quark_wilson), intent(inout) :: this integer, intent(in) :: id character(len=CHARLEN) :: str NULLIFY(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_std_quark_wilson), intent(inout) |
print out HMC quark even/odd preconditioned algorithm parameters on disyplay.
subroutine print_wilson_quark_hmc(this) ! ! print out HMC quark even/odd preconditioned algorithm parameters on disyplay. ! implicit none type(hmc_std_quark_wilson), intent(inout) :: this if (nodeid==0) then write(*,'("**** Quark # ",I3," HMC")')this%myid write(*,'("***** Physical parameters *****")') call print(this%phys_param) write(*,'("***** HMC Algorithm parameters *****")') 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_std_quark_wilson), intent(inout) |
Print out HMC quark even/odd-site preconditioned algorithm job statistics on display
subroutine print_stat_wilson_quark_hmc(this) ! ! Print out HMC quark even/odd-site preconditioned algorithm job statistics on display ! implicit none type(hmc_std_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")')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_std_quark_wilson), intent(inout) |
iout : | integer, intent(in) |
Read HMC quark even/odd-site preoconditioned algorithm parameters
subroutine read_wilson_quark_hmc(this,iout) ! ! Read HMC quark even/odd-site preoconditioned algorithm parameters ! use comlib implicit none type(hmc_std_quark_wilson), intent(inout) :: this integer, intent(in) :: iout allocate(this%y) call new(this%y) 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 return end subroutine
Subroutine : | |||
this : | type(hmc_std_quark_wilson),
intent(inout)
| ||
u : | type(vfield_gluon_wg), intent(in) |
Nothing is done
subroutine reweighting_wilson_quark_hmc(this,u) ! ! Nothing is done ! implicit none type(hmc_std_quark_wilson), intent(inout) :: this ! quark action/algorihtm parameters type(vfield_gluon_wg), intent(in) :: u return end subroutine
Subroutine : | |||
this : | type(hmc_std_quark_wilson), intent(in) | ||
iout : | integer, intent(in)
|
Save HMC Quark even/odd-site preoconditioned algorihtm parameters on measurement configuration
subroutine save_config_hmc_std_quark_wilson(this,iout) ! ! Save HMC Quark even/odd-site preoconditioned algorihtm parameters on measurement configuration ! implicit none type(hmc_std_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_std_quark_wilson), intent(inout) |
id : | integer, intent(in) |
set quark # id
subroutine set_id_hsqw(this,id) ! ! set quark # id ! implicit none type(hmc_std_quark_wilson), intent(inout) :: this integer, intent(in) :: id this%myid = id call set_id(this%phys_param,id) return end subroutine