Class field_gauge_class
In: FieldGaugeClass/field_gauge_class.F90
lattice_class timer_class error_class comlib sf_lattice_class field_gauge_class dot/f_154.png

This defines even/odd site orderd gauge fields and related routines.

Operations on fields are in linear algebra and are expressed with subroutine form. The subroutine name contains the fllowing strings expressing operations for the input arguments.

assign
substitution or data copy
accum
accumulation to first argument
conjg
take Hermitian conjugate
mult
multiplication of scalar on field or multiplication of two gauge fields
add
addition of two fields
sub
subtractoin of two fields
shift
shift cordinate of the field
abs2
squared norm of a field.

etc...


The subroutine name can contains several operations to do composite opertions. The order of subrutine arguments and operations in subroutine name is matched as fllows.

assign_add(u,v,w)
u assign v add w : u <= v + w
assign_sub(u,v,w)
u assign v sub w : u <= v - w
accum_sub(u,v)
u accum u sub v : u <= u - v
assign_mult(u,v,w)
u assign v mult w : u <= v * w
assign_shift(u,v,mu)
u(n) assign shift(v,mu) : u(n) <= v(n+mu)

etc...

Methods

Included Modules

lattice_class timer_class error_class comlib sf_lattice_class

Public Instance methods

Function :
my_abs2 :real(DP)
pe :type(sfield_gluon_eo_wog), intent(in)

return squared norm of gauge field (momentum)

[Source]

function abs2_sfg_eo_wog(pe) result(my_abs2)
!
! return squared norm of gauge field (momentum)
!
  implicit none
  type(sfield_gluon_eo_wog), intent(in) :: pe
  real(DP) :: my_abs2
  integer :: ix,iy,iz,itb,ic,jc

  my_abs2 = 0.0_DP
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,itb,ic,jc) REDUCTION(+:my_abs2)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    do itb=1,NTH
      do jc=1,COL
      do ic=1,COL
        my_abs2 = my_abs2+ real(pe%s(itb,iz,iy,ix)%u(ic,jc))**2 +aimag(pe%s(itb,iz,iy,ix)%u(ic,jc))**2
      enddo
      enddo
    enddo
  enddo
  enddo
  enddo
#ifndef _singlePU
  call comlib_sumcast(my_abs2)
#endif

  return
end function
Subroutine :
u2e :type(sfield_gluon_eo_wg), intent(inout)
u1e :type(sfield_gluon_eo_wg), intent(in)

u2e <= u2e + u1e

[Source]

subroutine accum_add_sfg_eo_wg(u2e,u1e)
!
! u2e <= u2e + u1e
!
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: u2e
  type(sfield_gluon_eo_wg), intent(in)    :: u1e
  integer :: ix,iy,iz,ieoxyz,itb,ic,jc

#ifdef _DEBUG_
  if (u1e%ieo /= u2e%ieo) then
    call error_stop("error with u1e%ieo /= u2e%ieo in accum_add_sfg_eo_wg.")
  endif
#endif
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ic,jc)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+u1e%ieo,2)
    do itb=1-ieoxyz,NTH-ieoxyz
      do jc=1,COL
      do ic=1,COL
        u2e%s(itb,iz,iy,ix)%u(ic,jc)=u2e%s(itb,iz,iy,ix)%u(ic,jc) +u1e%s(itb,iz,iy,ix)%u(ic,jc)
      enddo
      enddo
    enddo
  enddo
  enddo
  enddo

  return
end subroutine
Subroutine :
p2e :type(sfield_gluon_eo_wog), intent(inout)
p1e :type(sfield_gluon_eo_wog), intent(in)
coef :complex(DP), intent(in)
 p2e <= p2e + p1e * coef

[Source]

subroutine accum_add_cmult_sfg_eo_wog(p2e,p1e,coef)
!
!  p2e <= p2e + p1e * coef
!
  implicit none
  type(sfield_gluon_eo_wog), intent(inout) :: p2e
  type(sfield_gluon_eo_wog), intent(in)    :: p1e
  complex(DP), intent(in) :: coef
  integer :: ix,iy,iz,itb

  p2e%ieo=p1e%ieo
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,itb)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    do itb=1,NTH
      p2e%s(itb,iz,iy,ix)%u(:,:)=p2e%s(itb,iz,iy,ix)%u(:,:) +coef*p1e%s(itb,iz,iy,ix)%u(:,:)
    enddo
  enddo
  enddo
  enddo

  return
end subroutine
Subroutine :
u :type(sfield_gluon_wg), intent(inout)
coef :real(DP), intent(in)
itt0 :integer, intent(in)
itt1 :integer, optional, intent(in)

multiply scalar at time t on SU(3) field

ue(t,x,y,z) = ue(t,x,y,z)*coef for t = itt0 or (itt0 and itt1:optional)

[Source]

subroutine accum_mult_at_sf_boundary_u(u,coef,itt0,itt1)
!
! multiply scalar at time t on SU(3) field
!
! ue(t,x,y,z) = ue(t,x,y,z)*coef for t = itt0 or (itt0 and itt1:optional)
!
  implicit none
  type(sfield_gluon_wg),    intent(inout) :: u
  real(DP),                 intent(in) :: coef
  integer,                  intent(in) :: itt0
  integer, optional,        intent(in) :: itt1
  if (.not.present(itt1)) then
    call accum_mult_at_sf_boundary_u_eo(u%eo(0),coef,itt0)
    call accum_mult_at_sf_boundary_u_eo(u%eo(1),coef,itt0)
  else
    call accum_mult_at_sf_boundary_u_eo(u%eo(0),coef,itt0,itt1)
    call accum_mult_at_sf_boundary_u_eo(u%eo(1),coef,itt0,itt1)
  endif
  return
end subroutine
Subroutine :
ue :type(sfield_gluon_eo_wg), intent(inout)
coef :real(DP), intent(in)
itt0 :integer, intent(in)
itt1 :integer, optional, intent(in)

multiply scalar at time t on SU(3) field

ue(t,x,y,z) = ue(t,x,y,z)*coef for t = itt0 or (itt0 and itt1:optional)

[Source]

subroutine accum_mult_at_sf_boundary_u_eo(ue,coef,itt0,itt1)
!
! multiply scalar at time t on SU(3) field
!
! ue(t,x,y,z) = ue(t,x,y,z)*coef for t = itt0 or (itt0 and itt1:optional)
!
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: ue
  real(DP),                    intent(in) :: coef
  integer,                     intent(in) :: itt0
  integer, optional,           intent(in) :: itt1
  integer :: ieo
  integer :: ix,iy,iz,ieoxyz,ith0,ith1,jtt0,jtt1

  ieo = ue%ieo
  if (.not.present(itt1)) then
    ith0 = itt0/2
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,jtt0)
    do ix=1,NX
    do iy=1,NY
    do iz=1,NZ
      ieoxyz=mod(ipeo+ix+iy+iz+ieo,2)
      jtt0 = 2*ith0 + ieoxyz
      if (jtt0 == itt0) then
        ue%s(ith0,iz,iy,ix)%u(:,:) = ue%s(ith0,iz,iy,ix)%u(:,:)*coef
      endif
    enddo
    enddo
    enddo
  else
    ith0 = itt0/2
    ith1 = itt1/2
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,jtt0,jtt1)
    do ix=1,NX
    do iy=1,NY
    do iz=1,NZ
      ieoxyz=mod(ipeo+ix+iy+iz+ieo,2)
      jtt0 = 2*ith0 + ieoxyz
      if (jtt0 == itt0) then
        ue%s(ith0,iz,iy,ix)%u(:,:) = ue%s(ith0,iz,iy,ix)%u(:,:)*coef
      endif
      jtt1 = 2*ith1 + ieoxyz
      if (jtt1 == itt1) then
        ue%s(ith1,iz,iy,ix)%u(:,:) = ue%s(ith1,iz,iy,ix)%u(:,:)*coef
      endif
    enddo
    enddo
    enddo
  endif
  return
end subroutine
Subroutine :
u2e :type(sfield_gluon_eo_wg), intent(inout)
u1e :type(sfield_gluon_eo_wg), intent(in)

u2e <= u2e + u1e

[Source]

subroutine accum_sub_sfg_eo_wg(u2e,u1e)
!
! u2e <= u2e + u1e
!
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: u2e
  type(sfield_gluon_eo_wg), intent(in)    :: u1e
  integer :: ix,iy,iz,ieoxyz,itb,ic,jc

#ifdef _DEBUG_
  if (u1e%ieo /= u2e%ieo) then
    call error_stop("error with u1e%ieo /= u2e%ieo in accum_sub_sfg_eo_wg.")
  endif
#endif
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ic,jc)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+u1e%ieo,2)
    do itb=1-ieoxyz,NTH-ieoxyz
      do jc=1,COL
      do ic=1,COL
        u2e%s(itb,iz,iy,ix)%u(ic,jc)=u2e%s(itb,iz,iy,ix)%u(ic,jc) -u1e%s(itb,iz,iy,ix)%u(ic,jc)
      enddo
      enddo
    enddo
  enddo
  enddo
  enddo

  return
end subroutine
Subroutine :
u2e :type(sfield_gluon_eo_wg), intent(inout)
u1e :type(sfield_gluon_eo_wg), intent(in)
 u2e <= u1e

[Source]

subroutine assign_sfg_eo_wg(u2e,u1e)
  !
  !  u2e <= u1e
  !
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: u2e
  type(sfield_gluon_eo_wg), intent(in)    :: u1e
  integer :: ix,iy,iz,ieoxyz,itb
  u2e%ieo=u1e%ieo
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+u2e%ieo,2)
    do itb=1-ieoxyz,NTH-ieoxyz
      u2e%s(itb,iz,iy,ix)=u1e%s(itb,iz,iy,ix)
    enddo
  enddo
  enddo
  enddo
  return
end subroutine
Subroutine :
u3e :type(sfield_gluon_eo_wg), intent(inout)
u1e :type(sfield_gluon_eo_wg), intent(in)
u2e :type(sfield_gluon_eo_wg), intent(in)

u3e <= u1e + u2e

[Source]

subroutine assign_add_sfg_eo_wg(u3e,u1e,u2e)
!
! u3e <= u1e + u2e
!
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: u3e
  type(sfield_gluon_eo_wg), intent(in)    :: u1e,u2e
  integer :: ix,iy,iz,ieoxyz,itb,ic,jc

#ifdef _DEBUG_
  if (u1e%ieo /= u2e%ieo) then
    call error_stop("error with u1e%ieo /= u2e%ieo in assign_add_sfg_eo_wg.")
  endif
#endif
  u3e%ieo=u1e%ieo
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ic,jc)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+u1e%ieo,2)
    do itb=1-ieoxyz,NTH-ieoxyz
      do jc=1,COL
      do ic=1,COL
        u3e%s(itb,iz,iy,ix)%u(ic,jc)=u1e%s(itb,iz,iy,ix)%u(ic,jc) +u2e%s(itb,iz,iy,ix)%u(ic,jc)
      enddo
      enddo
    enddo
  enddo
  enddo
  enddo

  return
end subroutine
Subroutine :
u2e :type(sfield_gluon_eo_wg), intent(inout)
u1e :type(sfield_gluon_eo_wg), intent(in)
 u2e <= u1e

[Source]

subroutine assign_conjg_sfg_eo_wg(u2e,u1e)
!
!  u2e <= u1e
!
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: u2e 
  type(sfield_gluon_eo_wg), intent(in)    :: u1e 
  integer :: ix,iy,iz,ieoxyz,itb,ic,jc

  u2e%ieo=u1e%ieo
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ic,jc)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+u1e%ieo,2)
    do itb=1-ieoxyz,NTH-ieoxyz
      do ic=1,COL
      do jc=1,COL
        u2e%s(itb,iz,iy,ix)%u(ic,jc)=conjg(u1e%s(itb,iz,iy,ix)%u(jc,ic))
      enddo
      enddo
    enddo
  enddo
  enddo
  enddo

  return
end subroutine
Subroutine :
u3e :type(sfield_gluon_eo_wg), intent(inout)
u1e :type(sfield_gluon_eo_wg), intent(in)
u2e :type(sfield_gluon_eo_wg), intent(in)
 u3e <= (u1e)^+ * u2e

[Source]

subroutine assign_conjg1_mult_sfg_eo_wg(u3e,u1e,u2e)
!
!  u3e <= (u1e)^+ * u2e
!
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: u3e
  type(sfield_gluon_eo_wg), intent(in)    :: u1e,u2e
  integer :: ix,iy,iz,ieoxyz,itb,ic,jc

#ifdef _DEBUG_
  if (u1e%ieo /= u2e%ieo) then
    call error_stop("error with u1e%ieo /= u2e%ieo in assign_conjg1_mult_sfg_eo_wg.")
  endif
#endif
  u3e%ieo=u1e%ieo
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ic,jc)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+u3e%ieo,2)
    do itb=1-ieoxyz,NTH-ieoxyz
      do jc=1,COL
      do ic=1,COL
        u3e%s(itb,iz,iy,ix)%u(ic,jc)=conjg(u1e%s(itb,iz,iy,ix)%u(1,ic)) *      u2e%s(itb,iz,iy,ix)%u(1,jc) +conjg(u1e%s(itb,iz,iy,ix)%u(2,ic)) *      u2e%s(itb,iz,iy,ix)%u(2,jc) +conjg(u1e%s(itb,iz,iy,ix)%u(3,ic)) *      u2e%s(itb,iz,iy,ix)%u(3,jc)
      enddo
      enddo
    enddo
  enddo
  enddo
  enddo

  return
end subroutine
Subroutine :
u3e :type(sfield_gluon_eo_wg), intent(inout)
u1e :type(sfield_gluon_eo_wg), intent(in)
u2e :type(sfield_gluon_eo_wg), intent(in)
 u3e <= u1e * (u2e)^+

[Source]

subroutine assign_conjg2_mult_sfg_eo_wg(u3e,u1e,u2e)
!
!  u3e <= u1e * (u2e)^+
!
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: u3e
  type(sfield_gluon_eo_wg), intent(in)    :: u1e,u2e
  integer :: ix,iy,iz,ieoxyz,itb,ic,jc

#ifdef _DEBUG_
  if (u1e%ieo /= u2e%ieo) then
    call error_stop("error with u1e%ieo /= u2e%ieo in assign_conjg2_mult_sfg_eo_wg.")
  endif
#endif
  u3e%ieo=u1e%ieo
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ic,jc)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+u3e%ieo,2)
    do itb=1-ieoxyz,NTH-ieoxyz
      do jc=1,COL
      do ic=1,COL
        u3e%s(itb,iz,iy,ix)%u(ic,jc)=      u1e%s(itb,iz,iy,ix)%u(ic,1) *conjg(u2e%s(itb,iz,iy,ix)%u(jc,1)) +      u1e%s(itb,iz,iy,ix)%u(ic,2) *conjg(u2e%s(itb,iz,iy,ix)%u(jc,2)) +      u1e%s(itb,iz,iy,ix)%u(ic,3) *conjg(u2e%s(itb,iz,iy,ix)%u(jc,3))
      enddo
      enddo
    enddo
  enddo
  enddo
  enddo

  return
end subroutine
Subroutine :
u1o :type(sfield_gluon_eo_wg), intent(inout)
u1e :type(sfield_gluon_eo_wg), intent(in)
mu :integer, intent(in)
 u1o(n) <= (u1e(n+mu))^+

[Source]

subroutine assign_conjg_shift1_sfg_eo_wg(u1o,u1e,mu)
!
!  u1o(n) <= (u1e(n+mu))^+
!
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: u1o
  type(sfield_gluon_eo_wg), intent(in)    :: u1e
  integer, intent(in) :: mu
  integer :: ix,iy,iz,ieoxyz,itb
  integer :: ix2,iy2,iz2,itb2,ic,jc

  u1o%ieo=1-u1e%ieo
  select case (mu)
  case (1:NDIM)
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ix2,iy2,iz2,itb2,ic,jc)
    do ix=1,NX
    do iy=1,NY
    do iz=1,NZ
      ieoxyz=mod(ipeo+ix+iy+iz+u1o%ieo,2)
      ix2=ix+isx(mu)
      iy2=iy+isy(mu)
      iz2=iz+isz(mu)
      do itb=1-ieoxyz,NTH-ieoxyz
        itb2=itb+ieoxyz*ist(mu)
        do jc=1,COL
        do ic=1,COL
          u1o%s(itb,iz,iy,ix)%u(ic,jc)=conjg(u1e%s(itb2,iz2,iy2,ix2)%u(jc,ic))
        enddo
        enddo
      enddo
    enddo
    enddo
    enddo
  case (-NDIM:-1)
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ix2,iy2,iz2,itb2,ic,jc)
    do ix=1,NX
    do iy=1,NY
    do iz=1,NZ
      ieoxyz=mod(ipeo+ix+iy+iz+u1o%ieo,2)
      ix2=ix-isx(-mu)
      iy2=iy-isy(-mu)
      iz2=iz-isz(-mu)
      do itb=1-ieoxyz,NTH-ieoxyz
        itb2=itb-(1-ieoxyz)*ist(-mu)
        do jc=1,COL
        do ic=1,COL
          u1o%s(itb,iz,iy,ix)%u(ic,jc)=conjg(u1e%s(itb2,iz2,iy2,ix2)%u(jc,ic))
        enddo
        enddo
      enddo
    enddo
    enddo
    enddo
  case default
    call error_stop("error with |mu|>NDIM in assign_conjg_shift1_sfg_eo_wg.")
  end select

  return
end subroutine
Subroutine :
u2e :type(sfield_gluon_eo_wg), intent(inout)
u1e :type(sfield_gluon_eo_wg), intent(in)
mu :integer, intent(in)
nu :integer, intent(in)

u2e(n) <= (u1e(n+mu+nu))^+ (mu /= nu)

[Source]

subroutine assign_conjg_shift2_sfg_eo_wg(u2e,u1e,mu,nu)
!
! u2e(n) <= (u1e(n+mu+nu))^+  (mu /= nu)
!
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: u2e
  type(sfield_gluon_eo_wg), intent(in)    :: u1e
  integer, intent(in) :: mu,nu
  integer :: ix,iy,iz,ieoxyz,itb
  integer :: ix2,iy2,iz2,itb2,ic,jc

  if (mu == nu) then
    call error_stop("error with mu == nu in assign_conjg_shift2_sfg_eo_wg.")
  endif
  if ((ABS(mu) > NDIM) .OR. (ABS(nu) > NDIM)) then
    call error_stop("error with |mu| > NDIM or |nu| > NDIM in assign_conjg_shift2_sfg_eo_wg.")
  endif

  u2e%ieo=u1e%ieo
  if ((mu < 0) .AND. (nu < 0)) then
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ix2,iy2,iz2,itb2,ic,jc)
    do ix=1,NX
    do iy=1,NY
    do iz=1,NZ
      ieoxyz=mod(ipeo+ix+iy+iz+u2e%ieo,2)
      ix2=ix-isx(-mu)-isx(-nu)
      iy2=iy-isy(-mu)-isy(-nu)
      iz2=iz-isz(-mu)-isz(-nu)
      do itb=1-ieoxyz,NTH-ieoxyz
        itb2=itb-(1-ieoxyz)*ist(-mu)-(1-ieoxyz)*ist(-nu)
        do jc=1,COL
        do ic=1,COL
          u2e%s(itb,iz,iy,ix)%u(ic,jc)=conjg(u1e%s(itb2,iz2,iy2,ix2)%u(jc,ic))
        enddo
        enddo
      enddo
    enddo
    enddo
    enddo
  else if ((mu > 0) .AND. (nu < 0)) then
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ix2,iy2,iz2,itb2,ic,jc)
    do ix=1,NX
    do iy=1,NY
    do iz=1,NZ
      ieoxyz=mod(ipeo+ix+iy+iz+u2e%ieo,2)
      ix2=ix+isx(mu)-isx(-nu)
      iy2=iy+isy(mu)-isy(-nu)
      iz2=iz+isz(mu)-isz(-nu)
      do itb=1-ieoxyz,NTH-ieoxyz
        itb2=itb+(  ieoxyz)*ist(mu)-(1-ieoxyz)*ist(-nu)
        do jc=1,COL
        do ic=1,COL
          u2e%s(itb,iz,iy,ix)%u(ic,jc)=conjg(u1e%s(itb2,iz2,iy2,ix2)%u(jc,ic))
        enddo
        enddo
      enddo
    enddo
    enddo
    enddo
  else if ((mu < 0) .AND. (nu > 0)) then
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ix2,iy2,iz2,itb2,ic,jc)
    do ix=1,NX
    do iy=1,NY
    do iz=1,NZ
      ieoxyz=mod(ipeo+ix+iy+iz+u2e%ieo,2)
      ix2=ix-isx(-mu)+isx(nu)
      iy2=iy-isy(-mu)+isy(nu)
      iz2=iz-isz(-mu)+isz(nu)
      do itb=1-ieoxyz,NTH-ieoxyz
        itb2=itb-(1-ieoxyz)*ist(-mu)+(  ieoxyz)*ist(nu)
        do jc=1,COL
        do ic=1,COL
          u2e%s(itb,iz,iy,ix)%u(ic,jc)=conjg(u1e%s(itb2,iz2,iy2,ix2)%u(jc,ic))
        enddo
        enddo
      enddo
    enddo
    enddo
    enddo
  else if ((mu > 0) .AND. (nu > 0)) then
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ix2,iy2,iz2,itb2,ic,jc)
    do ix=1,NX
    do iy=1,NY
    do iz=1,NZ
      ieoxyz=mod(ipeo+ix+iy+iz+u2e%ieo,2)
      ix2=ix+isx(mu)+isx(nu)
      iy2=iy+isy(mu)+isy(nu)
      iz2=iz+isz(mu)+isz(nu)
      do itb=1-ieoxyz,NTH-ieoxyz
        itb2=itb+(  ieoxyz)*ist(mu)+(  ieoxyz)*ist(nu)
        do jc=1,COL
        do ic=1,COL
          u2e%s(itb,iz,iy,ix)%u(ic,jc)=conjg(u1e%s(itb2,iz2,iy2,ix2)%u(jc,ic))
        enddo
        enddo
      enddo
    enddo
    enddo
    enddo
  endif

  return
end subroutine
Subroutine :
p :type(vfield_gluon_wog), intent(inout)
u :type(vfield_gluon_wg), intent(in)
BB :type(vfield_gluon_wg), intent(in)
 pe <= u1e * u2e

[Source]

subroutine assign_mult_vfg_wg_wog(p,u,BB)
!
!  pe <= u1e * u2e
!
  implicit none
  type(vfield_gluon_wog), intent(inout) :: p
  type(vfield_gluon_wg),  intent(in)    :: u,BB
  integer :: mu

  do mu=1,NDIM
    call assign_mult(p%eo(0)%mu(mu),u%eo(0)%mu(mu),BB%eo(0)%mu(mu))
    call assign_mult(p%eo(1)%mu(mu),u%eo(1)%mu(mu),BB%eo(1)%mu(mu))
  enddo

  return
end subroutine
Subroutine :
pe :type(sfield_gluon_eo_wog), intent(inout)
u1e :type(sfield_gluon_eo_wg), intent(in)
u2e :type(sfield_gluon_eo_wg), intent(in)
 pe <= u1e * u2e

[Source]

subroutine assign_mult_sfg_eo_wg_wog(pe,u1e,u2e)
!
!  pe <= u1e * u2e
!
  implicit none
  type(sfield_gluon_eo_wog), intent(inout) :: pe
  type(sfield_gluon_eo_wg),  intent(in)    :: u1e,u2e
  integer :: ix,iy,iz,ieoxyz,itb,ic,jc,itb0

#ifdef _DEBUG_
  if (u1e%ieo /= u2e%ieo) then
    call error_stop("error with u1e%ieo /= u2e%ieo in assign_mult_sfg_eo_wog.")
  endif
#endif
  pe%ieo=u1e%ieo
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,itb0,ic,jc)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+pe%ieo,2)
    do itb=1-ieoxyz,NTH-ieoxyz
      itb0 = itb+ieoxyz
      do jc=1,COL
      do ic=1,COL
        pe%s(itb0,iz,iy,ix)%u(ic,jc)=u1e%s(itb,iz,iy,ix)%u(ic,1) *u2e%s(itb,iz,iy,ix)%u(1,jc) +u1e%s(itb,iz,iy,ix)%u(ic,2) *u2e%s(itb,iz,iy,ix)%u(2,jc) +u1e%s(itb,iz,iy,ix)%u(ic,3) *u2e%s(itb,iz,iy,ix)%u(3,jc)
      enddo
      enddo
    enddo
  enddo
  enddo
  enddo

  return
end subroutine
Subroutine :
u2e :type(sfield_gluon_eo_wg), intent(inout)
u1e :type(sfield_gluon_eo_wg), intent(in)
rcoef :real(DP), intent(in)
 u2e <= u1e * rcoef

[Source]

subroutine assign_rmult_sfg_eo_wg(u2e,u1e,rcoef)
!
!  u2e <= u1e * rcoef
!
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: u2e
  type(sfield_gluon_eo_wg), intent(in)    :: u1e
  real(DP), intent(in) :: rcoef
  integer :: ix,iy,iz,ieoxyz,itb,ic,jc

  u2e%ieo=u1e%ieo
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ic,jc)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+u1e%ieo,2)
    do itb=1-ieoxyz,NTH-ieoxyz
      do jc=1,COL
      do ic=1,COL
        u2e%s(itb,iz,iy,ix)%u(ic,jc)=rcoef*u1e%s(itb,iz,iy,ix)%u(ic,jc)
      enddo
      enddo
    enddo
  enddo
  enddo
  enddo

  return
end subroutine
Subroutine :
u3e :type(sfield_gluon_eo_wg), intent(inout)
u1e :type(sfield_gluon_eo_wg), intent(in)
u2e :type(sfield_gluon_eo_wg), intent(in)
 u3e <= u1e * u2e

[Source]

subroutine assign_mult_sfg_eo_wg(u3e,u1e,u2e)
!
!  u3e <= u1e * u2e
!
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: u3e
  type(sfield_gluon_eo_wg), intent(in)    :: u1e,u2e
  integer :: ix,iy,iz,ieoxyz,itb,ic,jc

#ifdef _DEBUG_
  if (u1e%ieo /= u2e%ieo) then
    call error_stop("error with u1e%ieo /= u2e%ieo in assign_mult_sfg_eo_wg.")
  endif
#endif
  u3e%ieo=u1e%ieo
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ic,jc)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+u3e%ieo,2)
    do itb=1-ieoxyz,NTH-ieoxyz
      do jc=1,COL
      do ic=1,COL
        u3e%s(itb,iz,iy,ix)%u(ic,jc)=u1e%s(itb,iz,iy,ix)%u(ic,1) *u2e%s(itb,iz,iy,ix)%u(1,jc) +u1e%s(itb,iz,iy,ix)%u(ic,2) *u2e%s(itb,iz,iy,ix)%u(2,jc) +u1e%s(itb,iz,iy,ix)%u(ic,3) *u2e%s(itb,iz,iy,ix)%u(3,jc)
      enddo
      enddo
    enddo
  enddo
  enddo
  enddo

  return
end subroutine
Subroutine :
u1o :type(sfield_gluon_eo_wg), intent(inout)
u1e :type(sfield_gluon_eo_wg), intent(in)
mu :integer, intent(in)
 u1o(n) <= u1e(n+mu)

[Source]

subroutine assign_shift1_sfg_eo_wg(u1o,u1e,mu)
!
!  u1o(n) <= u1e(n+mu)
!
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: u1o
  type(sfield_gluon_eo_wg), intent(in)    :: u1e
  integer, intent(in) :: mu
  integer :: ix,iy,iz,ieoxyz,itb
  integer :: ix2,iy2,iz2,itb2,ic,jc

  u1o%ieo=1-u1e%ieo
  select case (mu)
  case (1:NDIM)
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ix2,iy2,iz2,itb2)
    do ix=1,NX
    do iy=1,NY
    do iz=1,NZ
      ieoxyz=mod(ipeo+ix+iy+iz+u1o%ieo,2)
      ix2=ix+isx(mu)
      iy2=iy+isy(mu)
      iz2=iz+isz(mu)
      do itb=1-ieoxyz,NTH-ieoxyz
        itb2=itb+ieoxyz*ist(mu)
        u1o%s(itb,iz,iy,ix)=u1e%s(itb2,iz2,iy2,ix2)
      enddo
    enddo
    enddo
    enddo
  case (-NDIM:-1)
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ix2,iy2,iz2,itb2)
    do ix=1,NX
    do iy=1,NY
    do iz=1,NZ
      ieoxyz=mod(ipeo+ix+iy+iz+u1o%ieo,2)
      ix2=ix-isx(-mu)
      iy2=iy-isy(-mu)
      iz2=iz-isz(-mu)
      do itb=1-ieoxyz,NTH-ieoxyz
        itb2=itb-(1-ieoxyz)*ist(-mu)
        u1o%s(itb,iz,iy,ix)=u1e%s(itb2,iz2,iy2,ix2)
      enddo
    enddo
    enddo
    enddo
  case default
    call error_stop("error with |mu|>NDIM in assign_shift1_sfg_eo_wg.")
  end select

  return
end subroutine
Subroutine :
u2e :type(sfield_gluon_eo_wg), intent(inout)
u1e :type(sfield_gluon_eo_wg), intent(in)
mu :integer, intent(in)
nu :integer, intent(in)

u2e(n) <= u1e(n+mu+nu) (mu /= nu)

[Source]

subroutine assign_shift2_sfg_eo_wg(u2e,u1e,mu,nu)
!
! u2e(n) <= u1e(n+mu+nu)  (mu /= nu)
!
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: u2e
  type(sfield_gluon_eo_wg), intent(in)    :: u1e
  integer, intent(in) :: mu,nu
  integer :: ix,iy,iz,ieoxyz,itb
  integer :: ix2,iy2,iz2,itb2,ic,jc

  if (mu == nu) then
    call error_stop("error with mu == nu in assign_shift2_sfg_eo_wg.")
  endif
  if ((ABS(mu) > NDIM) .OR. (ABS(nu) > NDIM)) then
    call error_stop("error with |mu| > NDIM or |nu| > NDIM in assign_shift2_sfg_eo_wg.")
  endif

  u2e%ieo=u1e%ieo
  if ((mu < 0) .AND. (nu < 0)) then
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ix2,iy2,iz2,itb2)
    do ix=1,NX
    do iy=1,NY
    do iz=1,NZ
      ieoxyz=mod(ipeo+ix+iy+iz+u2e%ieo,2)
      ix2=ix-isx(-mu)-isx(-nu)
      iy2=iy-isy(-mu)-isy(-nu)
      iz2=iz-isz(-mu)-isz(-nu)
      do itb=1-ieoxyz,NTH-ieoxyz
        itb2=itb-(1-ieoxyz)*ist(-mu)-(1-ieoxyz)*ist(-nu)
        u2e%s(itb,iz,iy,ix)=u1e%s(itb2,iz2,iy2,ix2)
      enddo
    enddo
    enddo
    enddo
  else if ((mu > 0) .AND. (nu < 0)) then
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ix2,iy2,iz2,itb2)
    do ix=1,NX
    do iy=1,NY
    do iz=1,NZ
      ieoxyz=mod(ipeo+ix+iy+iz+u2e%ieo,2)
      ix2=ix+isx(mu)-isx(-nu)
      iy2=iy+isy(mu)-isy(-nu)
      iz2=iz+isz(mu)-isz(-nu)
      do itb=1-ieoxyz,NTH-ieoxyz
        itb2=itb+(  ieoxyz)*ist(mu)-(1-ieoxyz)*ist(-nu)
        u2e%s(itb,iz,iy,ix)=u1e%s(itb2,iz2,iy2,ix2)
      enddo
    enddo
    enddo
    enddo
  else if ((mu < 0) .AND. (nu > 0)) then
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ix2,iy2,iz2,itb2)
    do ix=1,NX
    do iy=1,NY
    do iz=1,NZ
      ieoxyz=mod(ipeo+ix+iy+iz+u2e%ieo,2)
      ix2=ix-isx(-mu)+isx(nu)
      iy2=iy-isy(-mu)+isy(nu)
      iz2=iz-isz(-mu)+isz(nu)
      do itb=1-ieoxyz,NTH-ieoxyz
        itb2=itb-(1-ieoxyz)*ist(-mu)+(  ieoxyz)*ist(nu)
        u2e%s(itb,iz,iy,ix)=u1e%s(itb2,iz2,iy2,ix2)
      enddo
    enddo
    enddo
    enddo
  else if ((mu > 0) .AND. (nu > 0)) then
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ix2,iy2,iz2,itb2)
    do ix=1,NX
    do iy=1,NY
    do iz=1,NZ
      ieoxyz=mod(ipeo+ix+iy+iz+u2e%ieo,2)
      ix2=ix+isx(mu)+isx(nu)
      iy2=iy+isy(mu)+isy(nu)
      iz2=iz+isz(mu)+isz(nu)
      do itb=1-ieoxyz,NTH-ieoxyz
        itb2=itb+(  ieoxyz)*ist(mu)+(  ieoxyz)*ist(nu)
        u2e%s(itb,iz,iy,ix)=u1e%s(itb2,iz2,iy2,ix2)
      enddo
    enddo
    enddo
    enddo
  endif

  return
end subroutine
Subroutine :
u3e :type(sfield_gluon_eo_wg), intent(inout)
u1e :type(sfield_gluon_eo_wg), intent(in)
u2e :type(sfield_gluon_eo_wg), intent(in)

u3e = u1e - u2e

[Source]

subroutine assign_sub_sfg_eo_wg(u3e,u1e,u2e)
!
! u3e = u1e - u2e
!
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: u3e
  type(sfield_gluon_eo_wg), intent(in)    :: u1e,u2e
  integer :: ix,iy,iz,ieoxyz,itb,ic,jc

#ifdef _DEBUG_
  if (u1e%ieo /= u2e%ieo) then
    call error_stop("error with u1e%ieo /= u2e%ieo in assign_sub_sfg_eo_wg.")
  endif
#endif
  u3e%ieo=u1e%ieo
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ic,jc)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+u1e%ieo,2)
    do itb=1-ieoxyz,NTH-ieoxyz
      do jc=1,COL
      do ic=1,COL
        u3e%s(itb,iz,iy,ix)%u(ic,jc)=u1e%s(itb,iz,iy,ix)%u(ic,jc) -u2e%s(itb,iz,iy,ix)%u(ic,jc)
      enddo
      enddo
    enddo
  enddo
  enddo
  enddo

  return
end subroutine
Subroutine :
u :type(tfield_gluon_wg), intent(inout)

[Source]

subroutine clear_tfg_wg(u)
  implicit none
  type(tfield_gluon_wg), intent(inout) :: u
  call clear(u%eo(0))
  call clear(u%eo(1))
  return
end subroutine
Subroutine :
u :type(vfield_gluon_wg), intent(inout)

[Source]

subroutine clear_vfg_wg(u)
  implicit none
  type(vfield_gluon_wg), intent(inout) :: u
  call clear(u%eo(0))
  call clear(u%eo(1))
  return
end subroutine
Subroutine :
u :type(vfield_gluon_wog), intent(inout)

[Source]

subroutine clear_vfg_wog(u)
  implicit none
  type(vfield_gluon_wog), intent(inout) :: u
  call clear(u%eo(0))
  call clear(u%eo(1))
  return
end subroutine
Subroutine :
ue :type(sfield_gluon_eo_wg), intent(inout)

[Source]

subroutine clear_sfg_eo_wg(ue)
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: ue
  integer :: ix,iy,iz,ieoxyz,itb,ic,jc

!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ic,jc)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+ue%ieo,2)
    do itb=1-ieoxyz,NTH-ieoxyz
      do jc=1,COL
      do ic=1,COL
        ue%s(itb,iz,iy,ix)%u(ic,jc)=Z0
      enddo
      enddo
    enddo
  enddo
  enddo
  enddo

  return
end subroutine
Subroutine :
ue :type(sfield_gluon_eo_wog), intent(inout)

[Source]

subroutine clear_sfg_eo_wog(ue)
  implicit none
  type(sfield_gluon_eo_wog), intent(inout) :: ue
  integer :: ix,iy,iz,itb

!$OMP PARALLEL DO PRIVATE(ix,iy,iz,itb)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    do itb=1,NTH
      ue%s(itb,iz,iy,ix)%u(:,:)=Z0
    enddo
  enddo
  enddo
  enddo

  return
end subroutine
Subroutine :
ue :type(tfield_gluon_eo_wg), intent(inout)

[Source]

subroutine clear_tfg_eo_wg(ue)
  implicit none
  type(tfield_gluon_eo_wg), intent(inout) :: ue
  integer :: mu
  do mu=1,NDIM*(NDIM-1)/2
    call clear(ue%munu(mu))
  enddo
  return
end subroutine
Subroutine :
ue :type(vfield_gluon_eo_wg), intent(inout)

[Source]

subroutine clear_vfg_eo_wg(ue)
  implicit none
  type(vfield_gluon_eo_wg), intent(inout) :: ue
  integer :: mu
  do mu=1,NDIM
    call clear(ue%mu(mu))
  enddo
  return
end subroutine
Subroutine :
ue :type(vfield_gluon_eo_wog), intent(inout)

[Source]

subroutine clear_vfg_eo_wog(ue)
  implicit none
  type(vfield_gluon_eo_wog), intent(inout) :: ue
  integer :: mu
  do mu=1,NDIM
    call clear(ue%mu(mu))
  enddo
  return
end subroutine
Subroutine :
u :type(sfield_gluon_wg), intent(inout)

[Source]

subroutine copy_sfg(u)
  implicit none
  type(sfield_gluon_wg), intent(inout) :: u
  call copy_boundary(u%eo(0))
  call copy_boundary(u%eo(1))
  return
end subroutine
Subroutine :
u :type(tfield_gluon_wg), intent(inout)

[Source]

subroutine copy_tfg(u)
  implicit none
  type(tfield_gluon_wg), intent(inout) :: u
  call copy_boundary(u%eo(0))
  call copy_boundary(u%eo(1))
  return
end subroutine
Subroutine :
u :type(vfield_gluon_wg), intent(inout)

[Source]

subroutine copy_vfg(u)
  implicit none
  type(vfield_gluon_wg), intent(inout) :: u
  call copy_boundary(u%eo(0))
  call copy_boundary(u%eo(1))
  return
end subroutine
Subroutine :
ue :type(sfield_gluon_eo_wg), intent(inout)

Boundary copy for su(3) scalar field with periodic boundary condition.

 ue : even/odd site su3 scalar field (ieo=0/1)

[Source]

subroutine copy_sfg_eo(ue)
!
! Boundary copy for su(3) scalar field with periodic boundary condition.
!
!  ue : even/odd site su3 scalar field (ieo=0/1)
!
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: ue

  type(fields_gauge) :: fields
  integer :: ieoxyz,itb0,itb1
  integer :: ix,iy,iz,itb
#ifndef _singlePU
  integer :: ibuff,ic,jc
#endif

  if (.not.m_is_initialized) call new_fields_gauge(fields)

  call tic(copy_sfg_time)
#ifndef _singlePU
  call comlib_barrier
#endif

!****** T - boundaray
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb0,itb1)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+ue%ieo,2)
    itb0=   ieoxyz *NTH
    itb1=(1-ieoxyz)*NTH
    ue%s(itb0,iz,iy,ix)=ue%s(itb1,iz,iy,ix)
  enddo
  enddo
  enddo

#ifndef _singlePU
!----------------------------------------------------------
! copy without cornar
! set send buffer
!----------------------------------------------------------
!****** Z - boundaray
#if _NDIMZ == 1
!$OMP PARALLEL DO PRIVATE(ix,iy,itb)
  do ix=1,NX
  do iy=1,NY
    do itb=0,NTH
      ue%s(itb,  0,iy,ix)=ue%s(itb,NZ,iy,ix)
      ue%s(itb,NZ1,iy,ix)=ue%s(itb, 1,iy,ix)
    enddo
  enddo
  enddo
#else
  ibuff=0
  do ix=1,NX
  do iy=1,NY
    do itb=0,NTH
      do jc=1,COL
      do ic=1,COL
        ibuff=ibuff+1
        gbuffup(ibuff,1,3)=ue%s(itb,NZ,iy,ix)%u(ic,jc)
        gbuffdn(ibuff,1,3)=ue%s(itb, 1,iy,ix)%u(ic,jc)
      enddo
      enddo
    enddo
  enddo
  enddo
  if (ibuff.NE.GBUFF_SIZE(3)) then
    call error_stop(" Error : com-buff size in Z")
  endif
#endif

!****** Y - boundaray
#if _NDIMY == 1
!$OMP PARALLEL DO PRIVATE(ix,iz,itb)
  do ix=1,NX
  do iz=1,NZ
    do itb=0,NTH
      ue%s(itb,iz,  0,ix)=ue%s(itb,iz,NY,ix)
      ue%s(itb,iz,NY1,ix)=ue%s(itb,iz, 1,ix)
    enddo
  enddo
  enddo
#else
  ibuff=0
  do ix=1,NX
  do iz=1,NZ
    do itb=0,NTH
      do jc=1,COL
      do ic=1,COL
        ibuff=ibuff+1
        gbuffup(ibuff,1,2)=ue%s(itb,iz,NY,ix)%u(ic,jc)
        gbuffdn(ibuff,1,2)=ue%s(itb,iz, 1,ix)%u(ic,jc)
      enddo
      enddo
    enddo
  enddo
  enddo
  if (ibuff.NE.GBUFF_SIZE(2)) then
    call error_stop(" Error : com-buff size in Y")
  endif
#endif

!****** X - boundaray
#if _NDIMX == 1
!$OMP PARALLEL DO PRIVATE(iy,iz,itb)
  do iy=1,NY
  do iz=1,NZ
    do itb=0,NTH
      ue%s(itb,iz,iy,  0)=ue%s(itb,iz,iy,NX)
      ue%s(itb,iz,iy,NX1)=ue%s(itb,iz,iy, 1)
    enddo
  enddo
  enddo
#else
  ibuff=0
  do iy=1,NY
  do iz=1,NZ
    do itb=0,NTH
      do jc=1,COL
      do ic=1,COL
        ibuff=ibuff+1
        gbuffup(ibuff,1,1)=ue%s(itb,iz,iy,NX)%u(ic,jc)
        gbuffdn(ibuff,1,1)=ue%s(itb,iz,iy, 1)%u(ic,jc)
      enddo
      enddo
    enddo
  enddo
  enddo
  if (ibuff.NE.GBUFF_SIZE(1)) then
    call error_stop(" Error : com-buff size in X")
  endif
#endif

#if _NDIMZ != 1
  call comlib_sendrecv(idgsendup(3))
  call comlib_sendrecv(idgsenddn(3))
#endif
#if _NDIMY != 1
  call comlib_sendrecv(idgsendup(2))
  call comlib_sendrecv(idgsenddn(2))
#endif
#if _NDIMX != 1
  call comlib_sendrecv(idgsendup(1))
  call comlib_sendrecv(idgsenddn(1))
#endif

!****** Z - boundaray
#if _NDIMZ != 1
  ibuff=0
  do ix=1,NX
  do iy=1,NY
    do itb=0,NTH
      do jc=1,COL
      do ic=1,COL
        ibuff=ibuff+1
        ue%s(itb,  0,iy,ix)%u(ic,jc)=gbuffup(ibuff,2,3)
        ue%s(itb,NZ1,iy,ix)%u(ic,jc)=gbuffdn(ibuff,2,3)
      enddo
      enddo
    enddo
  enddo
  enddo
#endif

!****** Y - boundaray
#if _NDIMY != 1
  ibuff=0
  do ix=1,NX
  do iz=1,NZ
    do itb=0,NTH
      do jc=1,COL
      do ic=1,COL
        ibuff=ibuff+1
        ue%s(itb,iz,  0,ix)%u(ic,jc)=gbuffup(ibuff,2,2)
        ue%s(itb,iz,NY1,ix)%u(ic,jc)=gbuffdn(ibuff,2,2)
      enddo
      enddo
    enddo
  enddo
  enddo
#endif

!****** X - boundaray
#if _NDIMX != 1
  ibuff=0
  do iy=1,NY
  do iz=1,NZ
    do itb=0,NTH
      do jc=1,COL
      do ic=1,COL
        ibuff=ibuff+1
        ue%s(itb,iz,iy,  0)%u(ic,jc)=gbuffup(ibuff,2,1)
        ue%s(itb,iz,iy,NX1)%u(ic,jc)=gbuffdn(ibuff,2,1)
      enddo
      enddo
    enddo
  enddo
  enddo
#endif

!--------------------------------------------------------
! copy cornar
!--------------------------------------------------------
  call comlib_barrier

!****** T - boundaray
!****** t-x
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb0,itb1)
  do ix=0,NX1,NX1
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+ue%ieo,2)
    itb0=   ieoxyz *NTH
    itb1=(1-ieoxyz)*NTH
    ue%s(itb0,iz,iy,ix)=ue%s(itb1,iz,iy,ix)
  enddo
  enddo
  enddo
!****** t-y
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb0,itb1)
  do ix=1,NX
  do iy=0,NY1,NY1
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+ue%ieo,2)
    itb0=   ieoxyz *NTH
    itb1=(1-ieoxyz)*NTH
    ue%s(itb0,iz,iy,ix)=ue%s(itb1,iz,iy,ix)
  enddo
  enddo
  enddo
!****** t-z
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb0,itb1)
  do ix=1,NX
  do iy=1,NY
  do iz=0,NZ1,NZ1
    ieoxyz=mod(ipeo+ix+iy+iz+ue%ieo,2)
    itb0=   ieoxyz *NTH
    itb1=(1-ieoxyz)*NTH
    ue%s(itb0,iz,iy,ix)=ue%s(itb1,iz,iy,ix)
  enddo
  enddo
  enddo

!*********************** copy cornar
!*********************** send
!****** x-y
#if _NDIMX == 1
!$OMP PARALLEL DO PRIVATE(iz,itb)
    do iz=1,NZ
      do itb=0,NTH
        ue%s(itb,iz,  0,  0)=ue%s(itb,iz,  0,NX)
        ue%s(itb,iz,NY1,  0)=ue%s(itb,iz,NY1,NX)
        ue%s(itb,iz,  0,NX1)=ue%s(itb,iz,  0, 1)
        ue%s(itb,iz,NY1,NX1)=ue%s(itb,iz,NY1, 1)
      enddo
    enddo
#else
    ibuff=0
    do iz=1,NZ
      do itb=0,NTH
        do jc=1,COL
        do ic=1,COL
          ibuff=ibuff+2
          gbuffup(ibuff-1,1,4)=ue%s(itb,iz,  0,NX)%u(ic,jc)
          gbuffup(ibuff  ,1,4)=ue%s(itb,iz,NY1,NX)%u(ic,jc)
          gbuffdn(ibuff-1,1,4)=ue%s(itb,iz,  0, 1)%u(ic,jc)
          gbuffdn(ibuff  ,1,4)=ue%s(itb,iz,NY1, 1)%u(ic,jc)
        enddo
        enddo
      enddo
    enddo
    if (ibuff.NE.GBUFF_SIZE(4)) then
      call error_stop(" Error : com-buff size in (x-y)")
    endif
#endif
!****** x-z
#if _NDIMX == 1
!$OMP PARALLEL DO PRIVATE(iy,itb)
    do iy=1,NY
      do itb=0,NTH
        ue%s(itb,  0,iy,  0)=ue%s(itb,  0,iy,NX)
        ue%s(itb,NZ1,iy,  0)=ue%s(itb,NZ1,iy,NX)
        ue%s(itb,  0,iy,NX1)=ue%s(itb,  0,iy, 1)
        ue%s(itb,NZ1,iy,NX1)=ue%s(itb,NZ1,iy, 1)
      enddo
    enddo
#else
    ibuff=0
    do iy=1,NY
      do itb=0,NTH
        do jc=1,COL
        do ic=1,COL
          ibuff=ibuff+2
          gbuffup(ibuff-1,1,5)=ue%s(itb,  0,iy,NX)%u(ic,jc)
          gbuffup(ibuff  ,1,5)=ue%s(itb,NZ1,iy,NX)%u(ic,jc)
          gbuffdn(ibuff-1,1,5)=ue%s(itb,  0,iy, 1)%u(ic,jc)
          gbuffdn(ibuff  ,1,5)=ue%s(itb,NZ1,iy, 1)%u(ic,jc)
        enddo
        enddo
      enddo
    enddo
    if (ibuff.NE.GBUFF_SIZE(5)) then
      call error_stop(" Error : com-buff size in (x-z)")
    endif
#endif
!****** y-z
#if _NDIMY == 1
!$OMP PARALLEL DO PRIVATE(ix,itb)
    do ix=1,NX
      do itb=0,NTH
        ue%s(itb,  0,  0,ix)=ue%s(itb,  0,NY,ix)
        ue%s(itb,NZ1,  0,ix)=ue%s(itb,NZ1,NY,ix)
        ue%s(itb,  0,NY1,ix)=ue%s(itb,  0, 1,ix)
        ue%s(itb,NZ1,NY1,ix)=ue%s(itb,NZ1, 1,ix)
      enddo
    enddo
#else
    ibuff=0
    do ix=1,NX
      do itb=0,NTH
        do jc=1,COL
        do ic=1,COL
          ibuff=ibuff+2
          gbuffup(ibuff-1,1,6)=ue%s(itb,  0,NY,ix)%u(ic,jc)
          gbuffup(ibuff  ,1,6)=ue%s(itb,NZ1,NY,ix)%u(ic,jc)
          gbuffdn(ibuff-1,1,6)=ue%s(itb,  0, 1,ix)%u(ic,jc)
          gbuffdn(ibuff  ,1,6)=ue%s(itb,NZ1, 1,ix)%u(ic,jc)
        enddo
        enddo
      enddo
    enddo
    if (ibuff.NE.GBUFF_SIZE(6)) then
      call error_stop(" Error : com-buff size in (y-z)")
    endif
#endif

#if _NDIMX != 1
    call comlib_sendrecv(idgsendup(4))
    call comlib_sendrecv(idgsenddn(4))
    call comlib_sendrecv(idgsendup(5))
    call comlib_sendrecv(idgsenddn(5))
#endif
#if _NDIMY != 1
    call comlib_sendrecv(idgsendup(6))
    call comlib_sendrecv(idgsenddn(6))
#endif

!*********************** recv
!****** x-y
#if _NDIMX != 1
    ibuff=0
    do iz=1,NZ
      do itb=0,NTH
        do jc=1,COL
        do ic=1,COL
          ibuff=ibuff+2
          ue%s(itb,iz,  0,  0)%u(ic,jc)=gbuffup(ibuff-1,2,4)
          ue%s(itb,iz,NY1,  0)%u(ic,jc)=gbuffup(ibuff  ,2,4)
          ue%s(itb,iz,  0,NX1)%u(ic,jc)=gbuffdn(ibuff-1,2,4)
          ue%s(itb,iz,NY1,NX1)%u(ic,jc)=gbuffdn(ibuff  ,2,4)
        enddo
        enddo
      enddo
    enddo
!****** x-z
    ibuff=0
    do iy=1,NY
      do itb=0,NTH
        do jc=1,COL
        do ic=1,COL
          ibuff=ibuff+2
          ue%s(itb,  0,iy,  0)%u(ic,jc)=gbuffup(ibuff-1,2,5)
          ue%s(itb,NZ1,iy,  0)%u(ic,jc)=gbuffup(ibuff  ,2,5)
          ue%s(itb,  0,iy,NX1)%u(ic,jc)=gbuffdn(ibuff-1,2,5)
          ue%s(itb,NZ1,iy,NX1)%u(ic,jc)=gbuffdn(ibuff  ,2,5)
        enddo
        enddo
      enddo
    enddo
#endif
!****** y-z
#if _NDIMY != 1
    ibuff=0
    do ix=1,NX
      do itb=0,NTH
        do jc=1,COL
        do ic=1,COL
          ibuff=ibuff+2
          ue%s(itb,  0,  0,ix)%u(ic,jc)=gbuffup(ibuff-1,2,6)
          ue%s(itb,NZ1,  0,ix)%u(ic,jc)=gbuffup(ibuff  ,2,6)
          ue%s(itb,  0,NY1,ix)%u(ic,jc)=gbuffdn(ibuff-1,2,6)
          ue%s(itb,NZ1,NY1,ix)%u(ic,jc)=gbuffdn(ibuff  ,2,6)
        enddo
        enddo
      enddo
    enddo
#endif

#else

!************************************************************************
!****** Z - boundaray
!$OMP PARALLEL DO  &
!$OMP PRIVATE(ix,iy,itb)
  do ix = 1,NX
  do iy = 1,NY
    do itb = 0,NTH
      ue%s(itb,  0,iy,ix)=ue%s(itb,NZ,iy,ix)
      ue%s(itb,NZ1,iy,ix)=ue%s(itb, 1,iy,ix)
    enddo
  enddo
  enddo

!****** Y - boundaray
!$OMP PARALLEL DO  &
!$OMP PRIVATE(ix,iz,itb)
  do ix = 1,NX
  do iz = 0,NZ1
    do itb = 0,NTH
      ue%s(itb,iz,  0,ix)=ue%s(itb,iz,NY,ix)
      ue%s(itb,iz,NY1,ix)=ue%s(itb,iz, 1,ix)
    enddo
  enddo
  enddo

!****** X - boundaray
!$OMP PARALLEL DO  &
!$OMP COLLAPSE(2)  &
!$OMP PRIVATE(iy,iz,itb)
  do iy = 0,NY1
  do iz = 0,NZ1
    do itb = 0,NTH
      ue%s(itb,iz,iy,  0)=ue%s(itb,iz,iy,NX)
      ue%s(itb,iz,iy,NX1)=ue%s(itb,iz,iy, 1)
    enddo
  enddo
  enddo

#endif

  call toc(copy_sfg_time)

  return
end subroutine
Subroutine :
ue :type(tfield_gluon_eo_wg), intent(inout)

[Source]

subroutine copy_tfg_eo(ue)
  implicit none
  type(tfield_gluon_eo_wg), intent(inout) :: ue
  integer :: munu
  do munu=1,NDIM*(NDIM-1)/2
    call copy_boundary(ue%munu(munu))
  enddo
  return
end subroutine
Subroutine :
ue :type(vfield_gluon_eo_wg), intent(inout)

[Source]

subroutine copy_vfg_eo(ue)
  implicit none
  type(vfield_gluon_eo_wg), intent(inout) :: ue
  integer :: mu
  do mu=1,NDIM
    call copy_boundary(ue%mu(mu))
  enddo
#ifdef _SF
  call set_sf_boundary(ue)
#endif
  return
end subroutine
copy_sfg_time
Variable :
copy_sfg_time :type(timer), save
: holds gauge field boundary copy timing
Function :
norm :real(DP)
F :type(vfield_gluon_wog), intent(in)
: MD force

return averaged force norm

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

[Source]

function get_force_norm_wog(F) result(norm)
!
! return averaged force norm
!
!  norm = SQRT(abs2(F)/(NDIM*SPACETIMEVOLUME))
!
  implicit none
  type(vfield_gluon_wog), intent(in) :: F  ! MD force
  real(DP) :: norm
  integer :: mu

  norm = 0.0_DP
  do mu=1,NDIM
    norm = norm + abs2(F%eo(0)%mu(mu))
    norm = norm + abs2(F%eo(1)%mu(mu))
  enddo
  norm = sqrt(norm/(NDIM*NTT*NTZ*NTY*NTX))

  return
end function
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^+.

[Source]

subroutine make_two_links(u,wl34,wl98)
!
! 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^+.
! 
  implicit none
  type(vfield_gluon_wg), intent(in)    :: u     ! gauge field
  type(tfield_gluon_wg), intent(inout) :: wl34  ! even/odd site upper right 2-links
  type(tfield_gluon_wg), intent(inout) :: wl98  ! even/odd site lower left  2-links

#ifdef _FLAT_
  call make_two_links_eo_flat(u%eo(0),u%eo(1),wl34%eo(0),wl98%eo(0))
  call make_two_links_eo_flat(u%eo(1),u%eo(0),wl34%eo(1),wl98%eo(1))
#else
  call make_two_links_eo(u%eo(0),u%eo(1),wl34%eo(0),wl98%eo(0))
  call make_two_links_eo(u%eo(1),u%eo(0),wl34%eo(1),wl98%eo(1))
#endif

  return
end subroutine
Subroutine :
u :type(sfield_gluon_wg), intent(inout)

set even/odd index on the field

[Source]

subroutine new_sfg_wg(u)
  !
  ! set even/odd index on the field 
  !
  implicit none
  type(sfield_gluon_wg), intent(inout) :: u
  call new(u%eo(0),0)
  call new(u%eo(1),1)
  return
end subroutine
Subroutine :
u :type(sfield_gluon_wog), intent(inout)

set even/odd index on the field

[Source]

subroutine new_sfg_wog(u)
  !
  ! set even/odd index on the field 
  !
  implicit none
  type(sfield_gluon_wog), intent(inout) :: u
  call new(u%eo(0),0)
  call new(u%eo(1),1)
  return
end subroutine
Subroutine :
u :type(tfield_gluon_wg), intent(inout)

set even/odd index on the field

[Source]

subroutine new_tfg_wg(u)
  !
  ! set even/odd index on the field 
  !
  implicit none
  type(tfield_gluon_wg), intent(inout) :: u
  call new(u%eo(0),0)
  call new(u%eo(1),1)
  return
end subroutine
Subroutine :
u :type(vfield_gluon_wg), intent(inout)

set even/odd index on the field

[Source]

subroutine new_vfg_wg(u)
  !
  ! set even/odd index on the field 
  !
  implicit none
  type(vfield_gluon_wg), intent(inout) :: u
  call new(u%eo(0),0)
  call new(u%eo(1),1)
  return
end subroutine
Subroutine :
u :type(vfield_gluon_wog), intent(inout)

set even/odd index on the field

[Source]

subroutine new_vfg_wog(u)
  !
  ! set even/odd index on the field 
  !
  implicit none
  type(vfield_gluon_wog), intent(inout) :: u
  call new(u%eo(0),0)
  call new(u%eo(1),1)
  return
end subroutine
Subroutine :
ue :type(sfield_gluon_eo_wg), intent(inout)
ieo :integer, intent(in)

set even/odd index on the field

[Source]

subroutine new_sfg_eo_wg(ue,ieo)
!
! set even/odd index on the field 
!
  implicit none
  type(sfield_gluon_eo_wg), intent(inout) :: ue
  integer, intent(in) :: ieo
  ue%ieo=ieo
  return
end subroutine
Subroutine :
ue :type(sfield_gluon_eo_wog), intent(inout)
ieo :integer, intent(in)

set even/odd index on the field

[Source]

subroutine new_sfg_eo_wog(ue,ieo)
!
! set even/odd index on the field 
!
  implicit none
  type(sfield_gluon_eo_wog), intent(inout) :: ue
  integer, intent(in) :: ieo
  ue%ieo=ieo
  return
end subroutine
Subroutine :
ue :type(tfield_gluon_eo_wg), intent(inout)
ieo :integer, intent(in)

set even/odd index on the field

[Source]

subroutine new_tfg_eo_wg(ue,ieo)
  !
  ! set even/odd index on the field 
  !
  implicit none
  type(tfield_gluon_eo_wg), intent(inout) :: ue
  integer, intent(in) :: ieo
  integer :: mu
  ue%ieo=ieo
  do mu=1,NDIM*(NDIM-1)/2
    call new(ue%munu(mu),ieo)
  enddo
  return
end subroutine
Subroutine :
ue :type(vfield_gluon_eo_wg), intent(inout)
ieo :integer, intent(in)

set even/odd index on the field

[Source]

subroutine new_vfg_eo_wg(ue,ieo)
  !
  ! set even/odd index on the field 
  !
  implicit none
  type(vfield_gluon_eo_wg), intent(inout) :: ue
  integer, intent(in) :: ieo
  integer :: mu
  ue%ieo=ieo
  do mu=1,NDIM
    call new(ue%mu(mu),ieo)
  enddo
  return
end subroutine
Subroutine :
ue :type(vfield_gluon_eo_wog), intent(inout)
ieo :integer, intent(in)

set even/odd index on the field

[Source]

subroutine new_vfg_eo_wog(ue,ieo)
  !
  ! set even/odd index on the field 
  !
  implicit none
  type(vfield_gluon_eo_wog), intent(inout) :: ue
  integer, intent(in) :: ieo
  integer :: mu
  ue%ieo=ieo
  do mu=1,NDIM
    call new(ue%mu(mu),ieo)
  enddo
  return
end subroutine
Subroutine :
u :type(vfield_gluon_wg), intent(in)
: gauge field
plq :real(DP), intent(out)
: averaged plaquett

Calculate space-time-direction averaged plaquette

[Source]

subroutine plaquette_u(u,plq)
!
! Calculate space-time-direction averaged plaquette
!
  use comlib
  implicit none
  type(vfield_gluon_wg), intent(in)  :: u    ! gauge field
  real(DP),              intent(out) :: plq  ! averaged plaquett
  real(DP) :: plqe,plqo

  call plaquette_u_eo(u%eo(0),u%eo(1),plqe)
  call plaquette_u_eo(u%eo(1),u%eo(0),plqo)

  plq=plqe+plqo

  return
end subroutine
Subroutine :

[Source]

subroutine print_copy_sfg_statistics()
  implicit none
  integer :: mu
  if (.not.m_is_initialized) return
  if (0 == nodeid) then
    do mu=1,NDIM*(NDIM-1)/2
      call comlib_print_statistics(idgsendup(mu))
      call comlib_print_statistics(idgsenddn(mu))
    enddo
  endif
  return
end subroutine
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)

[Source]

subroutine set_gaussian_noise_p(p)
!
! \set Gaussian noise on canonical momentum (su(3) Lie algebra) of gauge field (SU(3) Lie group)
!
  implicit none
  type(vfield_gluon_wog), intent(inout) :: p   !  gauge momentum field

  integer :: COLADJ
  parameter (COLADJ=COL**2-1)
  real(DP) :: He(COLADJ,0:NTH),Ho(COLADJ,0:NTH)
  complex(DP) :: SU(COL,COL,COLADJ)
  real(DP) :: pi2
  real(DP) :: th(COLADJ*NTH),yr(COLADJ*NTH),sth(COLADJ*NTH)
  integer :: ix,iy,iz,itb,ieoxyz,ic,jc,mu,icad,itb0
  integer :: ictb,ir
#ifdef _PARALLELCHECK
  integer :: itx,ity,itz,itt,ipx,ipy,ipz,ieo
  real(DP) :: HH(COLADJ)
  real(DP) :: pth(COLADJ),pyr(COLADJ)
#endif

  pi2 = PI*2
  call new(p)

  do icad = 1,COLADJ
    do jc = 1,COL
    do ic = 1,COL
      SU(ic,jc,icad) = z0
    enddo
    enddo
  enddo

!****************************
!* set su(3) generator 
!* (times sqrt(2))
!****************************
  SU(1,2,1) =   z1/sqrt(2.0_DP)
  SU(2,1,1) =   z1/sqrt(2.0_DP)
  SU(1,2,2) =  -zi/sqrt(2.0_DP)
  SU(2,1,2) =   zi/sqrt(2.0_DP)
  SU(1,1,3) =   z1/sqrt(2.0_DP)
  SU(2,2,3) =  -z1/sqrt(2.0_DP)
  SU(1,3,4) =   z1/sqrt(2.0_DP)
  SU(3,1,4) =   z1/sqrt(2.0_DP)
  SU(1,3,5) =  -zi/sqrt(2.0_DP)
  SU(3,1,5) =   zi/sqrt(2.0_DP)
  SU(2,3,6) =   z1/sqrt(2.0_DP)
  SU(3,2,6) =   z1/sqrt(2.0_DP)
  SU(2,3,7) =  -zi/sqrt(2.0_DP)
  SU(3,2,7) =   zi/sqrt(2.0_DP)
  SU(1,1,8) =   z1/sqrt(6.0_DP)
  SU(2,2,8) =   z1/sqrt(6.0_DP)
  SU(3,3,8) =-2*z1/sqrt(6.0_DP)

#ifdef _PARALLELCHECK
!***********************************
!* programmed independent on 
!* parallelization size
!***********************************
  do mu=1,NDIM
  do itx=1,NTX
  do ity=1,NTY
  do itz=1,NTZ
  do itt=1,NTT
    ieo = mod(itx+ity+itz+itt,2)
    ipx = (itx-1)/NX
    ipy = (ity-1)/NY
    ipz = (itz-1)/NZ
    ix  = mod(itx-1,NX)+1
    iy  = mod(ity-1,NY)+1
    iz  = mod(itz-1,NZ)+1
    itb = itt/2

    call get(g_rand,pyr)
    call get(g_rand,pth)

    if ( (ipx == ipsite(1)) .and. (ipy == ipsite(2)) .and. (ipz == ipsite(3)) ) then

      ieoxyz = mod(ipeo+iz+iy+ix+ieo,2)
      itb0 = itb + ieoxyz

      do icad=1,COLADJ
        pth(icad) = pth(icad)*pi2
        pyr(icad) = SQRT(-2.0_DP*LOG(pyr(icad)))
        pth(icad) =              SIN(pth(icad))
         HH(icad) = pyr(icad)*pth(icad)
      enddo
      do jc = 1,COL
      do ic = 1,COL
        p%eo(ieo)%mu(mu)%s(itb0,iz,iy,ix)%u(ic,jc)=HH(1)*SU(ic,jc,1) +HH(2)*SU(ic,jc,2) +HH(3)*SU(ic,jc,3) +HH(4)*SU(ic,jc,4) +HH(5)*SU(ic,jc,5) +HH(6)*SU(ic,jc,6) +HH(7)*SU(ic,jc,7) +HH(8)*SU(ic,jc,8)
      enddo
      enddo
    endif
  enddo
  enddo
  enddo
  enddo
  enddo
#else

!********************************
!* set initial momentum
!* with parallel rundom number
!* generator (M series).
!*
!* Outer (ix,iy,iz) loop is 
!* divided by NPULOCAL to achieve
!* better efficiency for element
!* parallelization.
!********************************
  do mu = 1,NDIM
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz,2)

!********************
!* even sites
!********************
      call get(g_rand,yr)
      call get(g_rand,th)
      do ictb=1,COLADJ*NTH
         yr(ictb) = SQRT(-2.0_DP*LOG(yr(ictb)))
        sth(ictb) = SIN(th(ictb)*pi2)
      enddo
      ir=1
      do itb=1-ieoxyz,NTH-ieoxyz
        do ic = 1,COLADJ
          He(ic,itb) = yr(ir)*sth(ir)
          ir=ir+1
        enddo
      enddo

!********************
!* odd sites
!********************
      call get(g_rand,yr)
      call get(g_rand,th)
      do ictb=1,COLADJ*NTH
         yr(ictb) = SQRT(-2.0_DP*LOG(yr(ictb)))
        sth(ictb) = SIN(th(ictb)*pi2)
      enddo
      ir=1
      do itb=ieoxyz,NTH+ieoxyz-1
        do ic = 1,COLADJ
          Ho(ic,itb) = yr(ir)*sth(ir)
          ir=ir+1
        enddo
      enddo

!********************
!* even sites
!********************
      do itb=1-ieoxyz,NTH-ieoxyz
        itb0 = itb + ieoxyz
        do jc = 1,COL
        do ic = 1,COL
          p%eo(0)%mu(mu)%s(itb0,iz,iy,ix)%u(ic,jc)=He(1,itb)*SU(ic,jc,1) +He(2,itb)*SU(ic,jc,2) +He(3,itb)*SU(ic,jc,3) +He(4,itb)*SU(ic,jc,4) +He(5,itb)*SU(ic,jc,5) +He(6,itb)*SU(ic,jc,6) +He(7,itb)*SU(ic,jc,7) +He(8,itb)*SU(ic,jc,8)
        enddo
        enddo
      enddo

!********************
!* odd sites
!********************
      do itb=ieoxyz,NTH+ieoxyz-1
        itb0 = 1+itb - ieoxyz
        do jc = 1,COL
        do ic = 1,COL
          p%eo(1)%mu(mu)%s(itb0,iz,iy,ix)%u(ic,jc)=Ho(1,itb)*SU(ic,jc,1) +Ho(2,itb)*SU(ic,jc,2) +Ho(3,itb)*SU(ic,jc,3) +Ho(4,itb)*SU(ic,jc,4) +Ho(5,itb)*SU(ic,jc,5) +Ho(6,itb)*SU(ic,jc,6) +Ho(7,itb)*SU(ic,jc,7) +Ho(8,itb)*SU(ic,jc,8)
        enddo
        enddo
      enddo
  enddo
  enddo
  enddo ! end of do ix,iy,iz
  enddo ! end of do mu
#endif

#ifdef _SF
  call set_sf_boundary(p)
#endif

  return
end subroutine
Subroutine :
u :type(vfield_gluon_wg), intent(inout)
: gauge field

Initialize gauge field to unit matrix (cold config)

[Source]

subroutine set_identity_u(u)
!
! Initialize gauge field to unit matrix (cold config)
!
  implicit none
  type(vfield_gluon_wg), intent(inout) :: u ! gauge field

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

!************************************************************
!* Set initial gauge link
!************************************************************
  call new(u)
  do mu = 1,NDIM
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,itb,ic,jc)
  do ix = 0,NX1
  do iy = 0,NY1
  do iz = 0,NZ1
  do itb = 0,NTH
    do jc = 1,COL
    do ic = 1,COL
      u%eo(0)%mu(mu)%s(itb,iz,iy,ix)%u(ic,jc) = Z0
      u%eo(1)%mu(mu)%s(itb,iz,iy,ix)%u(ic,jc) = Z0
    enddo
    enddo
    do ic = 1,COL
      u%eo(0)%mu(mu)%s(itb,iz,iy,ix)%u(ic,ic) = Z1
      u%eo(1)%mu(mu)%s(itb,iz,iy,ix)%u(ic,ic) = Z1
    enddo
  enddo
  enddo
  enddo
  enddo
  enddo

  call copy_boundary(u)
#ifdef _SF
  call set_sf_boundary(u)
#endif

  return
end subroutine
Subroutine :
p :type(vfield_gluon_wog), intent(inout)

Set Schroedinger Functional boundary on gauge link momentum (spatial). (set zero at itt=1 for spatial link momentum)

 p : gauge link momentum

[Source]

subroutine set_sf_boundary_p(p)
!
! Set Schroedinger Functional boundary on gauge link momentum (spatial).
! (set zero at itt=1 for spatial link momentum)
!
!  p : gauge link momentum
!
  implicit none
  type(vfield_gluon_wog),  intent(inout) :: p
  call set_sf_boundary(p%eo(0))
  call set_sf_boundary(p%eo(1))
  return
end subroutine
Subroutine :
pe :type(vfield_gluon_eo_wog), intent(inout)

Set Schroedinger Functional boundary on gauge link momentum (spatial). (set zero at itt=1 for spatial link momentum)

 pe : gauge link momentum (even/odd sites)

[Source]

subroutine set_sf_boundary_p_eo(pe)
!
! Set Schroedinger Functional boundary on gauge link momentum (spatial).
! (set zero at itt=1 for spatial link momentum)
!
!  pe : gauge link momentum (even/odd sites)
!
  implicit none
  type(vfield_gluon_eo_wog), intent(inout) :: pe
  integer :: ieo,mu,ix,iy,iz,ieoxyz
  ieo = pe%ieo
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,mu)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+ieo,2)
!========================================================
! SF boundary memo (wo ghost sites)
! do itb=1-ieoxyz,NTH-ieoxyz
! itb0 = itb + ieoxyz
! ieoxyz = 0 => itb = 1..NTH => it = 2,4,6,..,NTT
!               itb = 0   => it = 0  (outside of physical sites)
!
! ieoxyz = 1 => itb = 0,..NTH-1 => it = 1,3,5,...,NTT-1
!               itb = NTH => it = NTT+1 (outside of physical sites)
!               itb = 0   => it = 1 => itb0 = 1 
!========================================================
    if (1 == ieoxyz) then
      !==========
      ! t = 1
      !==========
      do mu=1,NDIM-1 ! spatial only
        pe%mu(mu)%s(  1,iz,iy,ix)%u(:,:) = Z0
      enddo
    endif
  enddo
  enddo
  enddo  ! end of do ix,iy,iz

  return
end subroutine
Subroutine :
u :type(vfield_gluon_wg), intent(inout)
conjugate :logical, optional, intent(in)

Set Schroedinger Functional boundary on gauge link.

     this : sf boundary condition
        u : link field

conjugate : option : if true, set conjugate boundary condition

u_k(t=1 ,x,y,z) = U0_k <= constant SU(3) u_k(t=NTT+1,x,y,z) = U1_k <= constant SU(3)

data on itt=0 ghost sites are not touched.

[Source]

subroutine set_sf_boundary_u(u,conjugate)
!
! Set Schroedinger Functional boundary on gauge link.
!
!      this : sf boundary condition
!         u : link field
! conjugate : option : if true, set conjugate boundary condition
!
! u_k(t=1    ,x,y,z) = U0_k <= constant SU(3)
! u_k(t=NTT+1,x,y,z) = U1_k <= constant SU(3)
!
! data on itt=0 ghost sites are not touched.
!
  implicit none
  type(vfield_gluon_wg),    intent(inout) :: u
  logical, optional,        intent(in) :: conjugate
  logical :: is_conjg
  if (.not.present(conjugate)) then
    is_conjg = .false.
  else
    is_conjg = conjugate
  endif
  call set_sf_boundary(u%eo(0),is_conjg)
  call set_sf_boundary(u%eo(1),is_conjg)
  return
end subroutine
Subroutine :
ue :type(vfield_gluon_eo_wg), intent(inout)
conjugate :logical, optional, intent(in)

Set Schroedinger Functional boundary on gauge link.

       ue : link field (even/odd sites)

conjugate : option : if true, set conjugate boundary condition

u_k(t=1 ,x,y,z) = U0_k <= constant SU(3) u_k(t=NTT+1,x,y,z) = U1_k <= constant SU(3)

data on itt=0 ghost sites are not touched.

[Source]

subroutine set_sf_boundary_u_eo(ue,conjugate)
!
! Set Schroedinger Functional boundary on gauge link.
!
!        ue : link field (even/odd sites)
! conjugate : option : if true, set conjugate boundary condition
!
! u_k(t=1    ,x,y,z) = U0_k <= constant SU(3)
! u_k(t=NTT+1,x,y,z) = U1_k <= constant SU(3)
!
! data on itt=0 ghost sites are not touched.
!
  use sf_lattice_class
  implicit none
  type(vfield_gluon_eo_wg), intent(inout) :: ue
  logical, optional,        intent(in) :: conjugate
  logical :: is_conjg
  type(sf_lattice_world) :: sf_lattice
  complex(DP) :: phase0(COL,NDIM-1),phase1(COL,NDIM-1)
  integer :: mu,ieo,ic,jc,ix,iy,iz,ieoxyz

  if (.not.is_initialized(sf_lattice)) call new(sf_lattice)

  if (.not.present(conjugate)) then
    is_conjg = .false.
  else
    is_conjg = conjugate
  endif

  phase0(:,:) = get_phase0(sf_lattice)
  phase1(:,:) = get_phase1(sf_lattice)
  if (is_conjg) then
    phase0(:,:) = conjg(phase0(:,:))
    phase1(:,:) = conjg(phase1(:,:))
  endif

  ieo = ue%ieo

!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,mu,ic)
  do ix=0,NX1
  do iy=0,NY1
  do iz=0,NZ1
    ieoxyz=mod(ipeo+ix+iy+iz+ieo,2)
    if (1 == ieoxyz) then
      do mu=1,NDIM-1 ! spatial only
        !============
        ! t = 1
        !============
        ue%mu(mu)%s(  0,iz,iy,ix)%u(:,:) = Z0
        do ic=1,COL
          ue%mu(mu)%s(  0,iz,iy,ix)%u(ic,ic) = phase0(ic,mu)
        enddo
        !============
        ! t = NTT+1
        !============
        ue%mu(mu)%s(NTH,iz,iy,ix)%u(:,:) = Z0
        do ic=1,COL
          ue%mu(mu)%s(NTH,iz,iy,ix)%u(ic,ic) = phase1(ic,mu)
        enddo
      enddo
    endif
  enddo
  enddo
  enddo  ! end of do ix,iy,iz

  return
end subroutine
sfield_gluon_eo_wg
Derived Type :
s(0:NTH,0:NZ1,0:NY1,0:NX1) :type(su3fm)
ieo :integer
idummy(3) :integer

su(3) matrix scalar field on even/odd sites with ghost sits used for gauge vector/tensor fields.

sfield_gluon_eo_wog
Derived Type :
s(NTH,NZ,NY,NX) :type(su3fm)
ieo :integer
idummy(3) :integer

su(3) matrix scalar field on even/odd sites without ghost sits. used for gauge mementum vector field.

sfield_gluon_wg
Derived Type :
eo(0:1) :type(sfield_gluon_eo_wg)

su(3) matrix scalar field with ghost sits

sfield_gluon_wog
Derived Type :
eo(0:1) :type(sfield_gluon_eo_wog)

su(3) matrix scalar field without ghost sits

su3fm
Derived Type :
u(COL,COL) :complex(DP)

su(3) matrix

tfield_gluon_eo_wg
Derived Type :
munu(NDIM*(NDIM-1)/2) :type(sfield_gluon_eo_wg)
ieo :integer
idummy(3) :integer

su(3) matrix symmetric tensor field on even/odd sites with ghost sits

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

su(3) matrix symmetric tensor field with ghost sits

Function :
my_trace :complex(DP)
ue :type(sfield_gluon_eo_wg), intent(in)

return trace of gauge fields (even/odd sites only)

[Source]

function trace_sfg_eo_wg(ue) result(my_trace)
!
! return trace of gauge fields (even/odd sites only)
!
  use comlib
  implicit none
  type(sfield_gluon_eo_wg), intent(in) :: ue
  complex(DP) :: my_trace
  integer :: ix,iy,iz,ieoxyz,itb,ic

  my_trace = Z0
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb) REDUCTION(+:my_trace)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    ieoxyz=mod(ipeo+ix+iy+iz+ue%ieo,2)
    do itb=1-ieoxyz,NTH-ieoxyz
      do ic=1,COL
        my_trace = my_trace+ue%s(itb,iz,iy,ix)%u(ic,ic)
      enddo
    enddo
  enddo
  enddo
  enddo
#ifndef _singlePU
  call comlib_sumcast(my_trace)
#endif

  return
end function
Subroutine :
p :type(vfield_gluon_wog), intent(inout)
 p <= (p - p^+) - tr[p - p^+]/3

[Source]

subroutine traceless_ahermite_vfg_wog(p)
!
!  p <= (p - p^+) - tr[p - p^+]/3
!
  implicit none
  type(vfield_gluon_wog), intent(inout) :: p
  integer :: mu
  do mu=1,NDIM
    call traceless_antihermite(p%eo(0)%mu(mu))
    call traceless_antihermite(p%eo(1)%mu(mu))
  enddo
  return
end subroutine
Subroutine :
pe :type(sfield_gluon_eo_wog), intent(inout)
 pe <= (pe - pe^+) - tr[pe - pe^+]/3

[Source]

subroutine traceless_ahermite_sfg_eo_wog(pe)
!
!  pe <= (pe - pe^+) - tr[pe - pe^+]/3
!
  implicit none
  type(sfield_gluon_eo_wog), intent(inout) :: pe
  complex(DP) :: trace,ah(COL,COL)
  integer :: ix,iy,iz,itb,ic,jc

!$OMP PARALLEL DO PRIVATE(ix,iy,iz,itb,ic,jc,trace,ah)
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
    do itb=1,NTH
      trace = ( pe%s(itb,iz,iy,ix)%u(1,1) +pe%s(itb,iz,iy,ix)%u(2,2) +pe%s(itb,iz,iy,ix)%u(3,3))/COL

      do ic=1,COL
        pe%s(itb,iz,iy,ix)%u(ic,ic)=pe%s(itb,iz,iy,ix)%u(ic,ic) - trace
      enddo

      do jc=1,COL
      do ic=1,COL
        ah(ic,jc) = pe%s(itb,iz,iy,ix)%u(ic,jc) -conjg(pe%s(itb,iz,iy,ix)%u(jc,ic))
      enddo
      enddo

      do jc=1,COL
      do ic=1,COL
        pe%s(itb,iz,iy,ix)%u(ic,jc) = ah(ic,jc)
      enddo
      enddo

    enddo
  enddo
  enddo
  enddo

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

[Source]

subroutine update_p_wog(p,dt,F)
!
! 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!)
!
  implicit none
  type(vfield_gluon_wog), intent(inout) :: p    ! gauge momentum field
  real(DP),               intent(in)    :: dt   ! MD time step
  type(vfield_gluon_wog), intent(in)    :: F    ! gauge force field

  call update_p_eo_wog(p%eo(0),dt,F%eo(0))
  call update_p_eo_wog(p%eo(1),dt,F%eo(1))

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

[Source]

subroutine update_u(u,dt,p)
!
! Update gauge field using canonical momentum :
!
!    u_mu(dt)=exp(i p_mu(dt')dt) u_mu(0)
!
  implicit none
  type(vfield_gluon_wg),  intent(inout) :: u   ! gauge link field
  real(DP),               intent(in)    :: dt  ! MD step size
  type(vfield_gluon_wog), intent(in)    :: p   ! gauge momentum field

  call update_u_eo(u%eo(0),dt,p%eo(0))
  call update_u_eo(u%eo(1),dt,p%eo(1))

  return
end subroutine
vfield_gluon_eo_wg
Derived Type :
mu(NDIM) :type(sfield_gluon_eo_wg)
ieo :integer
idummy(3) :integer

su(3) matrix vector field on even/odd sites with ghost sits

vfield_gluon_eo_wog
Derived Type :
mu(NDIM) :type(sfield_gluon_eo_wog)
ieo :integer
idummy(3) :integer

su(3) matrix vector field on even/odd sites without ghost sits

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

su(3) matrix vector field with ghost sits

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

su(3) matrix vector field without ghost sits