Class | field_gauge_class |
In: |
FieldGaugeClass/field_gauge_class.F90
|
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.
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.
etc...
Function : | |
my_abs2 : | real(DP) |
pe : | type(sfield_gluon_eo_wog), intent(in) |
return squared norm of gauge field (momentum)
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
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
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)
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)
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
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
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
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
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
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)^+
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))^+
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)
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
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
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
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
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)
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)
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
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) |
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) |
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) |
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) |
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) |
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) |
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) |
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) |
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) |
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) |
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) |
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)
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) |
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) |
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
Function : | |||
norm : | real(DP) | ||
F : | type(vfield_gluon_wog),
intent(in)
|
return averaged force norm
norm = SQRT(abs2(F)/(NDIM*SPACETIMEVOLUME))
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)
| ||
wl34 : | type(tfield_gluon_wg),
intent(inout)
| ||
wl98 : | type(tfield_gluon_wg),
intent(inout)
|
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^+.
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
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
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
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
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
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
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
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
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
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
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)
| ||
plq : | real(DP), intent(out)
|
Calculate space-time-direction averaged plaquette
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 : |
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)
|
set Gaussian noise on canonical momentum (su(3) Lie algebra) of gauge field (SU(3) Lie group)
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)
|
Initialize gauge field to unit matrix (cold config)
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
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)
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.
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.
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
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.
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.
Derived Type : | |
eo(0:1) : | type(sfield_gluon_eo_wg) |
su(3) matrix scalar field with ghost sits
Derived Type : | |
eo(0:1) : | type(sfield_gluon_eo_wog) |
su(3) matrix scalar field without ghost sits
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
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)
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
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
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)
| ||
dt : | real(DP), intent(in)
| ||
F : | type(vfield_gluon_wog),
intent(in)
|
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!)
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)
| ||
dt : | real(DP), intent(in)
| ||
p : | type(vfield_gluon_wog),
intent(in)
|
Update gauge field using canonical momentum :
u_mu(dt)=exp(i p_mu(dt')dt) u_mu(0)
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
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
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
Derived Type : | |
eo(0:1) : | type(vfield_gluon_eo_wg) |
su(3) matrix vector field with ghost sits
Derived Type : | |
eo(0:1) : | type(vfield_gluon_eo_wog) |
su(3) matrix vector field without ghost sits