Class hmc_gluon_class
In: HmcGluonClass/hmc_gluon_class.F90
lattice_class error_class logfile_class hmc_logfile_class action_base_class field_gauge_class comlib file_tools_class print_status_class hmc_status_class hmc_gluon_class dot/f_256.png

Definition of HMC gluon action

This action contains RG action (1x1-Plaquette + 1x2-Rectangular).

$Id: hmc_gluon_class.F90,v 1.24 2011/09/25 06:12:47 ishikawa Exp $

Methods

Included Modules

lattice_class error_class logfile_class hmc_logfile_class action_base_class field_gauge_class comlib file_tools_class print_status_class hmc_status_class

Public Instance methods

clear( u )
Subroutine :
u :type(tfield_gluon_wg), intent(inout)

Original external subprogram is field_gauge_class#clear

clear( u )
Subroutine :
u :type(vfield_gluon_wg), intent(inout)

Original external subprogram is field_gauge_class#clear

clear( u )
Subroutine :
u :type(vfield_gluon_wog), intent(inout)

Original external subprogram is field_gauge_class#clear

clear( ue )
Subroutine :
ue :type(sfield_gluon_eo_wg), intent(inout)

Original external subprogram is field_gauge_class#clear

clear( ue )
Subroutine :
ue :type(sfield_gluon_eo_wog), intent(inout)

Original external subprogram is field_gauge_class#clear

clear( ue )
Subroutine :
ue :type(tfield_gluon_eo_wg), intent(inout)

Original external subprogram is field_gauge_class#clear

clear( ue )
Subroutine :
ue :type(vfield_gluon_eo_wg), intent(inout)

Original external subprogram is field_gauge_class#clear

clear( ue )
Subroutine :
ue :type(vfield_gluon_eo_wog), intent(inout)

Original external subprogram is field_gauge_class#clear

copy_sfg_time
Variable :
copy_sfg_time :type(timer), save
: holds gauge field boundary copy timing

Original external subprogram is field_gauge_class#copy_sfg_time

Subroutine :
this :type(hmc_gluon), intent(inout)

delete gluon action

[Source]

subroutine delete_gluon(this)
!
! delete gluon action
!
  implicit none
  type(hmc_gluon), intent(inout) :: this
  call delete(m_reunit_log)
  return
end subroutine
Function :
diff :real(DP)
p0 :type(vfield_gluon_wog), intent(in)
p1 :type(vfield_gluon_wog), intent(in)

compute difference norm of two momentum fields (for reversibility check)

 delta_p = (|p0 - p1|^2/(COL^2*NDIM*NTX*NTY*NTZ*NTT))^(1/2)
      p0 : momentum 1
      p1 : momentum 2

[Source]

function delta_p(p0,p1) result(diff)
!
! compute difference norm of two momentum fields (for reversibility check)
!
!  delta_p = (|p0 - p1|^2/(COL^2*NDIM*NTX*NTY*NTZ*NTT))^(1/2)
!       p0 : momentum 1
!       p1 : momentum 2
!
  use comlib
  implicit none
  real(DP) :: diff
  type(vfield_gluon_wog), intent(in) :: p0,p1

  integer :: ix,iy,iz,itb,ieoxyz,ieo,ic,jc,mu

  diff = 0.0_DP
  do ieo=0,1
  do mu=1,NDIM
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    do itb=1,NTH
      do jc=1,COL
      do ic=1,COL
        diff = diff + real( p0%eo(ieo)%mu(mu)%s(itb,iz,iy,ix)%u(ic,jc) -p0%eo(ieo)%mu(mu)%s(itb,iz,iy,ix)%u(ic,jc))**2 + aimag( p0%eo(ieo)%mu(mu)%s(itb,iz,iy,ix)%u(ic,jc) -p1%eo(ieo)%mu(mu)%s(itb,iz,iy,ix)%u(ic,jc))**2
      enddo
      enddo
    enddo
  enddo
  enddo
  enddo ! end of do ix,iy,iz
  enddo ! end of do mu
  enddo
#ifndef _singlePU
  call comlib_sumcast(diff)
#endif

  diff = sqrt(diff/(COL**2*NDIM)/(NTT*NTZ*NTY*NTX))

  return
end function
Function :
diff :real(DP)
u0 :type(vfield_gluon_wg), intent(in)
u1 :type(vfield_gluon_wg), intent(in)

compute diferrence norm of two gauge fields (for reversibility check)

 delta_u = (|u0 - u1|^2/(COL^2*NDIM*NTX*NTY*NTZ*NTT))^(1/2)
      u0 : link 1
      u1 : link 2

[Source]

function delta_u(u0,u1) result(diff)
!
! compute diferrence norm of two gauge fields (for reversibility check)
!
!  delta_u = (|u0 - u1|^2/(COL^2*NDIM*NTX*NTY*NTZ*NTT))^(1/2)
!       u0 : link 1
!       u1 : link 2
!
  use comlib
  implicit none
  real(DP) :: diff
  type(vfield_gluon_wg), intent(in) :: u0,u1

  integer :: ix,iy,iz,itb,ieoxyz,ieo,ic,jc,mu

  diff = 0.0_DP
  do ieo=0,1
  do mu=1,NDIM
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+iz+iy+ix+ieo,2)
    do itb=1-ieoxyz,NTH-ieoxyz
      do jc=1,COL
      do ic=1,COL
        diff = diff + real( u0%eo(ieo)%mu(mu)%s(itb,iz,iy,ix)%u(ic,jc) -u0%eo(ieo)%mu(mu)%s(itb,iz,iy,ix)%u(ic,jc))**2 + aimag( u0%eo(ieo)%mu(mu)%s(itb,iz,iy,ix)%u(ic,jc) -u1%eo(ieo)%mu(mu)%s(itb,iz,iy,ix)%u(ic,jc))**2
      enddo
      enddo
    enddo
  enddo
  enddo
  enddo ! end of do ix,iy,iz
  enddo ! end of do mu
  enddo
#ifndef _singlePU
  call comlib_sumcast(diff)
#endif

  diff = sqrt(diff/(COL**2*NDIM)/(NTT*NTZ*NTY*NTX))

  return
end function
Subroutine :
this :type(hmc_gluon), intent(in)
: gauge action parameters
BB :type(vfield_gluon_wg), intent(inout)
: gauge pre-force field (accumulated)
idepth :integer, intent(in)
: MD depth for recursive MD integrator
u :type(vfield_gluon_wg), intent(in)
: gauge link fields
wl34 :type(tfield_gluon_wg), intent(in)
: 2-link field prepared by make_two_links
wl98 :type(tfield_gluon_wg), intent(in)
: 2-link field prepared by make_two_links

Calculate MD Force from gauge action.

[Source]

subroutine force_gluon(this,BB,idepth,u,wl34,wl98)
!
! Calculate MD Force from gauge action.
!
  implicit none
  type(hmc_gluon),       intent(in)    :: this     ! gauge action parameters
  type(vfield_gluon_wg), intent(inout) :: BB       ! gauge pre-force field (accumulated)
  integer,               intent(in)    :: idepth   ! MD depth for recursive MD integrator
  type(vfield_gluon_wg), intent(in)    :: u        ! gauge link fields
  type(tfield_gluon_wg), intent(in)    :: wl34     ! 2-link field prepared by make_two_links
  type(tfield_gluon_wg), intent(in)    :: wl98     ! 2-link field prepared by make_two_links

  if (idepth == this%depth) then
#ifdef _FLAT_
    call force_g_eo_flat(this,BB%eo(0),u%eo(0),u%eo(1),wl34%eo(0),wl34%eo(1),wl98%eo(0),wl98%eo(1))
    call force_g_eo_flat(this,BB%eo(1),u%eo(1),u%eo(0),wl34%eo(1),wl34%eo(0),wl98%eo(1),wl98%eo(0))
#else
    call force_gluon_eo(this,BB%eo(0),u%eo(0),u%eo(1),wl34%eo(0),wl34%eo(1),wl98%eo(0),wl98%eo(1))
    call force_gluon_eo(this,BB%eo(1),u%eo(1),u%eo(0),wl34%eo(1),wl34%eo(0),wl98%eo(1),wl98%eo(0))
#endif
  endif

  return
end subroutine
get_force_norm( F ) result(norm)
Function :
norm :real(DP)
F :type(vfield_gluon_wog), intent(in)
: MD force

return averaged force norm

 norm = SQRT(abs2(F)/(NDIM*SPACETIMEVOLUME))

Original external subprogram is field_gauge_class#get_force_norm

get_id( this ) result(id)
Function :
id :integer
this :class(action_parameters), intent(in)

Original external subprogram is action_base_class#get_id

gluon_rg
Derived Type :
beta :real(DP)
: beta
cRG(0:3) :real(DP)
: RG-gauge improvement coefficients
cRG0(0:3) :real(DP)
: RG-gauge improvement coefficients (original input)
u0 :real(DP)
: tadpole improvement coefficients
gluon_rg_sf
Derived Type :
ct :real(DP)
ct_rg :real(DP)
Subroutine :
this :type(hmc_gluon), intent(inout)
: HMC gauge action
ifi :integer, intent(in)
: 0 for initail momentum is generated, 1 for final, 2 for final
u :type(vfield_gluon_wg), intent(inout)
: gauge field
p :type(vfield_gluon_wog), intent(inout)
: momentum field
SG :real(DP), intent(out)
: gauge action value
SP :real(DP), intent(out)
: gauge kinetic term value
wl34 :type(tfield_gluon_wg), intent(in)
: 2-links field prepared by make_two_links
wl98 :type(tfield_gluon_wg), intent(in)
: 2-links field prepared by make_two_links

Calculate Action value for gluon action (with Plaquette + Rectangle)

 SG = beta/6 Sum[ c0   Re[Tr[1-Upl_{mu,nu}(n)]]
                 +c1 2 Re[Tr[1-Urt_{mu,nu}(n)]] ,
                            n=VOL, mu=1..NDIM, nu=1..NDIM, mu=!=nu]

[Source]

subroutine hamil_gluon(this,ifi,u,p,SG,SP,wl34,wl98)
!
! Calculate Action value for gluon action (with Plaquette + Rectangle)
! 
!  SG = beta/6 Sum[ c0   Re[Tr[1-Upl_{mu,nu}(n)]]
!                  +c1 2 Re[Tr[1-Urt_{mu,nu}(n)]] ,
!                             n=VOL, mu=1..NDIM, nu=1..NDIM, mu=!=nu]
! 
  use hmc_status_class
  implicit none
  type(hmc_gluon),        intent(inout) :: this  ! HMC gauge action
  integer,                intent(in)    :: ifi   ! 0 for initail momentum is generated, 1 for final, 2 for final
  type(vfield_gluon_wg),  intent(inout) :: u   ! gauge field
  type(vfield_gluon_wog), intent(inout) :: p   ! momentum field
  real(DP),               intent(out)   :: SG  ! gauge action value
  real(DP),               intent(out)   :: SP  ! gauge kinetic term value
  type(tfield_gluon_wg),  intent(in)    :: wl34 ! 2-links field prepared by make_two_links
  type(tfield_gluon_wg),  intent(in)    :: wl98 ! 2-links field prepared by make_two_links

  type(wloops_obj)   :: WLPe,WLPo
  real(DP)           :: SGe,SPe,SGo,SPo
  character(CHARLEN) :: str
  type(hmc_status) :: status

  if (ifi==0) call set_gaussian_noise(p) ! update p with Heatbath

  if (ifi==1) then
    call reunit_u_gluon(u)
    call copy_boundary(u)
  endif

  call hamil_gluon_eo(this,SGe,SPe,WLPe,p%eo(0),u%eo(1),wl34%eo(0),wl34%eo(1),wl98%eo(0))
  call hamil_gluon_eo(this,SGo,SPo,WLPo,p%eo(1),u%eo(0),wl34%eo(1),wl34%eo(0),wl98%eo(1))

  SG=SGe+SGo
  SP=SPe+SPo
  
  this%WLP%PLQ=WLPe%PLQ+WLPo%PLQ
  this%WLP%RCT=WLPe%RCT+WLPo%RCT

  write(str,'("@",I8,I2," Gauge: ID:",I3,"   SU:",E24.16,"   SP:",E24.16,"  PLQ:",E24.16)') get_trajectory_number(status),ifi,this%myid,SG,SP,this%WLP%PLQ
  call print_log_action(status,str)

  return
end subroutine
hmc_gluon
Derived Type :
file_name :character(CHARLEN), public
: filename for action parameters
gluon :type(gluon_rg_sf), public
gluon :type(gluon_rg), public
gluon :type(gluon_rg_sf), public
gluon :type(gluon_rg), public
WLP :type(wloops_obj), public

Gluon action parameters

make_two_links( u, wl34, wl98 )
Subroutine :
u :type(vfield_gluon_wg), intent(in)
: gauge field
wl34 :type(tfield_gluon_wg), intent(inout)
: even/odd site upper right 2-links
wl98 :type(tfield_gluon_wg), intent(inout)
: even/odd site lower left 2-links

Prepare 2-links field (tensor) for force contribution

    4---<---3     8       n
            |     |
      wl34  ^     ^  wl98
            |     |
    n       2     9---<---6    for mu < nu

For mu < nu cases, use wl34^+, wl98^+.

Original external subprogram is field_gauge_class#make_two_links

Subroutine :
this :type(hmc_gluon), intent(inout)
id :integer, intent(in)

Initialize/create gluon action

[Source]

subroutine new_gluon(this,id)
!
! Initialize/create gluon action
!
  implicit none
  type(hmc_gluon), intent(inout) :: this
  integer,         intent(in)    :: id
  this%file_name = REPEAT(" ",LEN(this%file_name))
  this%gluon%beta      = 0.0_DP
  this%gluon%cRG(:)    = 0.0_DP
  this%gluon%cRG0(:)   = 0.0_DP
  this%gluon%u0        = 1.0_DP
#ifdef _SF
  this%gluon%ct        = 0.0_DP
  this%gluon%ct_rg     = 0.0_DP
#endif
  this%myid = id
  call new(m_reunit_log,TRIM(m_reunit_log_fname))
  return
end subroutine
plaquette( u, plq )
Subroutine :
u :type(vfield_gluon_wg), intent(in)
: gauge field
plq :real(DP), intent(out)
: averaged plaquett

Calculate space-time-direction averaged plaquette

Original external subprogram is field_gauge_class#plaquette

Subroutine :
this :type(hmc_gluon), intent(in)

print out gluon action parameters on display

[Source]

subroutine print_gluon(this)
!
! print out gluon action parameters on display
!
  implicit none
  type(hmc_gluon), intent(in) :: this

  if (nodeid==0) then
    write(*,'(16X,"           Beta :",F12.6)') this%gluon%beta
  
#ifdef _RGGAUGE
    write(*,'(16X," RG-gauge  cRG0 :",E24.16)')this%gluon%cRG(0)
    write(*,'(16X,"           cRG1 :",E24.16)')this%gluon%cRG0(1)
    write(*,'(16X,"             u0 :",E24.16)')this%gluon%u0
    write(*,'(16X,"      cRG1/u0^2 :",E24.16)')this%gluon%cRG(1)
#endif

#ifdef _SF
    write(*,'(16X," SF boundary ct :",F12.6)') this%gluon%ct
#ifdef _RGGAUGE
    write(*,'(16X,"          ct_RG :",F12.6)') this%gluon%ct_rg
#endif
#endif

    write(*,'(16X," MD Force depth :",I3)')    this%depth
  endif

  return
end subroutine
print_copy_sfg_statistics( )
Subroutine :

Original external subprogram is field_gauge_class#print_copy_sfg_statistics

Subroutine :
this :type(hmc_gluon), intent(inout)

Read gluon action parameters (couplings etc.) form formatted file

[Source]

subroutine read_gluon(this)
!
! Read gluon action parameters (couplings etc.) form formatted file
!
  use comlib
  use file_tools_class
  use print_status_class
  implicit none
  type(hmc_gluon), intent(inout) :: this
  integer :: iout
  character(CHARLEN) :: char

  if (nodeid == 0) then
    write(*,'(80("-"))')
  endif

  iout = search_free_file_unit()

  if (nodeid==0) then
  open(iout,file=this%file_name,status='old',form='formatted')
  read(iout,*)this%gluon%beta
#ifdef _RGGAUGE
  read(iout,*)this%gluon%cRG0(1)
  read(iout,*)this%gluon%u0
#endif
#ifdef _SF
  read(iout,*)this%gluon%ct
#ifdef _RGGAUGE
  read(iout,*)this%gluon%ct_rg
#endif
#endif
  read(iout,*)this%depth
  close(iout)
  endif

#ifndef _singlePU
  call comlib_bcast(this%gluon%beta,0)
  call comlib_bcast(this%gluon%cRG0(1),0)
  call comlib_bcast(this%gluon%u0,0)
  call comlib_bcast(this%depth,0)
#ifdef _SF
  call comlib_bcast(this%gluon%ct,0)
  call comlib_bcast(this%gluon%ct_rg,0)
#endif
#endif

#ifdef _RGGAUGE
  this%gluon%cRG(0)=1.0d0-8.0d0*this%gluon%cRG0(1)
  this%gluon%cRG(1)=this%gluon%cRG0(1)/this%gluon%u0**2
#endif

  write(char,'(1X,A," for gauge action : read OK ",I3)') TRIM(this%file_name),nodeid
  call print_status(char)

  return
end subroutine
Subroutine :
u :type(vfield_gluon_wg), intent(inout)
reunit_log :type(logfile), optional, intent(inout)
hmc_reunit_log :type(hmc_logfile), optional, intent(inout)

Re-unitarize gauge field including ghost sites

 u : gauge link

[Source]

subroutine reunit_u_gluon(u,reunit_log,hmc_reunit_log)
!
! Re-unitarize gauge field including ghost sites
!
!  u : gauge link
!
  use comlib
  implicit none
  type(vfield_gluon_wg),       intent(inout) :: u
  type(logfile),     optional, intent(inout) :: reunit_log
  type(hmc_logfile), optional, intent(inout) :: hmc_reunit_log
  if (present(reunit_log)) then
    call reunit_u_gluon_eo(u%eo(0),reunit_log=reunit_log)
    call reunit_u_gluon_eo(u%eo(1),reunit_log=reunit_log)
  else if (present(hmc_reunit_log)) then
    call reunit_u_gluon_eo(u%eo(0),hmc_reunit_log=hmc_reunit_log)
    call reunit_u_gluon_eo(u%eo(1),hmc_reunit_log=hmc_reunit_log)
  else
    call reunit_u_gluon_eo(u%eo(0))
    call reunit_u_gluon_eo(u%eo(1))
  endif
  return
end subroutine
Subroutine :
this :type(hmc_gluon), intent(in)
iout :integer, intent(in)

Save gauge action parameters (couplings etc.) to measurement config file in unformatted form.

[Source]

subroutine save_config_gluon(this,iout)
!
! Save gauge action parameters (couplings etc.) to measurement config file 
! in unformatted form.
!
  implicit none
  type(hmc_gluon), intent(in) :: this
  integer, intent(in) :: iout
  write(iout)this%gluon%beta
#ifdef _RGGAUGE
  write(iout)this%gluon%cRG0(1)
  write(iout)this%gluon%u0
#endif
#ifdef _SF
  write(iout)this%gluon%ct
#ifdef _RGGAUGE
  write(iout)this%gluon%ct_rg
#endif
#endif
  return
end subroutine
set_id( this, id )
Subroutine :
this :class(action_parameters), intent(inout)
id :integer, intent(in)

Original external subprogram is action_base_class#set_id

set_identity( u )
Subroutine :
u :type(vfield_gluon_wg), intent(inout)
: gauge field

Initialize gauge field to unit matrix (cold config)

Original external subprogram is field_gauge_class#set_identity

tfield_gluon_wg
Derived Type :
eo(0:1) :type(tfield_gluon_eo_wg)

su(3) matrix symmetric tensor field with ghost sits

Original external subprogram is field_gauge_class#tfield_gluon_wg

update_p( p, dt, F )
Subroutine :
p :type(vfield_gluon_wog), intent(inout)
: gauge momentum field
dt :real(DP), intent(in)
: MD time step
F :type(vfield_gluon_wog), intent(in)
: gauge force field

update gauge momentum field using force F

 p_mu(n) = p_mu(n) - i dt F_mu(n)

 Force <= U *(D Hamil / D U) (<= fix me, or check!)

Original external subprogram is field_gauge_class#update_p

update_u( u, dt, p )
Subroutine :
u :type(vfield_gluon_wg), intent(inout)
: gauge link field
dt :real(DP), intent(in)
: MD step size
p :type(vfield_gluon_wog), intent(in)
: gauge momentum field

Update gauge field using canonical momentum :

   u_mu(dt)=exp(i p_mu(dt')dt) u_mu(0)

Original external subprogram is field_gauge_class#update_u

vfield_gluon_wg
Derived Type :
eo(0:1) :type(vfield_gluon_eo_wg)

su(3) matrix vector field with ghost sits

Original external subprogram is field_gauge_class#vfield_gluon_wg

vfield_gluon_wog
Derived Type :
eo(0:1) :type(vfield_gluon_eo_wog)

su(3) matrix vector field without ghost sits

Original external subprogram is field_gauge_class#vfield_gluon_wog

wloops_obj
Derived Type :
PLQ =0.0_DP :real(DP)
: PLQ = Sum[Re[1x1],vol,direction]/direction/Vol
RCT =0.0_DP :real(DP)
: RCT = Sum[Re[2x1+1x2],vol,direction]/direction/Vol
    = 2 Sum[Re[1x2],vol,direction]/direction/Vol

contains Wilson loop/plaquette values