Class hmc_std_quark_wilson_class
In: HmcStdQuarkWilsonClass/hmc_std_quark_wilson_class.F90
lattice_class error_class hmc_logfile_class solver_parameter_class chrolog_class hmc_status_class field_gauge_class field_fermion_class quark_clover_class quark_wilson_class comlib gpu_solver_class solver_class counter_class timer_class hmc_std_quark_wilson_class dot/f_265.png

Definition of wilson/clover quark Standard(even/odd preconditioned) HMC algorithm

Version

 $Id: hmc_std_quark_wilson_class.F90,v 1.32 2011/09/28 07:16:16 ishikawa Exp $

Methods

Included Modules

lattice_class error_class hmc_logfile_class solver_parameter_class chrolog_class hmc_status_class field_gauge_class field_fermion_class quark_clover_class quark_wilson_class comlib gpu_solver_class solver_class counter_class timer_class

Public Instance methods

HMC_ALGORITHM
Constant :
HMC_ALGORITHM = 1 :integer, parameter
Subroutine :
this :type(quark_clover), intent(in)
: quark parameters
this :type(quark_wilson), intent(in)
: quark parameters
iout :integer, intent(in)
: solver logfile unit number
err :real(DP), intent(in)
: solver stopping condition
iiter :integer, intent(out)
: total iteration count after convergence
xe :type(field_quark_eo_wg), intent(inout)
: even/odd site solution field
be :type(field_quark_eo_wg), intent(in)
: even/odd site source(rhs) field
u :type(vfield_gluon_wg), intent(in)
: gauge field
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

[Source]

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
assign_mult_eoprec_wd( this, yde, ye, u )
Subroutine :
this :type(quark_clover), intent(in)
: clover quark parameters and inverse clover term
yde :type(field_quark_eo_wg), intent(inout)
: even/odd site fermion field
ye :type(field_quark_eo_wg), intent(inout)
: even/odd site fermion field
u :type(vfield_gluon_wg), intent(in)
: gauge field

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

assign_mult_eoprec_wd( this, yde, ye, u )
Subroutine :
this :type(quark_wilson), intent(in)
: wilson quark parameter
yde :type(field_quark_eo_wg), intent(inout)
: even/odd fermion field
ye :type(field_quark_eo_wg), intent(inout)
: even/odd fermion field / boundary copied
u :type(vfield_gluon_wg), intent(in)
: gauge field

Multiply even/odd-site preconditioned Dirac operator :

       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

[Source]

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)
: quark action/algorithm parameters
BB :type(vfield_gluon_wg), intent(inout)
: pre-force contribution from hopping part
idepth :integer, intent(in)
: MD depth for recursive MD integrator
u :type(vfield_gluon_wg), intent(in)
: gauge field

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

[Source]

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)
: quark action/algorithm parameters
ifi :integer, intent(in)
: target switch : 0 for initial generating PF field, 1 for final
u :type(vfield_gluon_wg), intent(in)
: gauge field
action :real(DP), intent(out)
: output action value

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)

[Source]

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

[Source]

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

[Source]

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
hmc_std_quark_wilson
Derived Type :
phys_param :type(quark_clover), public
phys_param :type(quark_wilson), public
XX => NULL() :type(tfield_gluon_wg), public, pointer
: pre-force contribution from clvoer term
phys_param :type(quark_clover), public
phys_param :type(quark_wilson), public

ifdef _CLOVER

Subroutine :
this :type(hmc_std_quark_wilson), intent(inout)
: quark action/algorihtm parameters
u1 :type(vfield_gluon_wg), intent(in)
: MD final gauge field
u0 :type(vfield_gluon_wg), intent(in)
: MD inital gauge field
ihit :integer, intent(inout)
: 1 for accept, 0 for reject

Nothing is done

[Source]

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

[Source]

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.

[Source]

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

[Source]

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

[Source]

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)
: quark action/algorihtm parameters
u :type(vfield_gluon_wg), intent(in)

Nothing is done

[Source]

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)
: configuration file unit id

Save HMC Quark even/odd-site preoconditioned algorihtm parameters on measurement configuration

[Source]

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

[Source]

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