Class hmc_poly_quark_wilson_class
In: HmcPolyQuarkWilsonClass/hmc_poly_quark_wilson_class.F90
counter_class hmc_logfile_class error_class solver_parameter_class lattice_class metropolis_test_class hmc_status_class field_gauge_class field_fermion_class quark_clover_class quark_wilson_class comlib random_class timer_class solver_class hmc_poly_quark_wilson_class dot/f_134.png

Define HMC quark Polynomial HMC algorithm (even/odd-site preconditioned)

Version

$Id: hmc_poly_quark_wilson_class.F90,v 1.26 2011/06/21 03:29:27 ishikawa Exp $

Methods

Included Modules

counter_class hmc_logfile_class error_class solver_parameter_class lattice_class metropolis_test_class hmc_status_class field_gauge_class field_fermion_class quark_clover_class quark_wilson_class comlib random_class timer_class solver_class

Public Instance methods

PHMC_ALGORITHM
Constant :
PHMC_ALGORITHM =0 :integer, parameter
Subroutine :
this :type(hmc_poly_quark_wilson), intent(inout)

[Source]

subroutine delete_wilson_quark_phmc(this)
  implicit none
  type(hmc_poly_quark_wilson), intent(inout) :: this
  deallocate(this%y)
  call delete(this%phys_param)
  call delete(this%solver_log)
  call delete(this%gmp_status)
  NULLIFY(this%y)
  return
end subroutine
Subroutine :
this :type(hmc_poly_quark_wilson), intent(inout)
: quark action/algorithm parameters
BB :type(vfield_gluon_wg), intent(inout)
: pre-force contribution (dot{u})
idepth :integer, intent(in)
: MD depth for recuresive MD integrator
u :type(vfield_gluon_wg), intent(in)
: gauge field

Calc MD force from hoping part of the even/odd-site preconditioned Polynomial HMC action.

  Action : Sphmc
                           v
    Sphm : PHMC action, |P[Doo]yo|^2  (odd site only)
                 v
           where Doo = 1oo - kappa^2 Fo^-1 Moe Fe^-1 Meo

[Source]

subroutine force_wilson_quark_phmc(this,BB,idepth,u)
!
! Calc MD force from hoping part of the even/odd-site preconditioned Polynomial HMC action.
!
!   Action : Sphmc
!                            v
!     Sphm : PHMC action, |P[Doo]yo|^2  (odd site only)
!                  v
!            where Doo = 1oo - kappa^2 Fo^-1 Moe Fe^-1 Meo
! 
 
  implicit none
  type(hmc_poly_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
  type(field_quark_wg),    allocatable :: fy(:),fx,fz
  type(field_quark_eo_wg), allocatable :: xo,wo
#ifdef _CLOVER
  real(DP) :: dummy
  real(DP) :: csw
  complex(DP) :: zickd8
#endif
  integer :: imb,nmb
  complex(DP) :: ccoef
  complex(DP), allocatable :: zcoef(:)
  integer :: nf
  real(DP) :: kappa


#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

  nf    = get_nflavor(this%phys_param)
  kappa =   get_kappa(this%phys_param)

  nmb=this%npoly/(3-nf)

  allocate(zcoef(0:nmb))
  allocate(fy(0:nmb),fx,fz,xo,wo)
  do imb=0,nmb
    call new(fy(imb)%eo(0),0)
    call new(fy(imb)%eo(1),1)
  enddo

  kappa2 = kappa**2
  fcoef  = -kappa2
  do imb=1,nmb
    zcoef(imb)=-(this%zcoef(imb)/this%zcoef(imb-1))*kappa2
  enddo
#ifdef _CLOVER
    csw  =     get_csw(this%phys_param)
  zickd8 = -(zi*csw*kappa**3)/8.0_DP
#endif

!======================================================================
! Prepare gauge force contributions from polynomial pseudo ferimon
!======================================================================
  !===========================
  !  fyo(nmb) <= yo
  !===========================
  call assign(fy(nmb)%eo(1),this%y%eo(1))

  do imb=nmb,1,-1
    !*************************************************************
    !*                              v
    !*  fyo^(j-1) = fyo^(N) + zcoef Moo fyo^(j), for j=N,N-1,...,1
    !*              _
    !*  fye^(j-1) = Meo fyo^(j-1)
    !*
    !*************************************************************
    call assign_mult_hop(this%phys_param,fy(imb  )%eo(0),fy(imb)%eo(1),u)
    call assign_mult_hop(this%phys_param,fy(imb-1)%eo(1),fy(imb)%eo(0),u)
    !===========================================
    !  fyo(j-1) <= fyo(j-1)*zeof(j) + fyo(N)
    !===========================================
    call accum_mult_add(fy(imb-1)%eo(1),zcoef(imb),fy(nmb)%eo(1))
  enddo ! end of do imb

!***************************************
!* xo := F^-1oo gamma_5 fyo^(0)
!***************************************
#ifdef _CLOVER
  call assign_mult_gamma5(xo,fy(0)%eo(1))
  call accum_mult_clover_term(xo,this%phys_param)
#else
  call assign_mult_gamma5(xo,fy(0)%eo(1))
#endif

!************************************************************************
!* Force loop
!************************************************************************
  do imb=1,nmb ! polynomial loop
!*************************************************************
!*
!* xo =  (c_1/c_0)^* xo    (for imb=1)
!*                           v
!* xo =  (c_imb/c_{imb-1})^* Moo xo    (for imb=2,3,...,N)
!*
!* fxo = gamma_5 xo
!*               _
!* fxe = gamma_5 Meo xo
!*
!*************************************************************
    select case (imb)
    case (1)
      ccoef=conjg(this%zcoef(1)/this%zcoef(0))
    case (2:)
      ccoef=conjg(zcoef(imb))
    end select
    !====================
    !  xo <= xo * ccoef
    !====================
    call accum_mult(xo,ccoef)

    call assign_mult_hop(this%phys_param,this%y%eo(0),xo,u)

    call assign_mult_gamma5(fx%eo(1),xo)
    call assign_mult_gamma5(fx%eo(0),this%y%eo(0))

    call copy_boundary(fx%eo(1))
    call copy_boundary(fx%eo(0))

    if (imb < nmb) then
      call assign_mult_hop(this%phys_param,xo,this%y%eo(0),u)
    endif

!======================================
! force from hopping matrix
!======================================
    call force_hmc_hopping(BB,fcoef,fx,fy(imb))

#ifdef _CLOVER
    call assign(fz%eo(0),fy(imb)%eo(0))
    call assign_mult_hop(this%phys_param,fz%eo(1),fz%eo(0),u)
!======================================
! force from clover term in hopping
!======================================
    call force_hmc_clover(this%XX,zickd8,fx,fz)
#endif

  enddo ! end of do imb

  deallocate(zcoef)
  deallocate(fx,fy,fz,xo,wo)
  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 :
this :type(hmc_poly_quark_wilson), intent(inout)
: quark action/algorithm paramters
ifi :integer, intent(in)
: target switch,
  • 0 for initial with generating PF field,
  • 1 for final action,
  • 2 for revcheck
u :type(vfield_gluon_wg), intent(in)
: gauge field
action :real(DP), intent(out)
: action value

Calculate Hamiltonian for PHMC action

Action :

            v
 Sphmc = |P[Doo] yo|^2 - 2 log[det[F]] for Nf=2, (yo:lives on odd site only)

   or
            v
       = |T[Doo] yo|^2 - 1 log[det[F]] for Nf=1, (yo:lives on odd site only)

  where

  P[x] = (1/c_0)Sum[c_i x^i,i=0,npoly],

  and

  T[x] = (1/d_0)Sum[d_i x^i,i=0,npoly/2].

[Source]

subroutine hamil_wilson_quark_phmc(this,ifi,u,action)
!
! Calculate Hamiltonian for PHMC action
!
! Action :
!             v 
!  Sphmc = |P[Doo] yo|^2 - 2 log[det[F]] for Nf=2, (yo:lives on odd site only)
!
!    or
!             v
!        = |T[Doo] yo|^2 - 1 log[det[F]] for Nf=1, (yo:lives on odd site only)
!
!   where 
! 
!   P[x] = (1/c_0)Sum[c_i x^i,i=0,npoly],
! 
!   and
! 
!   T[x] = (1/d_0)Sum[d_i x^i,i=0,npoly/2].
! 
  use comlib
  implicit none
  type(hmc_poly_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 :: wo,ro
  integer :: imb,nmb,iiter
  integer :: nf
  real(DP) :: Spf,Sclv
  real(DP) :: err,  kappa
  complex(DP), allocatable :: zcoef(:)
  character(CHARLEN) :: str
  type(hmc_status) :: status

  nf    = get_nflavor(this%phys_param)
  kappa =   get_kappa(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


  select case(ifi)
  case(0) ! INITIAL

    write(str,'("# HAMIL SOLVER")')
    call print(this%solver_log,TRIM(str))

    nmb = this%npoly/(3-nf)
    allocate(zcoef(0:nmb))
    allocate(wo,ro)
    call new(wo,1)
    call new(ro,1)
    call set_gaussian_noise(wo)
    !=============================
    ! calc. pseudo fermion action
    !=============================
    Spf=abs2(wo)

    !=============================
    !            v
    !  Solve : P[Doo] yo = wo
    ! 
    ! Using preconditioning. Solve
    !            v      
    !          W[Doo] ro = wo
    !            v
    !   =>  yo = Doo ro
    !=============================
    err = this%update_q_inv%tol
    call assign_inv_mult_poly_prec(this,err,iiter,ro,wo,u)
    call inc(this%update_q_inv,iiter)

    call assign_mult_eoprec_wd(this%phys_param,this%y%eo(1),ro,u)

    !=================================================
    ! 
    !  For Nf=1,
    ! 
    !   yo <= T[D]^* yo
    !                                            v                            v
    !      = (1/d^*_0) [ d^*_0 + d^*_1 (-kappa^2 Moo) + ... + d^*_n (-kappa^2 Moo)^n ] yo
    !                                  v                             v
    !      = [ 1 - kappa^2 (d_1/d_0)^* Moo [ 1 - kappa^2 (d_2/d_1)^* Moo [ ....
    !                   v                                 v
    !             ..... Moo [ 1 - kappa^2 (d_n/d_{n-1})^* Moo ]...] yo
    !=================================================
    if (1 == nf) then

      do imb=1,this%npoly/2
        zcoef(imb)=-conjg((this%zcoef(imb)/this%zcoef(imb-1)))*(kappa**2)
      enddo

      call assign(wo,this%y%eo(1))

      do imb=this%npoly/2,1,-1
        call assign_mult_hop(this%phys_param,this%y%eo(0),          wo,u)
        call assign_mult_hop(this%phys_param,          wo,this%y%eo(0),u)
        call accum_mult_add(wo,zcoef(imb),this%y%eo(1))
      enddo ! end of do imb

      call assign(this%y%eo(1),wo)
    endif

    deallocate(zcoef,wo,ro)

  case(1:2) ! FINAL : REVCHECK 
    nmb=this%npoly/(3-nf)
    allocate(zcoef(0:nmb))
    allocate(wo)
    call new(wo,1)

    !====================
    !        v
    ! wo = P[Doo] yo
    !
    ! or
    !        v
    ! wo = T[Doo] yo
    !====================
    do imb=1,nmb
      zcoef(imb)=-(this%zcoef(imb)/this%zcoef(imb-1))*(kappa**2)
    enddo

    call assign(wo,this%y%eo(1))

    do imb=nmb,1,-1
      call assign_mult_hop(this%phys_param,this%y%eo(0),          wo,u)
      call assign_mult_hop(this%phys_param,          wo,this%y%eo(0),u)
      call accum_mult_add(wo,zcoef(imb),this%y%eo(1))
    enddo

    Spf=abs2(wo)
    deallocate(zcoef,wo)
  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_poly_quark_wilson), intent(in)

inquire whether this algorithm has GMP

[Source]

function has_gmp_wilson_quark_phmc(this) result(has_gmp)
!
! inquire whether this algorithm has GMP
!
  implicit none
  type(hmc_poly_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_poly_quark_wilson), intent(in)

inquire whether this algorithm has reweighting

[Source]

function has_reweighting_wilson_quark_phmc(this) result(has_rew)
!
! inquire whether this algorithm has reweighting
!
  implicit none
  type(hmc_poly_quark_wilson), intent(in) :: this
  logical :: has_rew
  has_rew = this%has_reweighting
  return
end function
hmc_poly_quark_wilson
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

HMC quark even/odd sitem preconditioned Polynomical HMC algorithm action

Subroutine :
this :type(hmc_poly_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

Do global Metropolis test for PHMC action (correction factor)

[Source]

subroutine metropolis_wilson_quark_phmc(this,u1,u0,ihit)
!
! Do global Metropolis test for PHMC action (correction factor)
!
  implicit none
  type(hmc_poly_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
  integer :: nf

  nf = get_nflavor(this%phys_param)

  if (this%switch_metropolis_test == SW_OFF) return
  if (ihit == 0) return
  select case (nf)
  case (1)
    call metropolis_nf1_phmc(this,u1,u0,ihit)
  case (2)
    call metropolis_nf2_phmc(this,u1,u0,ihit)
  case default
    call error_stop(" # of flavor exceeds 2.")
  end select

  return
end subroutine
Subroutine :
this :type(hmc_poly_quark_wilson), intent(inout)
id :integer, intent(in)

[Source]

subroutine new_wilson_quark_phmc(this,id)
  implicit none
  type(hmc_poly_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%update_q_inv)
  call new(this%gmp_inv)

  write(str,'(I2.2)')id
  this%solver_log_fname = TRIM(this%solver_log_fname)//TRIM(str)
  this%gmp_log_fname    = TRIM(this%gmp_log_fname)//TRIM(str)

  call new(this%solver_log,this%solver_log_fname)
  call new(this%gmp_status,this%gmp_log_fname,this%gmp_test_name,this%gmp_weight_name)

  return
end subroutine
Subroutine :
this :type(hmc_poly_quark_wilson), intent(inout)

print out quark PHMC algorithm parameters on display

[Source]

subroutine print_wilson_quark_phmc(this)
!
! print out quark PHMC algorithm parameters on display
!
  implicit none
  type(hmc_poly_quark_wilson), intent(inout) :: this
  integer :: imb

  if (nodeid==0) then

  write(*,'("**** Quark # ",I3," PHMC")')this%myid

  write(*,'("***** Physical parameters *****")')
  call print(this%phys_param)
  write(*,'("***** PHMC Algorithm parameters *****")')
  write(*,'("    Heat-Bath Stopping Condition:",E16.8," (pseudo-fermion)")')this%update_q_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

  select case (this%switch_metropolis_test)
  case (0)
    write(*,'("              Global Metropolis : No ")')
  case (1)
    write(*,'("              Global Metropolis : Yes")')
    write(*,'("  Metropolis Stopping Condition :",E16.8)') this%gmp_inv%tol
  case default
    call error_stop("Metropolis test switch error")
  end select

  write(*,'(14X," Polynomial Order :",I4)') this%npoly

  select case (this%ptype)
  case (0)
    write(*,'(16X,"Polynomial type : Chebyshev")')
  case (1)
    write(*,'(16X,"Polynomial type : Adopted")')
  case default
    call error_stop("Polynomial type error")
  end select

  write(*,'("        Polynomial Coefficients :")')
  write(*,'(4X,7X,"c_i (real)")')
  write(*,'(80("-"))')
  do imb=0,this%npoly
    write(*,'(I4,E24.16)') imb,this%pcoef(imb)
  enddo
  write(*,'("  Force Polynomial Coefficients :")')
  write(*,'(4X,7X,"d_i (complex)")')
  write(*,'(80("-"))')
  do imb=0,this%npoly/(3-this%nflavor)
    write(*,'(I4,2E24.16)') imb,this%zcoef(imb)
  enddo
  write(*,'(80("-"))')

  endif

  return
end subroutine
Subroutine :
this :type(hmc_poly_quark_wilson), intent(inout)

print out quark PHMC algorithm job statistics on display

[Source]

subroutine print_stat_wilson_quark_phmc(this)
!
! print out quark PHMC algorithm job statistics on display
!
  implicit none
  type(hmc_poly_quark_wilson), intent(inout) :: this
  character(CHARLEN) :: name
  real(DP) :: hb_ave,gmp_ave

  hb_ave = get_average(this%update_q_inv)
  if (nodeid==0) then
  write(*,'("**** Quark # ",I3," PHMC")')this%myid

  call print(this%phys_param)

  write(*,'(" PHMC algorithm ")')
  write(*,'(" Heat-Bath for Pseudofermion BiCGStab iter")')
  write(*,'("                       (averaged) :",F12.4)') hb_ave
  endif

  if (SW_ON == this%switch_metropolis_test) then

    call print_statistics(this%gmp_status)

    if (this%nflavor == 2) then
      gmp_ave = get_average(this%gmp_inv)
      if (nodeid==0) then
      write(*,'(" Global Metropolis BiCGStab iter")')
      write(*,'("                       (averaged) :",F12.4)')gmp_ave
      endif
    endif
  endif

  return
end subroutine
Subroutine :
this :type(hmc_poly_quark_wilson), intent(inout)
iout :integer, intent(in)

Read PHMC quark algorithm parameters and Initialize algorithm parameters(ex. polynomial coeffs)

[Source]

subroutine read_wilson_quark_phmc(this,iout)
!
! Read PHMC quark algorithm parameters and Initialize algorithm parameters(ex. polynomial coeffs)
!
  use comlib
  implicit none
  type(hmc_poly_quark_wilson), intent(inout) :: this
  integer, intent(in) :: iout
  real(DP) :: rtmp0
  integer :: imb

  allocate(this%y)

  call new(this%update_q_inv)
  call new(this%gmp_inv)

  call read(this%phys_param,iout)
  if (nodeid==0) then
  read(iout,*)
  read(iout,*)this%update_q_inv%tol
  read(iout,*)this%gmp_inv%tol
  read(iout,*)this%depth_hopping
  read(iout,*)this%depth_clvtrlog
  read(iout,*)this%ptype
  read(iout,*)this%npoly
  read(iout,*)this%switch_metropolis_test
  endif
#ifndef _singlePU
  call comlib_bcast(this%update_q_inv%tol,0)
  call comlib_bcast(this%gmp_inv%tol,0)
  call comlib_bcast(this%depth_hopping,0)
  call comlib_bcast(this%depth_clvtrlog,0)
  call comlib_bcast(this%ptype,0)
  call comlib_bcast(this%npoly,0)
  call comlib_bcast(this%switch_metropolis_test,0)
#endif

  allocate(this%pcoef(0:this%npoly))
  allocate(this%zcoef(0:this%npoly))

  this%nflavor = get_nflavor(this%phys_param)

  select case (this%ptype)
  case (0) ! Chebyshev (c_i=(-1)^i)
    this%pcoef(0)=1.0_DP
    do imb=1,this%npoly
      this%pcoef(imb)=-this%pcoef(imb-1)
    enddo
    select case (this%nflavor)
    case (1) ! for 1 flavor
      call sqrtpoly_chev(this%npoly,this%pcoef,this%zcoef)
    case (2) ! for 2 flavor
      do imb=0,this%npoly
        rtmp0=this%pcoef(imb)
        this%zcoef(imb)=cmplx(rtmp0,0.0_DP,kind=DP)
      enddo
    case default
      call error_stop(" Error in quark param file (nflavor)")
    end select

  case (1) ! Adopted
    if (nodeid==0) then
    do imb=0,this%npoly
      read(iout,*)this%pcoef(imb)
    enddo
    endif
#ifndef _singlePU
    do imb=0,this%npoly
      call comlib_bcast(this%pcoef(imb),0)
    enddo
#endif
    select case (this%nflavor)
    case (1) ! for 1 flavor
      call sqrtpoly(this%npoly,this%pcoef,this%zcoef)
    case (2) ! for 2 flavor
      do imb=0,this%npoly
        rtmp0=this%pcoef(imb)
        this%zcoef(imb)=cmplx(rtmp0,0.0_DP,kind=DP)
      enddo
    case default
      call error_stop(" Error in quark param file (nflavor)")
    end select
  case default
    call error_stop(" Error in quark param file (ptype)")
  end select

  this%nsqrt=2

  return
end subroutine
Subroutine :
this :type(hmc_poly_quark_wilson), intent(inout)
: quark action/algorihtm parameters
u :type(vfield_gluon_wg), intent(in)

Nothing is done

[Source]

subroutine reweighting_wilson_quark_phmc(this,u)
!
! Nothing is done
!
  implicit none
  type(hmc_poly_quark_wilson), intent(inout) :: this ! quark action/algorihtm parameters
  type(vfield_gluon_wg),       intent(in)    :: u
  return
end subroutine
Subroutine :
this :type(hmc_poly_quark_wilson), intent(in)
iout :integer, intent(in)

Save Quark parameter on measurement configuration

[Source]

subroutine save_config_hmc_poly_quark_wilson(this,iout)
!
! Save Quark parameter on measurement configuration
!
  implicit none
  type(hmc_poly_quark_wilson), intent(in) :: this
  integer, intent(in) :: iout
  call save_config(this%phys_param,iout)
  return
end subroutine
Subroutine :
this :type(hmc_poly_quark_wilson), intent(inout)
id :integer, intent(in)

set quark # id

[Source]

subroutine set_id_hpqw(this,id)
!
! set quark # id
!
  implicit none
  type(hmc_poly_quark_wilson), intent(inout) :: this
  integer, intent(in) :: id
  this%myid = id
  call set_id(this%phys_param,id)
  return
end subroutine