Class quark_clover_class
In: QuarkCloverClass/quark_clover_class.F90
lattice_class error_class hmc_logfile_class hmc_status_class field_gauge_class field_fermion_class field_clover_class quark_wilson_class comlib counter_class quark_clover_class dot/f_142.png

Define Wilson-Clover quark action parameters and related routines

version

$Id: quark_clover_class.F90,v 1.25 2011/09/25 06:13:11 ishikawa Exp $

Methods

Included Modules

lattice_class error_class hmc_logfile_class hmc_status_class field_gauge_class field_fermion_class field_clover_class quark_wilson_class comlib counter_class

Public Instance methods

Subroutine :
ye :type(field_quark_eo_wg), intent(inout)
: even/odd site fermion field (overwritten)
this :type(quark_clover), intent(in)
: clover quark paramters

Multiply inverse clover term matrix on fermion field (overwritten) (even/odd sites only)

   ye <= Fee \ ye

[Source]

subroutine accum_mult_clover_term_eo(ye,this)
!
! Multiply inverse clover term matrix on fermion field (overwritten) (even/odd sites only)
!
!    ye <= Fee \ ye
!
  implicit none
  type(field_quark_eo_wg), intent(inout) :: ye   ! even/odd site fermion field (overwritten)
  type(quark_clover),      intent(in)    :: this ! clover quark paramters

  complex(DP):: a1,a2,a3,a4,a5,a6
  complex(DP):: b1,b2,b3,b4,b5,b6
  complex(DP):: p1,p2,p3,p4,p5,p6
  complex(DP):: m1,m2,m3,m4,m5,m6
  complex(DP):: y1,y2,y3,y4,y5,y6,y7,y8,y9,ya,yb,yc
  integer :: ix,iy,iz,itb,ieoxyz,itb0
  integer :: ieo

  ieo = ye%ieo

!$OMP PARALLEL DO &
!$OMP PRIVATE(ix,iy,iz,ieoxyz,itb,itb0, &
!$OMP         y1,y2,y3,y4,y5,y6,y7,y8,y9,ya,yb,yc, &
!$OMP         p1,p2,p3,p4,p5,p6,m1,m2,m3,m4,m5,m6, &
!$OMP         a1,a2,a3,a4,a5,a6,b1,b2,b3,b4,b5,b6)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+ieo,2)
#ifdef _SF
!===============================
! Set SF boundary condition
! itt = 0,1 and NTT+1 is zero
!===============================
    ye%s(0,iz,iy,ix)%y(:,:) = Z0
    do itb=1,NTH-ieoxyz
#else
    do itb=1-ieoxyz,NTH-ieoxyz
#endif
      itb0 = itb+ieoxyz

#define _YE_   ye%s(itb,iz,iy,ix)
#define _YDE_  ye%s(itb,iz,iy,ix)
#define _FCLV_ this%fclinv
#include "mult_clover_term_eo_kernel.h90"

    enddo

#ifdef _SF
    if (ieoxyz==1) then
      ye%s(NTH,iz,iy,ix)%y(:,:) = Z0
    endif
#endif

  enddo
  enddo
  enddo ! end of do ix,iy,iz

  return
end subroutine
Subroutine :
qp1 :type(quark_clover), intent(inout)
qp2 :type(quark_clover), intent(in)

[Source]

subroutine assign_quark_clover(qp1,qp2)
  implicit none
  type(quark_clover), intent(inout) :: qp1
  type(quark_clover), intent(in)    :: qp2
  call assign(qp1%quark_wilson,qp2%quark_wilson)
  qp1%csw                   =  qp2%csw
  qp1%fclinv                => qp2%fclinv
  qp1%clover_term_log       =  qp2%clover_term_log
  qp1%clover_term_log_fname =  qp2%clover_term_log_fname
  return
end subroutine
Subroutine :
yde :type(field_quark_eo_wg), intent(inout)
: even/odd site fermion field
this :type(quark_clover), intent(in)
: clover quark parameters
ye :type(field_quark_eo_wg), intent(in)
: even/odd site fermion field

Multiply inverse clover term matrix on fermion field (even/odd sites only)

   yde <= Fee \ ye

[Source]

subroutine assign_mult_clover_term_eo(yde,this,ye)
!
! Multiply inverse clover term matrix on fermion field (even/odd sites only)
!
!    yde <= Fee \ ye
!
  implicit none
  type(field_quark_eo_wg), intent(inout) :: yde  ! even/odd site fermion field
  type(quark_clover),      intent(in)    :: this ! clover quark parameters
  type(field_quark_eo_wg), intent(in)    ::  ye  ! even/odd site fermion field

  complex(DP):: a1,a2,a3,a4,a5,a6
  complex(DP):: b1,b2,b3,b4,b5,b6
  complex(DP):: p1,p2,p3,p4,p5,p6
  complex(DP):: m1,m2,m3,m4,m5,m6
  complex(DP):: y1,y2,y3,y4,y5,y6,y7,y8,y9,ya,yb,yc
  integer :: ix,iy,iz,itb,ieoxyz,ieo,itb0

  ieo     = ye%ieo
  yde%ieo = ye%ieo

!***********************************************************
!* Calc:
!*        
!*  yde = F^-1ee ye
!*      
!***********************************************************
!$OMP PARALLEL DO &
!$OMP PRIVATE(ix,iy,iz,ieoxyz,itb,itb0, &
!$OMP         a1,a2,a3,a4,a5,a6,b1,b2,b3,b4,b5,b6, &
!$OMP         p1,p2,p3,p4,p5,p6,m1,m2,m3,m4,m5,m6, &
!$OMP         y1,y2,y3,y4,y5,y6,y7,y8,y9,ya,yb,yc)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+yde%ieo,2)
#ifdef _SF
!===============================
! Set SF boundary condition
! itt = 0,1 and NTT+1 is zero
!===============================
    yde%s(0,iz,iy,ix)%y(:,:) = Z0
    do itb=1,NTH-ieoxyz
#else
    do itb=1-ieoxyz,NTH-ieoxyz
#endif
      itb0 = itb+ieoxyz

#define _YE_    ye%s(itb,iz,iy,ix)
#define _YDE_  yde%s(itb,iz,iy,ix)
#define _FCLV_ this%fclinv
#include "mult_clover_term_eo_kernel.h90"

    enddo

#ifdef _SF
    if (ieoxyz==1) then
      yde%s(NTH,iz,iy,ix)%y(:,:) = Z0
    endif
#endif

  enddo
  enddo
  enddo ! end of do ix,iy,iz

  return
end subroutine
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

[Source]

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

  type(field_quark_eo_wg), allocatable :: Myo
  real(DP) :: kappa2,kappa

  allocate(Myo)
  kappa = get_kappa(this%quark_wilson)
  kappa2=-(kappa**2)

!***********************************
!*  Myo = F^-1oo Moe ye
!*  yde = F^-1ee Meo Myo = F^-1ee Meo F^-1oo Moe ye
!***********************************
  call assign_mult_hop(this,Myo, ye,u)
  call assign_mult_hop(this,yde,Myo,u)

!***********************************
!*  yde = ye - kappa^2 yde
!*      = ye - kappa^2 Meo Moe ye
!***********************************
  call accum_mult_add(yde,kappa2,ye)

  deallocate(Myo)

  return
end subroutine
Subroutine :
this :type(quark_clover), intent(in)
: clover quark parameter and inverse clover term
FMye :type(field_quark_eo_wg), intent(out)
: even/odd site fermion field
yo :type(field_quark_eo_wg), intent(inout)
: odd/even site fermion field
u :type(vfield_gluon_wg), intent(in)
: gauge field

Multiply clover term preconditioned Hopping operator (odd->even or even->odd hopping only)

  FMye <= Fee \ Meo yo

[Source]

subroutine assign_mult_clover_hopping_eo(this,FMye,yo,u)
!
! Multiply clover term preconditioned Hopping operator (odd->even or even->odd hopping only)
!
!   FMye <= Fee \ Meo yo
! 
  use counter_class
  implicit none
  type(quark_clover),      intent(in)    :: this  ! clover quark parameter and inverse clover term
  type(field_quark_eo_wg), intent(out)   :: FMye  ! even/odd site fermion field
  type(field_quark_eo_wg), intent(inout) ::   yo  ! odd/even site fermion field
  type(vfield_gluon_wg),   intent(in)    ::    u  ! gauge field

  call copy_boundary(yo)
  call mult_clv_hopping_tzyx_eo(FMye,yo,u,this%fclinv)
  call inc(mult_iter)

!  call assign_mult_hop(this%quark_wilson,FMye,yo,u)
!  call accum_mult_clover_term(FMye,this)

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

[Source]

subroutine delete_quark_clover(this)
  implicit none
  type(quark_clover), intent(inout) :: this
  call delete(this%quark_wilson)
  call delete(this%clover_term_log)
  return
end subroutine
Subroutine :

delete clover leaf term.

[Source]

subroutine delete_clover_leaf()
!
! delete clover leaf term.
!
  implicit none
  if (associated(m_ucl)) then
    deallocate(m_ucl)
  endif
  return
end subroutine
Subroutine :
this :type(quark_clover), intent(inout)

delete inverse clover term.

[Source]

subroutine delete_inverse_clover_term(this)
!
! delete inverse clover term.
!
  implicit none
  type(quark_clover), intent(inout) :: this
  if (associated(this%fclinv)) then
    deallocate(this%fclinv)
  endif
  return
end subroutine
Subroutine :
BB :type(vfield_gluon_wg), intent(inout)
: pre-force contirubution (dot{u}) accumulated.
XX :type(tfield_gluon_wg), intent(inout)
: pre-force contirubution (dot{sigma_{mu,nu}F_{mu,nu}) accumulated.
u :type(vfield_gluon_wg), intent(in)
: gauge field
wl34 :type(tfield_gluon_wg), intent(in)
: 2-link fields prepared by make_two_links
wl98 :type(tfield_gluon_wg), intent(in)
: 2-link fields prepared by make_two_links

Compute Force for clover term contributions.

  BB_{mu}(n) = Sum[vcl_{mu,nu}(n),nu=!mu]

        X4---<---3X       Xn       2X
         |       |         |       |
  vcl =  v       ^    -    ^       v
         |       |         |       |
        Xn       2X       X6---<---5X

  X : X_{mu,nu}(n) - X^+_{mu,nu}(n),
  where
      X_{mu,nu}(n) = \sigma_{mu,nu} (force contributions at site n)

[Source]

subroutine force_clover_terms(BB,XX,u,wl34,wl98)
!
! Compute Force for clover term contributions.
!
!   BB_{mu}(n) = Sum[vcl_{mu,nu}(n),nu=!mu]
!
!         X4---<---3X       Xn       2X
!          |       |         |       |         
!   vcl =  v       ^    -    ^       v
!          |       |         |       |         
!         Xn       2X       X6---<---5X
!
!   X : X_{mu,nu}(n) - X^+_{mu,nu}(n),
!   where
!       X_{mu,nu}(n) = \sigma_{mu,nu} (force contributions at site n)
!
  implicit none
  type(vfield_gluon_wg), intent(inout) :: BB      ! pre-force contirubution (\dot{u}) accumulated.
  type(tfield_gluon_wg), intent(inout) :: XX      ! pre-force contirubution (\dot{\sigma_{mu,nu}F_{mu,nu}) accumulated.
  type(vfield_gluon_wg), intent(in)    ::  u          ! gauge field
  type(tfield_gluon_wg), intent(in)    :: wl34,wl98   ! 2-link fields prepared by make_two_links

  integer :: ipl

  do ipl=1,NDIM*(NDIM-1)/2
    call copy_boundary(XX%eo(0)%munu(ipl))
    call copy_boundary(XX%eo(1)%munu(ipl))
  enddo

  call force_clover_terms_eo(BB%eo(0),XX%eo(0),XX%eo(1),u%eo(0),u%eo(1),wl34%eo(0),wl98%eo(1))
  call force_clover_terms_eo(BB%eo(1),XX%eo(1),XX%eo(0),u%eo(1),u%eo(0),wl34%eo(1),wl98%eo(0))

  return
end subroutine
Subroutine :
this :type(quark_clover), intent(in)
: quark parameters
XX :type(tfield_gluon_wg), intent(inout)
: pre-force contribution (dot{sigma_{mu,nu}F_{mu,nu}})

Calc force contiribution from log determinant of clover term derived from even/odd site preconditioning.

Action:

     SF = - Nf log[det[F]]

[Source]

subroutine force_clover_trlog(this,XX)
!
! Calc force contiribution from log determinant of clover term derived
! from even/odd site preconditioning.
!
! Action:
!
!      SF = - Nf log[det[F]]
!
  implicit none
  type(quark_clover),    intent(in)    :: this ! quark parameters
  type(tfield_gluon_wg), intent(inout) :: XX   ! pre-force contribution (\dot{\sigma_{mu,nu}F_{mu,nu}})

  call force_clover_trlog_eo(this,this%fclinv%eo(0),XX%eo(0))
  call force_clover_trlog_eo(this,this%fclinv%eo(1),XX%eo(1))

  return
end subroutine
Subroutine :
XX :type(tfield_gluon_wg), intent(inout)
: pre-force contribution for clvoer term (dot{T})
zickd8 :complex(DP), intent(in)
: force coefficient
fx :type(field_quark_wg), intent(in)
: pseudo fermion
fz :type(field_quark_wg), intent(in)
: pseudo fermion

Calc force from inverse clover term in even/odd preconditioned Wilson-Dirac(Clover) hopping operator

 XX(n) =  XX(n) + zickd8 * tr[ sigma_{mu,nu} fz(n) fx(n)^{+}]

where, sigma_{mu,nu} : Dirac rep.

[Source]

subroutine force_hmc_clover(XX,zickd8,fx,fz)
!
! Calc force from inverse clover term in even/odd preconditioned Wilson-Dirac(Clover) hopping operator
!                         
!  XX(n) =  XX(n) + zickd8 * tr[ sigma_{mu,nu} fz(n) fx(n)^{+}]
!                         
! where,   sigma_{mu,nu} : Dirac rep.
!
  implicit none
  type(tfield_gluon_wg), intent(inout) :: XX   ! pre-force contribution for clvoer term (\dot{T})
  complex(DP),           intent(in) :: zickd8  ! force coefficient
  type(field_quark_wg),  intent(in) ::  fx     ! pseudo fermion
  type(field_quark_wg),  intent(in) ::  fz     ! pseudo fermion

  call force_hmc_clover_eo(XX%eo(0),zickd8,fx%eo(0),fz%eo(0))
  call force_hmc_clover_eo(XX%eo(1),zickd8,fx%eo(1),fz%eo(1))

  return
end subroutine
Function :
csw :real(DP)
this :type(quark_clover), intent(in)

[Source]

function get_csw(this) result(csw)
  implicit none
  type(quark_clover), intent(in) :: this
  real(DP) :: csw
  csw = this%csw
  return
end function
get_id( this ) result(id)
Function :
id :integer
this :class(action_parameters), intent(in)

Original external subprogram is quark_wilson_class#get_id

Function :
kappa :real(DP)
this :type(quark_clover), intent(in)

[Source]

function get_kappa_quark_clover(this) result(kappa)
  implicit none
  type(quark_clover), intent(in) :: this
  real(DP) :: kappa
  kappa = get_kappa(this%quark_wilson)
  return
end function
Function :
nf :integer
this :type(quark_clover), intent(in)

[Source]

function get_nflavor_quark_clover(this) result(nf)
  implicit none
  type(quark_clover), intent(in) :: this
  integer :: nf
  nf = get_nflavor(this%quark_wilson)
  return
end function
Function :
ptr :type(field_clover_term), pointer
this :type(quark_clover), intent(in)

[Source]

function get_ptr_inverse_clover_term(this) result(ptr)
  implicit none
  type(quark_clover),       intent(in) :: this
  type(field_clover_term),  pointer    :: ptr
  ptr => this%fclinv
  return
end function
Subroutine :
this :type(quark_clover), intent(inout)
: quark parameters, and contains
 inverse clover term (in chiral rep., linear form)
action :real(DP), intent(out)
: action value

Calculate Hamiltonian of log determinant clover term. This also computes inverse clover term by side effect.

  Action :

    SF = - Nf log[det[F]]

[Source]

subroutine hamil_clover_trlog(this,action)
!
! Calculate Hamiltonian of log determinant clover term.
! This also computes inverse clover term by side effect.
!
!   Action :
!
!     SF = - Nf log[det[F]]
!
  implicit none
  type(quark_clover), intent(inout) :: this    !  quark parameters, and contains 
                                               !  inverse clover term (in chiral rep., linear form)
  real(DP),           intent(out)   :: action  ! action value
  real(DP) :: logdetfcl

!**************************************************
!* calc. inverse clover term and log[det[Fcl]]
!**************************************************
  call make_inverse_clover_term(this,1,logdetfcl)

  action=-logdetfcl*get_nflavor(this%quark_wilson)

  return
end subroutine
Subroutine :
u :type(vfield_gluon_wg), intent(in)
: gauge link
wl34 :type(tfield_gluon_wg), intent(in)
: 2-links prepared by make_two_links
wl98 :type(tfield_gluon_wg), intent(in)
: 2-links prepared by make_two_links

Calculate clover leaf (field strength) for Clover fermion

clover field strength is class/module variable and is calculated using the two-link fields.

[Source]

subroutine make_clover_leaf(u,wl34,wl98)
!
! Calculate clover leaf (field strength) for Clover fermion
!
! clover field strength is class/module variable and is calculated using the two-link fields.
!
  implicit none
  type(vfield_gluon_wg), intent(in) :: u         ! gauge link
  type(tfield_gluon_wg), intent(in) :: wl34,wl98 ! 2-links prepared by make_two_links

  integer :: iflag

  iflag = 0
  if (.not.associated(m_ucl)) then
    allocate(m_ucl)
  endif
  call make_clover_leaf_eo(m_ucl%eo(0),u%eo(0),u%eo(1),wl34%eo(0),wl34%eo(1),wl98%eo(0),wl98%eo(1))
  call make_clover_leaf_eo(m_ucl%eo(1),u%eo(1),u%eo(0),wl34%eo(1),wl34%eo(0),wl98%eo(1),wl98%eo(0))
  return
end subroutine
Subroutine :
this :type(quark_clover), intent(inout)
: quark parameter contains inverse clover term (in chiral rep. linear storage)
idet :integer, intent(in)
: if idet == 1 then logdetfcl is computed.
logdetfcl :real(DP), intent(out)
: log determinat of clover term.

Calculate inverse clover term matrix and log-determinant of clover term.

  fcl    = [ 1 - csw kappa/2 sigma_{mu,nu} F_{mu,nu}]
  fclinv =  fcl^-1
  log[det[fcl]]

This uses module variable, m_ucl : clover leaf field.

  [ 1 - csw kappa/2 sigma_{mu,nu} F_{mu,nu}]^-1  = diag(fclinv(1),fclinv(2)),

  where

    |   1   2   3   4   5   6  |
    |       7   8   9  10  11  |
    |          12  13  14  15  |
    |              16  17  18  | * 0.5
    |                  19  20  |
    |                      21  |

  are stored in fclinv(i).

 Note: 0.5 is multiplied to simplify conversion from chiral to Dirac rep.

[Source]

subroutine make_inverse_clover_term(this,idet,logdetfcl)
!
! Calculate inverse clover term matrix and log-determinant of clover term.
!
!   fcl    = [ 1 - csw kappa/2 sigma_{mu,nu} F_{mu,nu}]
!   fclinv =  fcl^-1
!   log[det[fcl]]
!
! This uses module variable, m_ucl : clover leaf field.
!
!   [ 1 - csw kappa/2 sigma_{mu,nu} F_{mu,nu}]^-1  = diag(fclinv(1),fclinv(2)),
!
!   where
!
!     |   1   2   3   4   5   6  |
!     |       7   8   9  10  11  |
!     |          12  13  14  15  |
!     |              16  17  18  | * 0.5
!     |                  19  20  |
!     |                      21  |
! 
!   are stored in fclinv(i).
!
!  Note: 0.5 is multiplied to simplify conversion from chiral to Dirac rep.
!
  use comlib
  implicit none
  type(quark_clover), intent(inout) :: this      ! quark parameter contains 
                                                 ! inverse clover term (in chiral rep. linear storage)
 
  integer,            intent(in)    :: idet      ! if idet == 1 then logdetfcl is computed.
  real(DP),           intent(out)   :: logdetfcl ! log determinat of clover term.

  type(f_c_m_eo), allocatable :: f(:),finv(:)
  real(DP) :: logdetf1cl(0:1)
  real(DP) :: logdetf2cl(0:1)
  integer :: ieo,iout
  real(DP) :: kappa,csw

   iout = get_file_unit(this%clover_term_log)
  kappa = get_kappa(this%quark_wilson)
    csw = this%csw

  logdetfcl=0.0_DP

  if (.not.associated(m_ucl)) then
    call error_stop("Clover leaf field is not created, Stop!")
  endif
  if (.not.associated(this%fclinv)) then
    allocate(this%fclinv)
  endif
  call new(this%fclinv)

  allocate(f(2),finv(2))

  do ieo=0,1

!*****************************
!* calc clover term matrix in chiral rep.
!*
!*  fcl = diag[f1cl,f2cl] = 1 - csw kappa/2 sigma_{mu,nu}F_{mu,nu}
!*
!* solve inverse matrix (even part)
!*
!*****************************
    call clover_f1f2(kappa,csw,f(1),f(2),m_ucl%eo(ieo))
    call _CLV_SOLVER(f(1)%m,iout,idet,finv(1)%m,logdetf1cl(ieo))
    call _CLV_SOLVER(f(2)%m,iout,idet,finv(2)%m,logdetf2cl(ieo))
    finv(1)%ieo=f(1)%ieo
    finv(2)%ieo=f(2)%ieo

!*****************************
!* chiral full matrix -> chiral linear vector
!*****************************
    call conv_matrix_to_linear_eo(finv(1),finv(2),this%fclinv%eo(ieo))

  enddo

  if (idet == 1) then
    logdetfcl=logdetf1cl(0)+logdetf2cl(0) +logdetf1cl(1)+logdetf2cl(1)
#ifndef _singlePU
    call comlib_sumcast(logdetfcl)
#endif
  endif

  deallocate(f,finv)

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

[Source]

subroutine new_quark_clover(this,id)
  implicit none
  type(quark_clover), intent(inout) :: this
  integer,            intent(in)    :: id
  character(len=CHARLEN) :: str
  call new(this%quark_wilson,id)
  this%action_name = ACTION_NAME
  write(str,'(I2.2)')id
  this%clover_term_log_fname = TRIM(this%clover_term_log_fname)//TRIM(str)
  call new(this%clover_term_log,this%clover_term_log_fname)
  return
end subroutine
Subroutine :
this :type(quark_clover), intent(inout)

[Source]

subroutine print_quark_clover(this)
  implicit none
  type(quark_clover), intent(inout) :: this
  call print(this%quark_wilson)
  write(*,'(14X,"              Csw :",F12.8)')this%csw
  return
end subroutine
quark_clover
Derived Type :

Clover quark action parameters (extends from quark_wilson)

quark_wilson
Derived Type :

Wilson quark type extended from the above type (extends from action_parameters)

Original external subprogram is quark_wilson_class#quark_wilson

Subroutine :
this :type(quark_clover), intent(inout)
iout :integer, intent(in)

[Source]

subroutine read_quark_clover(this,iout)
  use comlib
  implicit none
  type(quark_clover), intent(inout) :: this
  integer, intent(in) :: iout
  call read(this%quark_wilson,iout)
  if (0 == nodeid) then
    read(iout,*)this%csw
  endif
#ifndef _singlePU
  call comlib_bcast(this%csw,0)
#endif
  return
end subroutine
Subroutine :
this :type(quark_clover), intent(in)
iout :integer, intent(in)

[Source]

subroutine save_config_quark_clover(this,iout)
  implicit none
  type(quark_clover), intent(in) :: this
  integer, intent(in) :: iout
  call save_config(this%quark_wilson,iout)
  write(iout)this%csw
  return
end subroutine
Subroutine :
this :type(quark_clover), intent(inout)
csw :real(DP), intent(in)

[Source]

subroutine set_csw(this,csw)
  implicit none
  type(quark_clover), intent(inout) :: this
  real(DP), intent(in) :: csw
  this%csw = csw
  return
end subroutine
set_gaussian_noise( p )
Subroutine :
p :type(vfield_gluon_wog), intent(inout)
: gauge momentum field

set Gaussian noise on canonical momentum (su(3) Lie algebra) of gauge field (SU(3) Lie group)

Original external subprogram is field_gauge_class#set_gaussian_noise

set_gaussian_noise( y )
Subroutine :
y :type(field_quark_wg), intent(inout)

set Gaussian noise on y

Original external subprogram is field_fermion_class#set_gaussian_noise

set_gaussian_noise( ye )
Subroutine :
ye :type(field_quark_eo_wg), intent(inout)

set Gaussian noise on even/odd sites only

Original external subprogram is field_fermion_class#set_gaussian_noise

set_id( this, id )
Subroutine :
this :class(action_parameters), intent(inout)
id :integer, intent(in)

Original external subprogram is quark_wilson_class#set_id

Subroutine :
this :type(quark_clover), intent(inout)
kappa :real(DP), intent(in)

[Source]

subroutine set_kappa_quark_clover(this,kappa)
  implicit none
  type(quark_clover), intent(inout) :: this
  real(DP), intent(in) :: kappa
  call set_kappa(this%quark_wilson,kappa)
  return
end subroutine
Subroutine :
this :type(quark_clover), intent(inout)
nf :integer, intent(in)

[Source]

subroutine set_nflavor_quark_clover(this,nf)
  implicit none
  type(quark_clover), intent(inout) :: this
  integer,            intent(in)    :: nf
  call set_nflavor(this%quark_wilson,nf)
  return
end subroutine