Class | hmc_mp_quark_wilson_class |
In: |
HmcMpQuarkWilsonClass/hmc_mp_quark_wilson_class.F90
|
$Id: hmc_mp_quark_wilson_class.F90,v 1.3 2011/06/17 04:48:58 ishikawa Exp $
Subroutine : | |
this : | type(hmc_mp_quark_wilson), intent(inout) |
subroutine delete_mp_wq_hmc(this) implicit none type(hmc_mp_quark_wilson), intent(inout) :: this deallocate(this%y) NULLIFY(this%y) call delete(this%phys_param) call delete(this%solver_uv_log) call delete(this%solver_ir_log) return end subroutine
Subroutine : | |||
this : | type(hmc_mp_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 mass preconditioned pseudo fermion
Action : S_Q = S_Q(UV) + S_Q(IR) v S_Q(UV) : |(D(UV)oo)^-1 y(UV)o|^2 (odd site only) v v S_Q(IR) : |(D(IR)oo)^-1 D(UV)oo y(IR)o|^2 (odd site only) v v where Doo = 1 - Moo v Moo = kappa^2 Fo^-1 Moe Fe^-1 Meo
subroutine force_mp_wq_hmc(this,BB,idepth,u) ! ! Calc force from hopping part of mass preconditioned pseudo fermion ! ! Action : S_Q = S_Q(UV) + S_Q(IR) ! v ! S_Q(UV) : |(D(UV)oo)^-1 y(UV)o|^2 (odd site only) ! v v ! S_Q(IR) : |(D(IR)oo)^-1 D(UV)oo y(IR)o|^2 (odd site only) ! ! v v ! where Doo = 1 - Moo ! v ! Moo = kappa^2 Fo^-1 Moe Fe^-1 Meo ! ! implicit none type(hmc_mp_quark_wilson), intent(inout) :: this ! quark action parameters type(vfield_gluon_wg), intent(inout) :: BB ! pre-force contribution (\dot{u}) integer, intent(in) :: idepth ! MD depth for recursive MD integrator type(vfield_gluon_wg), intent(in) :: u ! gauge field #ifdef _CLOVER real(DP) :: dummy #endif #ifdef _CLOVER if (idepth == this%depth_hopping_uv .or. idepth == this%depth_hopping_ir .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_uv) then !================================================= ! UV part force !================================================= #ifdef _MD_DEPTH_ if (nodeid==0) write(*,'("depth=",I3," : force quark mprecd hmc hopping(UV)")')idepth #endif call force_uv_mp_wq_hmc(this,BB,u) endif if (idepth == this%depth_hopping_ir) then !================================================= ! IR part force !================================================= #ifdef _MD_DEPTH_ if (nodeid==0) write(*,'("depth=",I3," : force quark mprecd hmc hopping(IR)")')idepth #endif call force_ir_mp_wq_hmc(this,BB,u) endif #ifdef _CLOVER if (idepth == this%depth_clvtrlog) then !====================================== ! calc force from -Nf*log(det[F]) !====================================== #ifdef _MD_DEPTH_ if (nodeid==0) write(*,'("depth=",I3," : force quark std hmc trlog")')idepth #endif call force_clover_trlog(this%phys_param,this%XX) endif !================================= ! remove link to XX ! use XX for interfacing only !================================= this%XX => NULL() call delete_inverse_clover_term(this%phys_param) #endif return end subroutine
Subroutine : | |||
this : | type(hmc_mp_quark_wilson),
intent(inout)
| ||
ifi : | integer, intent(in)
| ||
u : | type(vfield_gluon_wg), intent(in)
| ||
action : | real(DP), intent(out)
|
Calculate Hamiltonian of Mass preconditioned action
Action : v v v Shmc = |D(IR)oo^-1 D(UV)oo y(IR)o|^2 + |D(UV)oo^-1 y(UV)o| - 2 log[det[F]] (yo : lives only on odd field)
subroutine hamil_mp_wq_hmc(this,ifi,u,action) ! ! Calculate Hamiltonian of Mass preconditioned action ! ! Action : ! v v v ! Shmc = |D(IR)oo^-1 D(UV)oo y(IR)o|^2 + |D(UV)oo^-1 y(UV)o| - 2 log[det[F]] ! ! (yo : lives only on odd field) ! ! use comlib use hmc_status_class implicit none type(hmc_mp_quark_wilson), intent(inout) :: this ! quark action parameters integer, intent(in) :: ifi ! targeg swithch, 0 for initial, 1 for final type(vfield_gluon_wg), intent(in) :: u ! gauge field real(DP), intent(out) :: action ! action real(DP) :: Spf(NUM_PF),Sclv type(field_quark_eo_wg), allocatable :: wo,ro character(CHARLEN) :: str real(DP) :: err, kappa_uv integer :: iiter,iout type(hmc_status) :: status #ifdef _CLOVER type(quark_clover), allocatable :: uv_quark_param #else type(quark_wilson), allocatable :: uv_quark_param #endif 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 !======================= ! Kappa_uv = Kappa*rho !======================= allocate(uv_quark_param) uv_quark_param = this%phys_param kappa_uv = get_kappa(this%phys_param) * this%rho call set_kappa(uv_quark_param,kappa_uv) select case(ifi) case(0) ! INITIAL if (this%ncron > 0) then allocate(this%cron_uv1) allocate(this%cron_uv2) allocate(this%cron_ir1) allocate(this%cron_ir2) call new(this%cron_uv1,NSIZE=COL*SPIN*NTH*NZ*NY*NX,NMAX=this%ncron) call new(this%cron_uv2,NSIZE=COL*SPIN*NTH*NZ*NY*NX,NMAX=this%ncron) call new(this%cron_ir1,NSIZE=COL*SPIN*NTH*NZ*NY*NX,NMAX=this%ncron) call new(this%cron_ir2,NSIZE=COL*SPIN*NTH*NZ*NY*NX,NMAX=this%ncron) endif !=================================================== ! UV part !=================================================== allocate(wo) call new(wo,1) call set_gaussian_noise(wo) Spf(UVPF) = abs2(wo) !======================== ! v ! y(UV)o = D(UV)oo wo !======================== call assign_mult_eoprec_wd(uv_quark_param,this%y(UVPF)%eo(1),wo,u) if (this%ncron > 0) then call pack(wo,this%cron_uv1%store_vec) call store(this%cron_uv1) endif !=================================================== ! IR part !=================================================== call set_gaussian_noise(wo) Spf(IRPF) = abs2(wo) !=============================== ! v v ! yo = D(UV)oo^-1 D(IR)oo wo !=============================== allocate(ro) call assign_mult_eoprec_wd(this%phys_param,ro,wo,u) err = this%hamil_uv_inv%tol iout = get_file_unit(this%solver_uv_log) call assign_inv_mult_eoprec_wd(uv_quark_param,iout,err,iiter,this%y(IRPF)%eo(1),ro,u) call inc(this%hamil_uv_inv,iiter) deallocate(wo,ro) case(1:2) ! FINAL:REVCHECK !=================================================== ! UV part !=================================================== !===================================== ! v ! solve : wo = D(UV)oo^-1 y(UV)o !===================================== allocate(wo) err = this%hamil_uv_inv%tol iout = get_file_unit(this%solver_uv_log) call assign_inv_mult_eoprec_wd(uv_quark_param,iout,err,iiter,wo,this%y(UVPF)%eo(1),u,this%cron_uv1) call inc(this%hamil_uv_inv,iiter) !====================================== ! v ! action = |D(UV)oo^-1 y(UV)o|^2 ! = |wo|^2 !====================================== Spf(UVPF) = abs2(wo) !=================================================== ! IR part !=================================================== !========================================== ! v ! solve : wo = D(IR)oo^-1 D(UV)oo y(IR)o !========================================== allocate(ro) call assign_mult_eoprec_wd(uv_quark_param,ro,this%y(IRPF)%eo(1),u) err = this%hamil_ir_inv%tol iout = get_file_unit(this%solver_ir_log) call assign_inv_mult_eoprec_wd(this%phys_param,iout,err,iiter,wo,ro,u) call inc(this%hamil_ir_inv,iiter) !====================================== ! v v ! action = |D(IR)oo^-1 D(UV)oo y(IR)o|^2 ! = |wo|^2 !====================================== Spf(IRPF) = abs2(wo) deallocate(wo,ro) if (this%ncron > 0) then call delete(this%cron_uv1) call delete(this%cron_uv2) call delete(this%cron_ir1) call delete(this%cron_ir2) deallocate(this%cron_uv1) deallocate(this%cron_uv2) deallocate(this%cron_ir1) deallocate(this%cron_ir2) endif end select deallocate(uv_quark_param) if (0 == this%depth_hopping_uv) Spf(UVPF) = 0.0_DP if (0 == this%depth_hopping_ir) Spf(IRPF) = 0.0_DP #ifdef _CLOVER call delete_inverse_clover_term(this%phys_param) #endif action=Spf(IRPF)+Spf(UVPF)+Sclv write(str,'("@",I8,I2," Quark: ID:",I3," SQPF(IR):",E24.16," SQPF(UV):",E24.16," SCLV:",E24.16)') get_trajectory_number(status),ifi,this%myid,Spf(IRPF),Spf(UVPF),Sclv call print_log_action(status,str) return end subroutine
Function : | |
has_gmp : | logical |
this : | type(hmc_mp_quark_wilson), intent(in) |
function has_gmp_mp_wq_hmc(this) result(has_gmp) implicit none type(hmc_mp_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_mp_quark_wilson), intent(in) |
function has_reweighting_mp_wq_hmc(this) result(has_rew) implicit none type(hmc_mp_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 |
Mass + even/odd site preconditioned quark action
Subroutine : | |||
this : | type(hmc_mp_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_mp_wq_hmc(this,u1,u0,ihit) ! ! Nothing is done ! implicit none type(hmc_mp_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_mp_quark_wilson), intent(inout) |
id : | integer, intent(in) |
subroutine new_mp_wq_hmc(this,id) implicit none type(hmc_mp_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_uv_inv) call new(this%force_ir_inv) call new(this%hamil_uv_inv) call new(this%hamil_ir_inv) this%depth_hopping_uv = 1 this%depth_hopping_ir = 1 this%depth_clvtrlog = 1 write(str,'(I2.2)')id this%solver_uv_log_fname = ADJUSTL(TRIM(this%solver_uv_log_fname)//TRIM(str)) this%solver_ir_log_fname = ADJUSTL(TRIM(this%solver_ir_log_fname)//TRIM(str)) call new(this%solver_uv_log,TRIM(this%solver_uv_log_fname)) call new(this%solver_ir_log,TRIM(this%solver_ir_log_fname)) return end subroutine
Subroutine : | |
this : | type(hmc_mp_quark_wilson), intent(inout) |
Print out Mass Preconditioned quark action/algorithm parameters
subroutine print_mp_wq_hmc(this) ! ! Print out Mass Preconditioned quark action/algorithm parameters ! implicit none type(hmc_mp_quark_wilson), intent(inout) :: this character(CHARLEN) :: str if (nodeid==0) then write(*,'("**** Quark # ",I3," MassPreconditioned HMC")')this%myid write(*,'("***** Physical parameters *****")') call print(this%phys_param) write(*,'("***** HMC Algorithm parameters *****")') write(*,'(" Mass ratio rho kappa_uv/kappa :",E16.8)')this%rho write(*,'(" UV Force Stopping Condition :",E16.8)')this%force_uv_inv%tol write(*,'(" IR Force Stopping Condition :",E16.8)')this%force_ir_inv%tol write(*,'(" UV Hamil Stopping Condition :",E16.8)')this%hamil_uv_inv%tol write(*,'(" IR Hamil Stopping Condition :",E16.8)')this%hamil_ir_inv%tol write(*,'(" force depth UV PSfermion part :",I3)')this%depth_hopping_uv write(*,'(" force depth IR PSfermion part :",I3)')this%depth_hopping_ir #ifdef _CLOVER write(*,'(" force depth CLV trlog part :",I3)')this%depth_clvtrlog #endif write(*,'(" Chronological History :",I3)')this%ncron write(*,'(80("-"))') endif if (this%rho > 1.0_DP) then write(str,'("Mass ratio should be less than one!, rho=",E24.15)')this%rho call error_stop(str) endif return end subroutine
Subroutine : | |
this : | type(hmc_mp_quark_wilson), intent(inout) |
Print out Mass preconditioned quark action/algorithm statistics on display
subroutine print_stat_mp_wq_hmc(this) ! ! Print out Mass preconditioned quark action/algorithm statistics on display ! implicit none type(hmc_mp_quark_wilson), intent(inout) :: this real(DP) :: f_uv_ave,f_ir_ave real(DP) :: h_uv_ave,h_ir_ave f_uv_ave = get_average(this%force_uv_inv) f_ir_ave = get_average(this%force_ir_inv) h_uv_ave = get_average(this%hamil_uv_inv) h_ir_ave = get_average(this%hamil_ir_inv) if (nodeid==0) then write(*,'("**** Quark # ",I3," MassPreconditioned HMC")')this%myid call print(this%phys_param) write(*,'(" MPHMC algorithm ")') write(*,'(" Mass ratio rho kappa_uv/kappa :",E16.8)')this%rho write(*,'(" UV Force BiCGStab iter")') write(*,'(" (averaged) :",F12.4)') f_uv_ave write(*,'(" IR Force BiCGStab iter")') write(*,'(" (averaged) :",F12.4)') f_ir_ave write(*,'(" UV Hamiltonian BiCGStab iter")') write(*,'(" (averaged) :",F12.4)') h_uv_ave write(*,'(" IR Hamiltonian BiCGStab iter")') write(*,'(" (averaged) :",F12.4)') h_ir_ave write(*,'(" Chronological History :",I3)')this%ncron endif return end subroutine
Subroutine : | |
this : | type(hmc_mp_quark_wilson), intent(inout) |
iout : | integer, intent(in) |
Read Mass preconditioned HMC quark algorithm parameters and initialize parameters
subroutine read_mp_wq_hmc(this,iout) ! ! Read Mass preconditioned HMC quark algorithm parameters ! and initialize parameters ! use comlib implicit none type(hmc_mp_quark_wilson), intent(inout) :: this integer, intent(in) :: iout integer :: ipf allocate(this%y(NUM_PF)) do ipf=1,NUM_PF call new(this%y(ipf)) enddo call new(this%hamil_uv_inv) call new(this%hamil_ir_inv) call new(this%force_uv_inv) call new(this%force_ir_inv) call read(this%phys_param,iout) if (nodeid==0) then read(iout,*) read(iout,*)this%rho read(iout,*)this%hamil_uv_inv%tol read(iout,*)this%force_uv_inv%tol read(iout,*)this%depth_hopping_uv read(iout,*)this%depth_hopping_ir #ifdef _CLOVER read(iout,*)this%depth_clvtrlog #endif read(iout,*)this%ncron endif #ifndef _singlePU call comlib_bcast(this%rho,0) call comlib_bcast(this%hamil_uv_inv%tol,0) call comlib_bcast(this%force_uv_inv%tol,0) call comlib_bcast(this%depth_hopping_uv,0) call comlib_bcast(this%depth_hopping_ir,0) #ifdef _CLOVER call comlib_bcast(this%depth_clvtrlog,0) #endif call comlib_bcast(this%ncron,0) #endif this%hamil_ir_inv%tol = this%hamil_uv_inv%tol this%force_ir_inv%tol = this%force_uv_inv%tol return end subroutine
Subroutine : | |||
this : | type(hmc_mp_quark_wilson),
intent(inout)
| ||
u : | type(vfield_gluon_wg), intent(in) |
Nothing is done
subroutine reweighting_mp_wq_hmc(this,u) ! ! Nothing is done ! implicit none type(hmc_mp_quark_wilson), intent(inout) :: this ! quark action/algorihtm parameters type(vfield_gluon_wg), intent(in) :: u return end subroutine
Subroutine : | |
this : | type(hmc_mp_quark_wilson), intent(in) |
iout : | integer, intent(in) |
Save Quark parameters on configuration
subroutine save_config_mp_wq_hmc(this,iout) ! ! Save Quark parameters on configuration ! implicit none type(hmc_mp_quark_wilson), intent(in) :: this integer, intent(in) :: iout call save_config(this%phys_param,iout) return end subroutine
Subroutine : | |
this : | type(hmc_mp_quark_wilson), intent(inout) |
id : | integer, intent(in) |
subroutine set_id_mp_wq_hmc(this,id) implicit none type(hmc_mp_quark_wilson), intent(inout) :: this integer, intent(in) :: id this%myid = id call set_id(this%phys_param,id) return end subroutine