Class | field_fermion_class |
In: |
FieldFermionClass/field_fermion_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.
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) |
q : | type(field_quark_wg), intent(in) |
return |q|^2
function abs2_fq_wg( q ) result(my_abs2) ! ! return |q|^2 ! implicit none type(field_quark_wg), intent(in) :: q real(DP) :: my_abs2 my_abs2 = abs2(q%eo(0)) + abs2(q%eo(1)) return end function
Function : | |
my_abs2 : | real(DP) |
qe : | type(field_quark_eo_wg), intent(in) |
return |qe|^2
function abs2_fq_eo_wg( qe ) result(my_abs2) ! ! return |qe|^2 ! implicit none type(field_quark_eo_wg), intent(in) :: qe real(DP) :: my_abs2 integer :: ix,iy,iz,itb,ic,is,ieoxyz my_abs2=0.0_DP !$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ic,is) REDUCTION(+:my_abs2) do ix=1,NX do iy=1,NY do iz=1,NZ ieoxyz=mod(ipeo+ix+iy+iz+qe%ieo,2) #ifdef _SF do itb=1,NTH-ieoxyz #else do itb=1-ieoxyz,NTH-ieoxyz #endif do is=1,SPIN do ic=1,COL my_abs2 = my_abs2 + real(qe%s(itb,iz,iy,ix)%y(ic,is))**2 +aimag(qe%s(itb,iz,iy,ix)%y(ic,is))**2 enddo enddo enddo enddo enddo enddo #ifndef _singlePU call comlib_sumcast(my_abs2) #endif return end function
Subroutine : | |
q2 : | type(field_quark_wg), intent(inout) |
q1 : | type(field_quark_wg), intent(in) |
q2 <= q2 + q1
subroutine accum_add_fq_wg(q2,q1) ! ! q2 <= q2 + q1 ! implicit none type(field_quark_wg), intent(inout) :: q2 type(field_quark_wg), intent(in) :: q1 call accum_add(q2%eo(0),q1%eo(0)) call accum_add(q2%eo(1),q1%eo(1)) return end subroutine
Subroutine : | |
q2e : | type(field_quark_eo_wg), intent(inout) |
q1e : | type(field_quark_eo_wg), intent(in) |
q2e <= q2e + q1e
subroutine accum_add_fq_eo_wg(q2e,q1e) ! ! q2e <= q2e + q1e ! implicit none type(field_quark_eo_wg), intent(inout) :: q2e type(field_quark_eo_wg), intent(in) :: q1e integer :: ix,iy,iz,ieoxyz,itb,ic,is #ifdef _DEBUG_ if (q1e%ieo /= q2e%ieo) then if (nodeid == 0) then write(*,'("error in wqf_eo_wg_accum_add(even/odd)")') endif stop endif #endif !$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+q2e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q2e%s(itb,iz,iy,ix)%y(:,:)=q2e%s(itb,iz,iy,ix)%y(:,:) +q1e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
q2 : | type(field_quark_wg), intent(inout) |
q1 : | type(field_quark_wg), intent(in) |
coef : | complex(DP), intent(in) |
q2 <= q2 + q1 * coef
subroutine accum_add_cmult_fq_wg(q2,q1,coef) ! ! q2 <= q2 + q1 * coef ! implicit none type(field_quark_wg), intent(inout) :: q2 type(field_quark_wg), intent(in) :: q1 complex(DP), intent(in) :: coef call accum_add_mult(q2%eo(0),q1%eo(0),coef) call accum_add_mult(q2%eo(1),q1%eo(1),coef) return end subroutine
Subroutine : | |
q2 : | type(field_quark_wg), intent(inout) |
q1 : | type(field_quark_wg), intent(in) |
coef : | real(DP), intent(in) |
q2 <= q2 + q1 * coef
subroutine accum_add_rmult_fq_wg(q2,q1,coef) ! ! q2 <= q2 + q1 * coef ! implicit none type(field_quark_wg), intent(inout) :: q2 type(field_quark_wg), intent(in) :: q1 real(DP), intent(in) :: coef call accum_add_mult(q2%eo(0),q1%eo(0),coef) call accum_add_mult(q2%eo(1),q1%eo(1),coef) return end subroutine
Subroutine : | |
q2e : | type(field_quark_eo_wg), intent(inout) |
q1e : | type(field_quark_eo_wg), intent(in) |
coef : | complex(DP), intent(in) |
q2e <= q2e + q1e * coef
subroutine accum_add_cmult_fq_eo_wg(q2e,q1e,coef) ! ! q2e <= q2e + q1e * coef ! implicit none type(field_quark_eo_wg), intent(inout) :: q2e type(field_quark_eo_wg), intent(in) :: q1e complex(DP), intent(in) :: coef integer :: ix,iy,iz,ieoxyz,itb,ic,is #ifdef _DEBUG_ if (q1e%ieo /= q2e%ieo) then if (nodeid == 0) then write(*,'("error in wqf_eo_wg_cmult1_accum_add(even/odd)")') endif stop endif #endif !$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+q2e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q2e%s(itb,iz,iy,ix)%y(:,:)=q2e%s(itb,iz,iy,ix)%y(:,:) +coef*q1e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
q2e : | type(field_quark_eo_wg), intent(inout) |
q1e : | type(field_quark_eo_wg), intent(in) |
coef : | real(DP), intent(in) |
q2e <= q2e + q1e * coef
subroutine accum_add_rmult_fq_eo_wg(q2e,q1e,coef) ! ! q2e <= q2e + q1e * coef ! implicit none type(field_quark_eo_wg), intent(inout) :: q2e type(field_quark_eo_wg), intent(in) :: q1e real(DP), intent(in) :: coef integer :: ix,iy,iz,ieoxyz,itb,ic,is #ifdef _DEBUG_ if (q1e%ieo /= q2e%ieo) then if (nodeid == 0) then write(*,'("error in wqf_eo_wg_rmult1_accum_add(even/odd)")') endif stop endif #endif !$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+q2e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q2e%s(itb,iz,iy,ix)%y(:,:)=q2e%s(itb,iz,iy,ix)%y(:,:) +coef*q1e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
q1 : | type(field_quark_wg), intent(inout) |
coef : | complex(DP), intent(in) |
q1 <= q1 * coef
subroutine accum_cmult_fq_wg(q1,coef) ! ! q1 <= q1 * coef ! implicit none type(field_quark_wg), intent(inout) :: q1 complex(DP), intent(in) :: coef call accum_mult(q1%eo(0),coef) call accum_mult(q1%eo(1),coef) return end subroutine
Subroutine : | |
q1 : | type(field_quark_wg), intent(inout) |
coef : | real(DP), intent(in) |
q1 <= q1 * coef
subroutine accum_rmult_fq_wg(q1,coef) ! ! q1 <= q1 * coef ! implicit none type(field_quark_wg), intent(inout) :: q1 real(DP), intent(in) :: coef call accum_mult(q1%eo(0),coef) call accum_mult(q1%eo(1),coef) return end subroutine
Subroutine : | |
q1e : | type(field_quark_eo_wg), intent(inout) |
coef : | complex(DP), intent(in) |
q1e <= q1e * coef
subroutine accum_cmult_fq_eo_wg(q1e,coef) ! ! q1e <= q1e * coef ! implicit none type(field_quark_eo_wg), intent(inout) :: q1e complex(DP), intent(in) :: coef integer :: ix,iy,iz,ieoxyz,itb,ic,is !$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+q1e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q1e%s(itb,iz,iy,ix)%y(:,:)=coef*q1e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
q1e : | type(field_quark_eo_wg), intent(inout) |
coef : | real(DP), intent(in) |
q1e <= q1e * coef
subroutine accum_rmult_fq_eo_wg(q1e,coef) ! ! q1e <= q1e * coef ! implicit none type(field_quark_eo_wg), intent(inout) :: q1e real(DP), intent(in) :: coef integer :: ix,iy,iz,ieoxyz,itb,ic,is !$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+q1e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q1e%s(itb,iz,iy,ix)%y(:,:)=coef*q1e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
q2 : | type(field_quark_wg), intent(inout) |
coef : | complex(DP), intent(in) |
q1 : | type(field_quark_wg), intent(in) |
q2 <= q2 * coef + q1
subroutine accum_cmult_add_fq_wg(q2,coef,q1) ! ! q2 <= q2 * coef + q1 ! implicit none type(field_quark_wg), intent(inout) :: q2 complex(DP), intent(in) :: coef type(field_quark_wg), intent(in) :: q1 call accum_mult_add(q2%eo(0),coef,q1%eo(0)) call accum_mult_add(q2%eo(1),coef,q1%eo(1)) return end subroutine
Subroutine : | |
q2 : | type(field_quark_wg), intent(inout) |
coef : | real(DP), intent(in) |
q1 : | type(field_quark_wg), intent(in) |
q2 <= q2 * coef + q1
subroutine accum_rmult_add_fq_wg(q2,coef,q1) ! ! q2 <= q2 * coef + q1 ! implicit none type(field_quark_wg), intent(inout) :: q2 real(DP), intent(in) :: coef type(field_quark_wg), intent(in) :: q1 call accum_mult_add(q2%eo(0),coef,q1%eo(0)) call accum_mult_add(q2%eo(1),coef,q1%eo(1)) return end subroutine
Subroutine : | |
q2e : | type(field_quark_eo_wg), intent(inout) |
coef : | complex(DP), intent(in) |
q1e : | type(field_quark_eo_wg), intent(in) |
q2e <= q2e * coef + q1e
subroutine accum_cmult_add_fq_eo_wg(q2e,coef,q1e) ! ! q2e <= q2e * coef + q1e ! implicit none type(field_quark_eo_wg), intent(inout) :: q2e complex(DP), intent(in) :: coef type(field_quark_eo_wg), intent(in) :: q1e integer :: ix,iy,iz,ieoxyz,itb,ic,is #ifdef _DEBUG_ if (q1e%ieo /= q2e%ieo) then if (nodeid == 0) then write(*,'("error in wqf_eo_wg_cmult2_accum_add(even/odd)")') endif stop endif #endif !$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+q2e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q2e%s(itb,iz,iy,ix)%y(:,:)=q1e%s(itb,iz,iy,ix)%y(:,:) +coef*q2e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
q2e : | type(field_quark_eo_wg), intent(inout) |
coef : | real(DP), intent(in) |
q1e : | type(field_quark_eo_wg), intent(in) |
q2e <= q2e * coef + q1e
subroutine accum_rmult_add_fq_eo_wg(q2e,coef,q1e) ! ! q2e <= q2e * coef + q1e ! implicit none type(field_quark_eo_wg), intent(inout) :: q2e real(DP), intent(in) :: coef type(field_quark_eo_wg), intent(in) :: q1e integer :: ix,iy,iz,ieoxyz,itb,ic,is #ifdef _DEBUG_ if (q1e%ieo /= q2e%ieo) then if (nodeid == 0) then write(*,'("error in wqf_eo_wg_rmult2_accum_add(even/odd)")') endif stop endif #endif !$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+q2e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q2e%s(itb,iz,iy,ix)%y(:,:)=q1e%s(itb,iz,iy,ix)%y(:,:) +coef*q2e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
y : | type(field_quark_wg), intent(inout) |
Multiply gamma_5 self (dirac rep.)
y <= gamma_5 y
subroutine accum_mult_gamma5(y) ! ! Multiply gamma_5 self (dirac rep.) ! ! y <= gamma_5 y ! implicit none type(field_quark_wg), intent(inout) :: y call accum_mult_gamma5_eo(y%eo(0)) call accum_mult_gamma5_eo(y%eo(1)) return end subroutine
Subroutine : | |||
ye : | type(field_quark_eo_wg),
intent(inout)
|
Multiply gamma_5 self (dirac rep.)
ye <= gamma_5 ye
subroutine accum_mult_gamma5_eo(ye) ! ! Multiply gamma_5 self (dirac rep.) ! ! ye <= gamma_5 ye ! implicit none type(field_quark_eo_wg), intent(inout) :: ye ! even/odd site fermion complex(DP) :: y1,y2,y3,y4 integer ix,iy,iz,itb,ic,ieoxyz !***************************** !* (0 0 1 0) !* gamma_5 = (0 0 0 1) !* (1 0 0 0) !* (0 1 0 0) !***************************** !$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ic,y1,y2,y3,y4) do ix=1,NX do iy=1,NY do iz=1,NZ ieoxyz=mod(ipeo+ix+iy+iz+ye%ieo,2) #ifdef _SF ye%s(0,iz,iy,ix)%y(:,:) = Z0 do itb=1,NTH-ieoxyz #else do itb=1-ieoxyz,NTH-ieoxyz #endif do ic=1,COL y1=ye%s(itb,iz,iy,ix)%y(ic,1) y2=ye%s(itb,iz,iy,ix)%y(ic,2) y3=ye%s(itb,iz,iy,ix)%y(ic,3) y4=ye%s(itb,iz,iy,ix)%y(ic,4) ye%s(itb,iz,iy,ix)%y(ic,1)=y3 ye%s(itb,iz,iy,ix)%y(ic,2)=y4 ye%s(itb,iz,iy,ix)%y(ic,3)=y1 ye%s(itb,iz,iy,ix)%y(ic,4)=y2 enddo enddo #ifdef _SF if (1==ieoxyz) then ye%s(NTH,iz,iy,ix)%y(:,:) = Z0 endif #endif enddo enddo enddo return end subroutine
Subroutine : | |
q2 : | type(field_quark_wg), intent(inout) |
q1 : | type(field_quark_wg), intent(in) |
q2 <= q2 - q1
subroutine accum_sub_fq_wg(q2,q1) ! ! q2 <= q2 - q1 ! implicit none type(field_quark_wg), intent(inout) :: q2 type(field_quark_wg), intent(in) :: q1 call accum_sub(q2%eo(0),q1%eo(0)) call accum_sub(q2%eo(1),q1%eo(1)) return end subroutine
Subroutine : | |
q2e : | type(field_quark_eo_wg), intent(inout) |
q1e : | type(field_quark_eo_wg), intent(in) |
q2e <= q2e - q1e
subroutine accum_sub_fq_eo_wg(q2e,q1e) ! ! q2e <= q2e - q1e ! implicit none type(field_quark_eo_wg), intent(inout) :: q2e type(field_quark_eo_wg), intent(in) :: q1e integer :: ix,iy,iz,ieoxyz,itb,ic,is #ifdef _DEBUG_ if (q1e%ieo /= q2e%ieo) then if (nodeid == 0) then write(*,'("error in wqf_eo_wg_accum_sub(even/odd)")') endif stop endif #endif !$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+q2e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q2e%s(itb,iz,iy,ix)%y(:,:)=q2e%s(itb,iz,iy,ix)%y(:,:) -q1e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
q2 : | type(field_quark_wg), intent(inout) |
q1 : | type(field_quark_wg), intent(in) |
q2 <= q1
subroutine assign_fq_wg(q2,q1) ! ! q2 <= q1 ! implicit none type(field_quark_wg), intent(inout) :: q2 type(field_quark_wg), intent(in) :: q1 call assign(q2%eo(0),q1%eo(0)) call assign(q2%eo(1),q1%eo(1)) return end subroutine
Subroutine : | |
q2e : | type(field_quark_eo_wg), intent(inout) |
q1e : | type(field_quark_eo_wg), intent(in) |
q2 <= q1
subroutine assign_fq_eo_wg(q2e,q1e) ! ! q2 <= q1 ! implicit none type(field_quark_eo_wg), intent(inout) :: q2e type(field_quark_eo_wg), intent(in) :: q1e integer :: ix,iy,iz,itb,ieoxyz q2e%ieo=q1e%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+q2e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q2e%s(itb,iz,iy,ix)%y(:,:)=q1e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
q3 : | type(field_quark_wg), intent(inout) |
q1 : | type(field_quark_wg), intent(in) |
q2 : | type(field_quark_wg), intent(in) |
q3 <= q1 + q2
subroutine assign_add_fq_wg(q3,q1,q2) ! ! q3 <= q1 + q2 ! implicit none type(field_quark_wg), intent(inout) :: q3 type(field_quark_wg), intent(in) :: q1,q2 call assign_add(q3%eo(0),q1%eo(0),q2%eo(0)) call assign_add(q3%eo(1),q1%eo(1),q2%eo(1)) return end subroutine
Subroutine : | |
q3e : | type(field_quark_eo_wg), intent(inout) |
q1e : | type(field_quark_eo_wg), intent(in) |
q2e : | type(field_quark_eo_wg), intent(in) |
q3e <= q1e + q2e
subroutine assign_add_fq_eo_wg(q3e,q1e,q2e) ! ! q3e <= q1e + q2e ! implicit none type(field_quark_eo_wg), intent(inout) :: q3e type(field_quark_eo_wg), intent(in) :: q1e,q2e integer :: ix,iy,iz,ieoxyz,itb,ic,is #ifdef _DEBUG_ if (q1e%ieo /= q2e%ieo) then if (nodeid == 0) then write(*,'("error in wqf_eo_wg_add(even/odd)")') endif stop endif #endif q3e%ieo=q1e%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+q3e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q3e%s(itb,iz,iy,ix)%y(:,:)=q1e%s(itb,iz,iy,ix)%y(:,:) +q2e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
q2 : | type(field_quark_wg), intent(inout) |
q1 : | type(field_quark_wg), intent(in) |
coef : | complex(DP), intent(in) |
q2 <= q1 * coef
subroutine assign_cmult_fq_wg(q2,q1,coef) ! ! q2 <= q1 * coef ! implicit none type(field_quark_wg), intent(inout) :: q2 type(field_quark_wg), intent(in) :: q1 complex(DP), intent(in) :: coef call assign_mult(q2%eo(0),q1%eo(0),coef) call assign_mult(q2%eo(1),q1%eo(1),coef) return end subroutine
Subroutine : | |
q2 : | type(field_quark_wg), intent(inout) |
q1 : | type(field_quark_wg), intent(in) |
coef : | real(DP), intent(in) |
q2 <= q1 * coef
subroutine assign_rmult_fq_wg(q2,q1,coef) ! ! q2 <= q1 * coef ! implicit none type(field_quark_wg), intent(inout) :: q2 type(field_quark_wg), intent(in) :: q1 real(DP), intent(in) :: coef call assign_mult(q2%eo(0),q1%eo(0),coef) call assign_mult(q2%eo(1),q1%eo(1),coef) return end subroutine
Subroutine : | |
q2e : | type(field_quark_eo_wg), intent(inout) |
q1e : | type(field_quark_eo_wg), intent(in) |
coef : | complex(DP), intent(in) |
q2e <= q1e * coef
subroutine assign_cmult_fq_eo_wg(q2e,q1e,coef) ! ! q2e <= q1e * coef ! implicit none type(field_quark_eo_wg), intent(inout) :: q2e type(field_quark_eo_wg), intent(in) :: q1e complex(DP), intent(in) :: coef integer :: ix,iy,iz,ieoxyz,itb,ic,is q2e%ieo=q1e%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+q2e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q2e%s(itb,iz,iy,ix)%y(:,:)=coef*q1e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
q2e : | type(field_quark_eo_wg), intent(inout) |
q1e : | type(field_quark_eo_wg), intent(in) |
coef : | real(DP), intent(in) |
q2e <= q1e * coef
subroutine assign_rmult_fq_eo_wg(q2e,q1e,coef) ! ! q2e <= q1e * coef ! implicit none type(field_quark_eo_wg), intent(inout) :: q2e type(field_quark_eo_wg), intent(in) :: q1e real(DP), intent(in) :: coef integer :: ix,iy,iz,ieoxyz,itb,ic,is q2e%ieo=q1e%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+q2e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q2e%s(itb,iz,iy,ix)%y(:,:)=coef*q1e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
q3 : | type(field_quark_wg), intent(inout) |
q1 : | type(field_quark_wg), intent(in) |
coef : | complex(DP), intent(in) |
q2 : | type(field_quark_wg), intent(in) |
q3 <= q1 * coef + q2
subroutine assign_cmult_add_fq_wg(q3,q1,coef,q2) ! ! q3 <= q1 * coef + q2 ! implicit none type(field_quark_wg), intent(inout) :: q3 type(field_quark_wg), intent(in) :: q1,q2 complex(DP), intent(in) :: coef call assign_mult_add(q3%eo(0),q1%eo(0),coef,q2%eo(0)) call assign_mult_add(q3%eo(1),q1%eo(1),coef,q2%eo(1)) return end subroutine
Subroutine : | |
q3 : | type(field_quark_wg), intent(inout) |
q1 : | type(field_quark_wg), intent(in) |
coef : | real(DP), intent(in) |
q2 : | type(field_quark_wg), intent(in) |
q3 <= q1 * coef + q2
subroutine assign_rmult_add_fq_wg(q3,q1,coef,q2) ! ! q3 <= q1 * coef + q2 ! implicit none type(field_quark_wg), intent(inout) :: q3 type(field_quark_wg), intent(in) :: q1,q2 real(DP), intent(in) :: coef call assign_mult_add(q3%eo(0),q1%eo(0),coef,q2%eo(0)) call assign_mult_add(q3%eo(1),q1%eo(1),coef,q2%eo(1)) return end subroutine
Subroutine : | |
q3e : | type(field_quark_eo_wg), intent(inout) |
q1e : | type(field_quark_eo_wg), intent(in) |
coef : | complex(DP), intent(in) |
q2e : | type(field_quark_eo_wg), intent(in) |
q3e <= q1e * coef + q2e
subroutine assign_cmult_add_fq_eo_wg(q3e,q1e,coef,q2e) ! ! q3e <= q1e * coef + q2e ! implicit none type(field_quark_eo_wg), intent(inout) :: q3e type(field_quark_eo_wg), intent(in) :: q1e,q2e complex(DP), intent(in) :: coef integer :: ix,iy,iz,ieoxyz,itb,ic,is #ifdef _DEBUG_ if (q1e%ieo /= q2e%ieo) then if (nodeid == 0) then write(*,'("error in wqf_eo_wg_cmult_add(even/odd)")') endif stop endif #endif q3e%ieo=q1e%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+q3e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q3e%s(itb,iz,iy,ix)%y(:,:)=coef*q1e%s(itb,iz,iy,ix)%y(:,:) +q2e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
q3e : | type(field_quark_eo_wg), intent(inout) |
q1e : | type(field_quark_eo_wg), intent(in) |
coef : | real(DP), intent(in) |
q2e : | type(field_quark_eo_wg), intent(in) |
q3e <= q1e * coef + q2e
subroutine assign_rmult_add_fq_eo_wg(q3e,q1e,coef,q2e) ! ! q3e <= q1e * coef + q2e ! implicit none type(field_quark_eo_wg), intent(inout) :: q3e type(field_quark_eo_wg), intent(in) :: q1e,q2e real(DP), intent(in) :: coef integer :: ix,iy,iz,ieoxyz,itb,ic,is #ifdef _DEBUG_ if (q1e%ieo /= q2e%ieo) then if (nodeid == 0) then write(*,'("error in wqf_eo_wg_rmult_add(even/odd)")') endif stop endif #endif q3e%ieo=q1e%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+q3e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q3e%s(itb,iz,iy,ix)%y(:,:)=coef*q1e%s(itb,iz,iy,ix)%y(:,:) +q2e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
gy : | type(field_quark_wg), intent(out) |
y : | type(field_quark_wg), intent(in) |
Multiply gamma_5 (dirac rep.)
gy <= gamma_5 y
subroutine assign_mult_gamma5(gy,y) ! ! Multiply gamma_5 (dirac rep.) ! ! gy <= gamma_5 y ! implicit none type(field_quark_wg), intent(in) :: y type(field_quark_wg), intent(out):: gy call assign_mult_gamma5_eo(gy%eo(0),y%eo(0)) call assign_mult_gamma5_eo(gy%eo(1),y%eo(1)) return end subroutine
Subroutine : | |||
gye : | type(field_quark_eo_wg),
intent(out)
| ||
ye : | type(field_quark_eo_wg),
intent(in)
|
Multiply gamma_5 (dirac rep.)
gye <= gamma_5 ye
subroutine assign_mult_gamma5_eo(gye,ye) ! ! Multiply gamma_5 (dirac rep.) ! ! gye <= gamma_5 ye ! implicit none type(field_quark_eo_wg), intent(in) :: ye ! even/odd site fermion type(field_quark_eo_wg), intent(out):: gye ! even/odd site fermion complex(DP) :: y1,y2,y3,y4 integer ix,iy,iz,itb,ic,ieoxyz gye%ieo=ye%ieo !***************************** !* (0 0 1 0) !* gamma_5 = (0 0 0 1) !* (1 0 0 0) !* (0 1 0 0) !***************************** !$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ic,y1,y2,y3,y4) do ix=1,NX do iy=1,NY do iz=1,NZ ieoxyz=mod(ipeo+ix+iy+iz+ye%ieo,2) #ifdef _SF gye%s(0,iz,iy,ix)%y(:,:) = Z0 do itb=1,NTH-ieoxyz #else do itb=1-ieoxyz,NTH-ieoxyz #endif do ic=1,COL y1=ye%s(itb,iz,iy,ix)%y(ic,1) y2=ye%s(itb,iz,iy,ix)%y(ic,2) y3=ye%s(itb,iz,iy,ix)%y(ic,3) y4=ye%s(itb,iz,iy,ix)%y(ic,4) gye%s(itb,iz,iy,ix)%y(ic,1)=y3 gye%s(itb,iz,iy,ix)%y(ic,2)=y4 gye%s(itb,iz,iy,ix)%y(ic,3)=y1 gye%s(itb,iz,iy,ix)%y(ic,4)=y2 enddo enddo #ifdef _SF if (1==ieoxyz) then gye%s(NTH,iz,iy,ix)%y(:,:) = Z0 endif #endif enddo enddo enddo return end subroutine
Subroutine : | |
q3 : | type(field_quark_wg), intent(inout) |
q1 : | type(field_quark_wg), intent(in) |
coef : | complex(DP), intent(in) |
q2 : | type(field_quark_wg), intent(in) |
q3 <= q1 * coef - q2
subroutine assign_cmult_sub_fq_wg(q3,q1,coef,q2) ! ! q3 <= q1 * coef - q2 ! implicit none type(field_quark_wg), intent(inout) :: q3 type(field_quark_wg), intent(in) :: q1,q2 complex(DP), intent(in) :: coef call assign_mult_sub(q3%eo(0),q1%eo(0),coef,q2%eo(0)) call assign_mult_sub(q3%eo(1),q1%eo(1),coef,q2%eo(1)) return end subroutine
Subroutine : | |
q3 : | type(field_quark_wg), intent(inout) |
q1 : | type(field_quark_wg), intent(in) |
coef : | real(DP), intent(in) |
q2 : | type(field_quark_wg), intent(in) |
q3 = q1 * coef - q2
subroutine assign_rmult_sub_fq_wg(q3,q1,coef,q2) ! ! q3 = q1 * coef - q2 ! implicit none type(field_quark_wg), intent(inout) :: q3 type(field_quark_wg), intent(in) :: q1,q2 real(DP), intent(in) :: coef call assign_mult_sub(q3%eo(0),q1%eo(0),coef,q2%eo(0)) call assign_mult_sub(q3%eo(1),q1%eo(1),coef,q2%eo(1)) return end subroutine
Subroutine : | |
q3e : | type(field_quark_eo_wg), intent(inout) |
q1e : | type(field_quark_eo_wg), intent(in) |
coef : | complex(DP), intent(in) |
q2e : | type(field_quark_eo_wg), intent(in) |
q3e <= q1e * coef - q2e
subroutine assign_cmult_sub_fq_eo_wg(q3e,q1e,coef,q2e) ! ! q3e <= q1e * coef - q2e ! implicit none type(field_quark_eo_wg), intent(inout) :: q3e type(field_quark_eo_wg), intent(in) :: q1e,q2e complex(DP), intent(in) :: coef integer :: ix,iy,iz,ieoxyz,itb,ic,is #ifdef _DEBUG_ if (q1e%ieo /= q2e%ieo) then if (nodeid == 0) then write(*,'("error in wqf_eo_wg_cmult_sub(even/odd)")') endif stop endif #endif q3e%ieo=q1e%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+q3e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q3e%s(itb,iz,iy,ix)%y(:,:)=coef*q1e%s(itb,iz,iy,ix)%y(:,:) -q2e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
q3e : | type(field_quark_eo_wg), intent(inout) |
q1e : | type(field_quark_eo_wg), intent(in) |
coef : | real(DP), intent(in) |
q2e : | type(field_quark_eo_wg), intent(in) |
q3e = q1e * coef - q2e
subroutine assign_rmult_sub_fq_eo_wg(q3e,q1e,coef,q2e) ! ! q3e = q1e * coef - q2e ! implicit none type(field_quark_eo_wg), intent(inout) :: q3e type(field_quark_eo_wg), intent(in) :: q1e,q2e real(DP), intent(in) :: coef integer :: ix,iy,iz,ieoxyz,itb,ic,is #ifdef _DEBUG_ if (q1e%ieo /= q2e%ieo) then if (nodeid == 0) then write(*,'("error in wqf_eo_wg_rmult_sub(even/odd)")') endif stop endif #endif q3e%ieo=q1e%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+q3e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q3e%s(itb,iz,iy,ix)%y(:,:)=coef*q1e%s(itb,iz,iy,ix)%y(:,:) -q2e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
q3 : | type(field_quark_wg), intent(inout) |
q1 : | type(field_quark_wg), intent(in) |
q2 : | type(field_quark_wg), intent(in) |
q3 <= q1 - q2
subroutine assign_sub_fq_wg(q3,q1,q2) ! ! q3 <= q1 - q2 ! implicit none type(field_quark_wg), intent(inout) :: q3 type(field_quark_wg), intent(in) :: q1,q2 call assign_sub(q3%eo(0),q1%eo(0),q2%eo(0)) call assign_sub(q3%eo(1),q1%eo(1),q2%eo(1)) return end subroutine
Subroutine : | |
q3e : | type(field_quark_eo_wg), intent(inout) |
q1e : | type(field_quark_eo_wg), intent(in) |
q2e : | type(field_quark_eo_wg), intent(in) |
q3e <= q1e - q2e
subroutine assign_sub_fq_eo_wg(q3e,q1e,q2e) ! ! q3e <= q1e - q2e ! implicit none type(field_quark_eo_wg), intent(inout) :: q3e type(field_quark_eo_wg), intent(in) :: q1e,q2e integer :: ix,iy,iz,ieoxyz,itb,ic,is #ifdef _DEBUG_ if (q1e%ieo /= q2e%ieo) then if (nodeid == 0) then write(*,'("error in wqf_eo_wg_sub(even/odd)")') endif stop endif #endif q3e%ieo=q1e%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+q3e%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz q3e%s(itb,iz,iy,ix)%y(:,:)=q1e%s(itb,iz,iy,ix)%y(:,:) -q2e%s(itb,iz,iy,ix)%y(:,:) enddo enddo enddo enddo return end subroutine
Subroutine : | |
q : | type(field_quark_wg), intent(inout) |
set zero on field
subroutine clear_fq_wg(q) ! ! set zero on field ! implicit none type(field_quark_wg), intent(inout) :: q call clear(q%eo(0)) call clear(q%eo(1)) return end subroutine
Subroutine : | |
q : | type(field_quark_wog), intent(inout) |
set zero on field
subroutine clear_fq_wog(q) ! ! set zero on field ! implicit none type(field_quark_wog), intent(inout) :: q call clear(q%eo(0)) call clear(q%eo(1)) return end subroutine
Subroutine : | |
qe : | type(field_quark_eo_wg), intent(inout) |
set zero on field
subroutine clear_fq_eo_wg(qe) ! ! set zero on field ! implicit none type(field_quark_eo_wg), intent(inout) :: qe integer :: ix,iy,iz,itb,ieoxyz !$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+qe%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz qe%s(itb,iz,iy,ix)%y(:,:)=Z0 enddo enddo enddo enddo return end subroutine
Subroutine : | |
qe : | type(field_quark_eo_wog), intent(inout) |
set zero on field
subroutine clear_fq_eo_wog(qe) ! ! set zero on field ! implicit none type(field_quark_eo_wog), intent(inout) :: qe 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 qe%s(itb,iz,iy,ix)%y(:,:)=Z0 enddo enddo enddo enddo return end subroutine
Subroutine : | |
y : | type(field_quark_wg), intent(inout) |
Boundary copy
subroutine copy_fq_wg(y) ! ! Boundary copy ! implicit none type(field_quark_wg), intent(inout) :: y type(fields_fermion) :: fields call copy_boundary(y%eo(0)) call copy_boundary(y%eo(1)) return end subroutine
Subroutine : | |
ye : | type(field_quark_eo_wg), intent(inout) |
Boundary copy on even/odd sites only
subroutine copy_fq_eo_wg(ye) ! ! Boundary copy on even/odd sites only ! implicit none type(field_quark_eo_wg), intent(inout) :: ye type(fields_fermion) :: fields integer :: ix,iy,iz,itb,ieoxyz,ic,is,itb0,itb1 #ifndef _singlePU integer :: ibuff #endif if (.not.m_is_initialized) call new_fields_fermion(fields) call tic(copy_fq_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+ye%ieo,2) #ifdef _SF !================================ ! itt=1 and NTT+1 is SF boundary ! itt=0 is unphysical ! set zero !================================ ye%s(0,iz,iy,ix)%y(:,:) = Z0 ! itt=0 or 1 if (1==ieoxyz) then ye%s(NTH,iz,iy,ix)%y(:,:) = Z0 ! itt=NTT+1 endif #else itb0= ieoxyz *NTH itb1=(1-ieoxyz)*NTH ye%s(itb0,iz,iy,ix)=ye%s(itb1,iz,iy,ix) #endif enddo enddo enddo !--------------------------------- ! copy without cornar ! set 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 ye%s(itb, 0,iy,ix)=ye%s(itb,NZ,iy,ix) ye%s(itb,NZ1,iy,ix)=ye%s(itb, 1,iy,ix) enddo enddo enddo #else ibuff=0 do ix=1,NX do iy=1,NY do itb=0,NTH do is=1,SPIN do ic=1,COL ibuff=ibuff+1 fbuffup(ibuff,1,3)=ye%s(itb,NZ,iy,ix)%y(ic,is) fbuffdn(ibuff,1,3)=ye%s(itb, 1,iy,ix)%y(ic,is) enddo enddo enddo enddo enddo if (ibuff.NE.FBUFF_SIZE(3)) then write(*,'(" Error : com-buff size in Z")') stop 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 ye%s(itb,iz, 0,ix)=ye%s(itb,iz,NY,ix) ye%s(itb,iz,NY1,ix)=ye%s(itb,iz, 1,ix) enddo enddo enddo #else ibuff=0 do ix=1,NX do iz=1,NZ do itb=0,NTH do is=1,SPIN do ic=1,COL ibuff=ibuff+1 fbuffup(ibuff,1,2)=ye%s(itb,iz,NY,ix)%y(ic,is) fbuffdn(ibuff,1,2)=ye%s(itb,iz, 1,ix)%y(ic,is) enddo enddo enddo enddo enddo if (ibuff.NE.FBUFF_SIZE(2)) then write(*,'(" Error : com-buff size in Y")') stop 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 ye%s(itb,iz,iy, 0)=ye%s(itb,iz,iy,NX) ye%s(itb,iz,iy,NX1)=ye%s(itb,iz,iy, 1) enddo enddo enddo #else ibuff=0 do iy=1,NY do iz=1,NZ do itb=0,NTH do is=1,SPIN do ic=1,COL ibuff=ibuff+1 fbuffup(ibuff,1,1)=ye%s(itb,iz,iy,NX)%y(ic,is) fbuffdn(ibuff,1,1)=ye%s(itb,iz,iy, 1)%y(ic,is) enddo enddo enddo enddo enddo if (ibuff.NE.FBUFF_SIZE(1)) then write(*,'(" Error : com-buff size in X")') stop endif #endif #if _NDIMZ != 1 call comlib_sendrecv(idfsendup(3)) call comlib_sendrecv(idfsenddn(3)) #endif #if _NDIMY != 1 call comlib_sendrecv(idfsendup(2)) call comlib_sendrecv(idfsenddn(2)) #endif #if _NDIMX != 1 call comlib_sendrecv(idfsendup(1)) call comlib_sendrecv(idfsenddn(1)) #endif !--------------------------------- ! receive !--------------------------------- !****** Z - boundaray #if _NDIMZ != 1 ibuff=0 do ix=1,NX do iy=1,NY do itb=0,NTH do is=1,SPIN do ic=1,COL ibuff=ibuff+1 ye%s(itb, 0,iy,ix)%y(ic,is)=fbuffup(ibuff,2,3) ye%s(itb,NZ1,iy,ix)%y(ic,is)=fbuffdn(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 is=1,SPIN do ic=1,COL ibuff=ibuff+1 ye%s(itb,iz, 0,ix)%y(ic,is)=fbuffup(ibuff,2,2) ye%s(itb,iz,NY1,ix)%y(ic,is)=fbuffdn(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 is=1,SPIN do ic=1,COL ibuff=ibuff+1 ye%s(itb,iz,iy, 0)%y(ic,is)=fbuffup(ibuff,2,1) ye%s(itb,iz,iy,NX1)%y(ic,is)=fbuffdn(ibuff,2,1) enddo enddo enddo enddo enddo #endif call toc(copy_fq_time) return end subroutine
Derived Type : | |||
s(0:NTH,0:NZ1,0:NY1,0:NX1) : | type(su3fv_spinor)
| ||
ieo : | integer
| ||
idummy(3) : | integer
|
quark field on even/odd sites with ghost sites
Derived Type : | |
s(NTH,NZ,NY,NX) : | type(su3fv_spinor) |
ieo : | integer |
idummy(3) : | integer |
quark field on even/odd sites without ghost sites
Derived Type : | |||
eo(0:1) : | type(field_quark_eo_wg)
|
quark field with ghost sites
Subroutine : | |||
BB : | type(vfield_gluon_wg), intent(inout)
| ||
fcoef : | real(DP), intent(in)
| ||
fx : | type(field_quark_wg),
intent(in)
| ||
fy : | type(field_quark_wg),
intent(in)
|
Calc MD force from hopping matrix
User should copy boundary of the input field (fx,fy) before calling this subroutine.
\[
BB_{\mu}(n) = BB_{\mu}(n) + \mathrm{fcoef}\times \mathrm{tr}\left[(1-\gamma_{\mu}) fy(n+\hat{\mu}) fx(n)^{\dag} + fx(n+\hat{\mu}) fy(n)^{\dag} (1+\gamma_{\mu})\right].
\] this comes from \[
\delta S = fx^{\dag}\delta M_{hop} fy + fy^{\dag}\delta M_{hop}^{\dag} fx,
\] with hopping matrix $ M_{hop} $ .
subroutine force_hmc_hopping(BB,fcoef,fx,fy) ! ! Calc MD force from hopping matrix ! ! User should \copy boundary of the input field (fx,fy) before calling this subroutine. ! !\[ ! BB_{\mu}(n) = BB_{\mu}(n) ! + \mathrm{fcoef}\times \mathrm{tr}\left[(1-\gamma_{\mu}) fy(n+\hat{\mu}) fx(n)^{\dag} ! + fx(n+\hat{\mu}) fy(n)^{\dag} (1+\gamma_{\mu})\right]. !\] ! this comes from !\[ ! \delta S = fx^{\dag}\delta M_{hop} fy + fy^{\dag}\delta M_{hop}^{\dag} fx, !\] ! with hopping matrix $ M_{hop} $ . ! use field_gauge_class implicit none type(vfield_gluon_wg), intent(inout) :: BB ! pre-force real(DP), intent(in) :: fcoef ! force coefficient type(field_quark_wg), intent(in) :: fx ! external fermion field type(field_quark_wg), intent(in) :: fy ! external fermion field call force_hmc_hopping_eo(BB%eo(0),fcoef,fx%eo(0),fx%eo(1),fy%eo(0),fy%eo(1)) call force_hmc_hopping_eo(BB%eo(1),fcoef,fx%eo(1),fx%eo(0),fy%eo(1),fy%eo(0)) return end subroutine
Subroutine : | |||
BBe : | type(vfield_gluon_eo_wg), intent(inout)
| ||
fcoef : | real(DP), intent(in)
| ||
fxe : | type(field_quark_eo_wg),
intent(in)
| ||
fxo : | type(field_quark_eo_wg),
intent(in)
| ||
fye : | type(field_quark_eo_wg),
intent(in)
| ||
fyo : | type(field_quark_eo_wg),
intent(in)
|
Calc MD force from hopping matrix (even/odd sites only)
User should copy boundary of the input field (fx,fy) before calling this subroutine.
\[
BB_{\mu}(n) = BB_{\mu}(n) + \mathrm{fcoef}\times \mathrm{tr}\left[(1-\gamma_{\mu}) fy(n+\hat{\mu}) fx(n)^{\dag} + fx(n+\hat{\mu}) fy(n)^{\dag} (1+\gamma_{\mu})\right].
\] this comes from \[
\delta S = fx^{\dag}\delta M_{hop} fy + fy^{\dag}\delta M_{hop}^{\dag} fx,
\] with hopping matrix $ M_{hop} $ .
subroutine force_hmc_hopping_eo(BBe,fcoef,fxe,fxo,fye,fyo) ! ! Calc MD force from hopping matrix (even/odd sites only) ! ! User should \copy boundary of the input field (fx,fy) before calling this subroutine. ! !\[ ! BB_{\mu}(n) = BB_{\mu}(n) ! + \mathrm{fcoef}\times \mathrm{tr}\left[(1-\gamma_{\mu}) fy(n+\hat{\mu}) fx(n)^{\dag} ! + fx(n+\hat{\mu}) fy(n)^{\dag} (1+\gamma_{\mu})\right]. !\] ! this comes from !\[ ! \delta S = fx^{\dag}\delta M_{hop} fy + fy^{\dag}\delta M_{hop}^{\dag} fx, !\] ! with hopping matrix $ M_{hop} $ . ! use field_gauge_class implicit none type(vfield_gluon_eo_wg), intent(inout) :: BBe ! even/odd site force contribution (\dot{u}) real(DP), intent(in) :: fcoef ! force coefficient type(field_quark_eo_wg), intent(in) :: fxe ! even/odd site pseudo fermion type(field_quark_eo_wg), intent(in) :: fxo ! odd/even site pseudo fermion type(field_quark_eo_wg), intent(in) :: fye ! even/odd site pseudo fermion type(field_quark_eo_wg), intent(in) :: fyo ! odd/even site pseudo fermion complex(DP) :: b1e(COL,SPIN),b2e(COL,SPIN) integer :: ix,iy,iz,itb,itb2,ieoxyz,ic,jc #define _X_ itb,iz,iy,ix #define _XTUP_ itb2,iz,iy,ix #define _XZUP_ itb,iz+1,iy,ix #define _XYUP_ itb,iz,iy+1,ix #define _XXUP_ itb,iz,iy,ix+1 #define MULTJ(x) cmplx(-aimag(x),real(x),kind=DP) if ( .not. ((BBe%mu(1)%ieo==fxe%ieo).and. (BBe%mu(2)%ieo==fxe%ieo).and. (BBe%mu(3)%ieo==fxe%ieo).and. (BBe%mu(4)%ieo==fxe%ieo).and. (fxe%ieo==fye%ieo).and. (fxo%ieo==fyo%ieo).and. (fxe%ieo/=fxo%ieo).and. (fye%ieo/=fyo%ieo)) ) then write(*,*)"error in force_hmc_hopping_eo (even/odd)" stop endif !************************************************************************ !* !* BB = BB + fcoef * tr[(1-gamma_mu) fy(n+mu) fx(n)^{+} !* + fx(n+mu) fy(n)^{+} (1+gamma_mu)] !* !* = BB + fcoef * tr[ b1(n) fx(n)^{+} + fx(n+mu) b2(n)^{+} ] !* !* where !* !* b1(n)=(1-gamma_mu)fy(n+mu) !* b2(n)=(1+gamma_mu)fy(n) !* !************************************************************************ !$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,itb2,b1e,b2e) do ix=1,NX do iy=1,NY do iz=1,NZ ieoxyz=mod(ipeo+ix+iy+iz+fxe%ieo,2) do itb=1-ieoxyz,NTH-ieoxyz !================ ! T-direction !================ itb2=itb+ieoxyz b1e(:,3) = fyo%s(_XTUP_)%y(:,3)*2.0_DP b1e(:,4) = fyo%s(_XTUP_)%y(:,4)*2.0_DP b2e(:,1) = fye%s(_X_ )%y(:,1)*2.0_DP b2e(:,2) = fye%s(_X_ )%y(:,2)*2.0_DP do jc=1,COL do ic=1,COL BBe%mu(4)%s(_X_)%u(ic,jc)=BBe%mu(4)%s(_X_)%u(ic,jc) +fcoef*( +conjg( fxe%s(_X_ )%y(jc,3))* b1e(ic,3) +conjg( fxe%s(_X_ )%y(jc,4))* b1e(ic,4) + fxo%s(_XTUP_)%y(ic,1) *conjg(b2e(jc,1)) + fxo%s(_XTUP_)%y(ic,2) *conjg(b2e(jc,2)) ) enddo enddo !================ ! Z-direction !================ b1e(:,1) = fyo%s(_XZUP_)%y(:,1) + MULTJ(fyo%s(_XZUP_)%y(:,3)) b1e(:,2) = fyo%s(_XZUP_)%y(:,2) - MULTJ(fyo%s(_XZUP_)%y(:,4)) b1e(:,3) = -MULTJ(b1e(:,1)) b1e(:,4) = +MULTJ(b1e(:,2)) b2e(:,1) = fye%s(_X_ )%y(:,1) - MULTJ(fye%s(_X_ )%y(:,3)) b2e(:,2) = fye%s(_X_ )%y(:,2) + MULTJ(fye%s(_X_ )%y(:,4)) b2e(:,3) = +MULTJ(b2e(:,1)) b2e(:,4) = -MULTJ(b2e(:,2)) do jc=1,COL do ic=1,COL BBe%mu(3)%s(_X_)%u(ic,jc)=BBe%mu(3)%s(_X_)%u(ic,jc) +fcoef*( +conjg( fxe%s(_X_ )%y(jc,1))* b1e(ic,1) +conjg( fxe%s(_X_ )%y(jc,2))* b1e(ic,2) +conjg( fxe%s(_X_ )%y(jc,3))* b1e(ic,3) +conjg( fxe%s(_X_ )%y(jc,4))* b1e(ic,4) + fxo%s(_XZUP_)%y(ic,1) *conjg(b2e(jc,1)) + fxo%s(_XZUP_)%y(ic,2) *conjg(b2e(jc,2)) + fxo%s(_XZUP_)%y(ic,3) *conjg(b2e(jc,3)) + fxo%s(_XZUP_)%y(ic,4) *conjg(b2e(jc,4)) ) enddo enddo !================ ! Y-direction !================ b1e(:,1) = fyo%s(_XYUP_)%y(:,1) + fyo%s(_XYUP_)%y(:,4) b1e(:,2) = fyo%s(_XYUP_)%y(:,2) - fyo%s(_XYUP_)%y(:,3) b1e(:,3) = -b1e(:,2) b1e(:,4) = +b1e(:,1) b2e(:,1) = fye%s(_X_ )%y(:,1) - fye%s(_X_ )%y(:,4) b2e(:,2) = fye%s(_X_ )%y(:,2) + fye%s(_X_ )%y(:,3) b2e(:,3) = +b2e(:,2) b2e(:,4) = -b2e(:,1) do jc=1,COL do ic=1,COL BBe%mu(2)%s(_X_)%u(ic,jc)=BBe%mu(2)%s(_X_)%u(ic,jc) +fcoef*( +conjg( fxe%s(_X_ )%y(jc,1))* b1e(ic,1) +conjg( fxe%s(_X_ )%y(jc,2))* b1e(ic,2) +conjg( fxe%s(_X_ )%y(jc,3))* b1e(ic,3) +conjg( fxe%s(_X_ )%y(jc,4))* b1e(ic,4) + fxo%s(_XYUP_)%y(ic,1) *conjg(b2e(jc,1)) + fxo%s(_XYUP_)%y(ic,2) *conjg(b2e(jc,2)) + fxo%s(_XYUP_)%y(ic,3) *conjg(b2e(jc,3)) + fxo%s(_XYUP_)%y(ic,4) *conjg(b2e(jc,4)) ) enddo enddo !================ ! X-direction !================ b1e(:,1) = fyo%s(_XXUP_)%y(:,1) + MULTJ(fyo%s(_XXUP_)%y(:,4)) b1e(:,2) = fyo%s(_XXUP_)%y(:,2) + MULTJ(fyo%s(_XXUP_)%y(:,3)) b1e(:,3) = -MULTJ(b1e(:,2)) b1e(:,4) = -MULTJ(b1e(:,1)) b2e(:,1) = fye%s(_X_ )%y(:,1) - MULTJ(fye%s(_X_ )%y(:,4)) b2e(:,2) = fye%s(_X_ )%y(:,2) - MULTJ(fye%s(_X_ )%y(:,3)) b2e(:,3) = +MULTJ(b2e(:,2)) b2e(:,4) = +MULTJ(b2e(:,1)) do jc=1,COL do ic=1,COL BBe%mu(1)%s(_X_)%u(ic,jc)=BBe%mu(1)%s(_X_)%u(ic,jc) +fcoef*( +conjg( fxe%s(_X_ )%y(jc,1))* b1e(ic,1) +conjg( fxe%s(_X_ )%y(jc,2))* b1e(ic,2) +conjg( fxe%s(_X_ )%y(jc,3))* b1e(ic,3) +conjg( fxe%s(_X_ )%y(jc,4))* b1e(ic,4) + fxo%s(_XXUP_)%y(ic,1) *conjg(b2e(jc,1)) + fxo%s(_XXUP_)%y(ic,2) *conjg(b2e(jc,2)) + fxo%s(_XXUP_)%y(ic,3) *conjg(b2e(jc,3)) + fxo%s(_XXUP_)%y(ic,4) *conjg(b2e(jc,4)) ) enddo enddo enddo enddo enddo enddo ! end of do ix,iy,iz #undef _X_ #undef _XTUP_ #undef _XZUP_ #undef _XYUP_ #undef _XXUP_ #undef MULTJ return end subroutine
Subroutine : | |||
yde : | type(field_quark_eo_wg),
intent(inout)
| ||
yo : | type(field_quark_eo_wg),
intent(in)
| ||
u : | type(vfield_gluon_wg), intent(in)
|
Multiply hopping matrix (odd->even/even->odd sites only)
User should copy boundary of the input field before calling this subroutine.
yde <= Meo yo
\[
M(n,m) = \sum_{\mu=1}^{4} \left[ (1-\gamma_{\mu})U_{\mu}(n)\delta_{n+\hat{\mu},m} +(1+\gamma_{\mu})U_{\mu}^{\dag}(m)\delta_{n-\hat{\mu},m}\right]
\]
subroutine mult_hopping_tzyx_eo(yde,yo,u) ! ! Multiply hopping matrix (odd->even/even->odd sites only) ! ! User should copy boundary of the input field before calling this subroutine. ! ! yde <= Meo yo ! !\[ ! M(n,m) = \sum_{\mu=1}^{4} \left[ (1-\gamma_{\mu})U_{\mu}(n)\delta_{n+\hat{\mu},m} ! +(1+\gamma_{\mu})U_{\mu}^{\dag}(m)\delta_{n-\hat{\mu},m}\right] !\] ! use field_gauge_class implicit none type(field_quark_eo_wg), intent(inout):: yde ! even/odd site fermion vector (output) type(field_quark_eo_wg), intent(in) :: yo ! odd/even site fermion vector (input) type(vfield_gluon_wg), intent(in) :: u ! gauge field type(su3fv_spinor) :: yt complex(DP) :: w(COL,SPIN/2),y(COL,SPIN/2) integer :: ix,iy,iz,itb,itbup,itbdn,ieoxyz,ieo,ioe ioe = yo%ieo ieo = 1-ioe yde%ieo = ieo !$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,itbup,itbdn,yt,w,y) do ix=1,NX do iy=1,NY do iz=1,NZ ieoxyz=mod(ipeo+ix+iy+iz+yde%ieo,2) itbup= +ieoxyz itbdn=-(1-ieoxyz) #ifdef _SF yde%s(0,iz,iy,ix)%y(:,:) = Z0 do itb=1,NTH-ieoxyz #else do itb=1-ieoxyz,NTH-ieoxyz #endif #include "mult_hop_kernel.h90" yde%s(itb,iz,iy,ix) = yt enddo #ifdef _SF if (1==ieoxyz) then yde%s(NTH,iz,iy,ix)%y(:,:) = Z0 endif #endif enddo enddo enddo return end subroutine
Subroutine : | |
q : | type(field_quark_wg), intent(inout) |
Initialize field
subroutine new_fq_wg(q) ! ! Initialize field ! implicit none type(field_quark_wg), intent(inout) :: q call new(q%eo(0),0) call new(q%eo(1),1) return end subroutine
Subroutine : | |
q : | type(field_quark_wog), intent(inout) |
Initialize field
subroutine new_fq_wog(q) ! ! Initialize field ! implicit none type(field_quark_wog), intent(inout) :: q call new(q%eo(0),0) call new(q%eo(1),1) return end subroutine
Subroutine : | |
qe : | type(field_quark_eo_wg), intent(inout) |
ieo : | integer, intent(in) |
Initialize field by setting the even/odd ness
subroutine new_fq_eo_wg(qe,ieo) ! ! Initialize field by setting the even/odd ness ! implicit none type(field_quark_eo_wg), intent(inout) :: qe integer, intent(in) :: ieo qe%ieo=ieo return end subroutine
Subroutine : | |
qe : | type(field_quark_eo_wog), intent(inout) |
ieo : | integer, intent(in) |
Initialize field by setting the even/odd ness
subroutine new_fq_eo_wog(qe,ieo) ! ! Initialize field by setting the even/odd ness ! implicit none type(field_quark_eo_wog), intent(inout) :: qe integer, intent(in) :: ieo qe%ieo=ieo return end subroutine
Subroutine : | |
y : | type(field_quark_wg), intent(in) |
v(COL*SPIN*NT*NZ*NY*NX) : | complex(DP), intent(out) |
\Pack (copy) fermion field to one dimensional complex array Ghost sites data are not copied. User should keep the even-odd-ness of the field in mind.
subroutine pack_fq_wg(y,v) ! ! \Pack (copy) fermion field to one dimensional complex array ! Ghost sites data are not copied. ! ! User should keep the even-odd-ness of the field in mind. ! implicit none type(field_quark_wg), intent(in) :: y complex(DP), intent(out) :: v(COL*SPIN*NT*NZ*NY*NX) integer, parameter :: NSITE = COL*SPIN*NT*NZ*NY*NX integer, parameter :: NHSITE = COL*SPIN*NTH*NZ*NY*NX call pack(y%eo(0),v(1:NHSITE)) call pack(y%eo(1),v(NHSITE+1:NSITE)) return end subroutine
Subroutine : | |
yeo : | type(field_quark_eo_wg), intent(in) |
vec(COL*SPIN*NTH*NZ*NY*NX) : | complex(DP), intent(out) |
\Pack (copy) fermion field to one dimensional complex array Ghost sites data are not copied. User should keep the even-odd-ness of the field in mind.
subroutine pack_fq_eo_wg(yeo,vec) ! ! \Pack (copy) fermion field to one dimensional complex array ! Ghost sites data are not copied. ! ! User should keep the even-odd-ness of the field in mind. ! implicit none type(field_quark_eo_wg), intent(in) :: yeo complex(DP), intent(out) :: vec(COL*SPIN*NTH*NZ*NY*NX) integer :: ix,iy,iz,itb,ic,is,ieoxyz,isite integer :: ix_base integer :: iy_base integer :: iz_base integer :: itb_base integer :: is_base !$OMP PARALLEL DO & !$OMP PRIVATE(ix,iy,iz,ieoxyz,itb,ic,is, & !$OMP ix_base,iy_base,iz_base,itb_base,is_base,isite) do ix=1,NX ix_base=(ix -1)*COL*SPIN*NTH*NZ*NY do iy=1,NY iy_base=(iy -1)*COL*SPIN*NTH*NZ +ix_base do iz=1,NZ iz_base=(iz -1)*COL*SPIN*NTH +iy_base ieoxyz=mod(ipeo+ix+iy+iz+yeo%ieo,2) #ifdef _SF if (1==ieoxyz) then !=============================== ! itt=1 is SF boundary ! set zero on fermion field !=============================== itb = 0 itb_base=(itb+ieoxyz-1)*COL*SPIN +iz_base do is=1,SPIN is_base=(is -1)*COL +itb_base do ic=1,COL isite= ic+is_base vec(isite) = Z0 enddo enddo endif do itb=1,NTH-ieoxyz itb_base=(itb+ieoxyz-1)*COL*SPIN +iz_base do is=1,SPIN is_base=(is -1)*COL +itb_base do ic=1,COL isite= ic+is_base vec(isite)=yeo%s(itb,iz,iy,ix)%y(ic,is) enddo enddo enddo #else do itb=1-ieoxyz,NTH-ieoxyz itb_base=(itb+ieoxyz-1)*COL*SPIN +iz_base do is=1,SPIN is_base=(is -1)*COL +itb_base do ic=1,COL isite= ic+is_base vec(isite)=yeo%s(itb,iz,iy,ix)%y(ic,is) enddo enddo enddo #endif enddo enddo enddo return end subroutine
Subroutine : |
subroutine print_copy_fq_statistics() implicit none integer :: mu if (.not.m_is_initialized) return if (0 == nodeid) then do mu=1,NDIM-1 call comlib_print_statistics(idfsendup(mu)) call comlib_print_statistics(idfsenddn(mu)) enddo endif return end subroutine
Function : | |
my_prod : | complex(DP) |
q1 : | type(field_quark_wg), intent(in) |
q2 : | type(field_quark_wg), intent(in) |
return inner-product : q1’ .dot. q2
function prod_fq_wg( q1, q2 ) result(my_prod) ! ! return inner-product : q1' .dot. q2 ! implicit none type(field_quark_wg), intent(in) :: q1,q2 complex(DP) :: my_prod my_prod = prod(q1%eo(0),q2%eo(0)) + prod(q1%eo(1),q2%eo(1)) return end function
Function : | |
my_prod : | complex(DP) |
q1e : | type(field_quark_eo_wg), intent(in) |
q2e : | type(field_quark_eo_wg), intent(in) |
return inner-product : q1e’ .dot. q2e
function prod_fq_eo_wg( q1e, q2e ) result(my_prod) ! ! return inner-product : q1e' .dot. q2e ! implicit none type(field_quark_eo_wg), intent(in) :: q1e,q2e complex(DP) :: my_prod integer :: ix,iy,iz,ieoxyz,itb,ic,is #ifdef _DEBUG_ if (q1e%ieo /= q2e%ieo) then write(*,'("qf_eo_wg_prod: my_prod error in ieo")') stop endif #endif my_prod = Z0 !$OMP PARALLEL DO PRIVATE(ix,iy,iz,ieoxyz,itb,ic,is) REDUCTION(+:my_prod) do ix=1,NX do iy=1,NY do iz=1,NZ ieoxyz=mod(ipeo+ix+iy+iz+q1e%ieo,2) #ifdef _SF do itb=1,NTH-ieoxyz #else do itb=1-ieoxyz,NTH-ieoxyz #endif do is=1,SPIN do ic=1,COL my_prod = my_prod + conjg(q1e%s(itb,iz,iy,ix)%y(ic,is)) *q2e%s(itb,iz,iy,ix)%y(ic,is) enddo enddo enddo enddo enddo enddo #ifndef _singlePU call comlib_sumcast(my_prod) #endif return end function
Subroutine : | |
y : | type(field_quark_wg), intent(inout) |
set Gaussian noise on y
subroutine set_gaussian_noise_fq(y) ! ! set Gaussian noise on y ! implicit none type(field_quark_wg), intent(inout) :: y call set_gaussian_noise(y%eo(0)) call set_gaussian_noise(y%eo(1)) return end subroutine
Subroutine : | |
ye : | type(field_quark_eo_wg), intent(inout) |
set Gaussian noise on even/odd sites only
subroutine set_gaussian_noise_fq_eo(ye) ! ! set Gaussian noise on even/odd sites only ! implicit none type(field_quark_eo_wg), intent(inout) :: ye integer, parameter :: NCSTB=COL*SPIN*NTH real(DP) :: pi2 real(DP) :: th(NCSTB), yr(NCSTB) real(DP) :: cth(NCSTB),sth(NCSTB) integer :: ix,iy,iz,itb,ieoxyz,is,icstb integer :: ir #ifdef _PARALLELCHECK integer :: ic,itx,ity,itz,itt,ieo integer :: ipx,ipy,ipz,ics real(DP) :: pyr(CLSP) real(DP) :: pth(CLSP) real(DP) :: pcth(CLSP) real(DP) :: psth(CLSP) #endif pi2=PI*2 #ifdef _PARALLELCHECK !*********************************** !* Programmed independent on !* parallelization size. !*********************************** 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)) .and. (ieo == ye%ieo ) ) then do is=1,SPIN do ic=1,COL ics = ic + (is-1)*COL pth(ics) = pth(ics)*pi2 pyr(ics) = SQRT(-LOG(pyr(ics))) pcth(ics) = COS(pth(ics)) psth(ics) = SIN(pth(ics)) enddo enddo do is=1,SPIN do ic=1,COL ics = ic + (is-1)*COL ye%s(itb,iz,iy,ix)%y(ic,is)=pyr(ics)* cmplx(pcth(ics),psth(ics),kind=DP) enddo enddo endif enddo enddo enddo enddo #else do ix=1,NX do iy=1,NY do iz=1,NZ ieoxyz=mod(ipeo+ix+iy+iz+ye%ieo,2) call get(g_rand,yr) call get(g_rand,th) do icstb=1,NCSTB th(icstb) = th(icstb)*pi2 yr(icstb) = sqrt(-log(yr(icstb))) cth(icstb) = cos(th(icstb)) sth(icstb) = sin(th(icstb)) enddo ir=1 do itb=1-ieoxyz,NTH-ieoxyz do is = 1,SPIN ye%s(itb,iz,iy,ix)%y(1,is)=yr(ir )*cmplx(cth(ir ),sth(ir ),kind=DP) ye%s(itb,iz,iy,ix)%y(2,is)=yr(ir+1)*cmplx(cth(ir+1),sth(ir+1),kind=DP) ye%s(itb,iz,iy,ix)%y(3,is)=yr(ir+2)*cmplx(cth(ir+2),sth(ir+2),kind=DP) ir=ir+3 enddo enddo enddo enddo enddo ! end of do ixyz #endif #ifdef _SF call set_sf_boundary(ye) #endif return end subroutine
Subroutine : | |
y : | type(field_quark_wg), intent(inout) |
Set Schroedinger Functional boundary for fermion vector (set zero at itt=0,1,NTT+1) y : fermion field
subroutine set_sf_boundary_y(y) ! ! Set Schroedinger Functional boundary for fermion vector ! (set zero at itt=0,1,NTT+1) ! ! y : fermion field ! implicit none type(field_quark_wg), intent(inout) :: y call set_sf_boundary(y%eo(0)) call set_sf_boundary(y%eo(1)) return end subroutine
Subroutine : | |
ye : | type(field_quark_eo_wg), intent(inout) |
Set Schroedinger Functional boundary for fermion vector (set zero at itt=0,1,NTT+1) ye : fermion field (even/odd sites)
subroutine set_sf_boundary_y_eo(ye) ! ! Set Schroedinger Functional boundary for fermion vector ! (set zero at itt=0,1,NTT+1) ! ! ye : fermion field (even/odd sites) ! use sf_lattice_class implicit none type(field_quark_eo_wg), intent(inout) :: ye type(sf_lattice_world) :: sf_lattice integer :: ix,iy,iz,ieoxyz,ic,is,ieo if (.not.is_initialized(sf_lattice)) call new(sf_lattice) ieo = ye%ieo do ix=1,NX do iy=1,NY do iz=1,NZ ieoxyz=mod(ipeo+ix+iy+iz+ieo,2) !======================================================== ! SF boundary memo ! do itb=1-ieoxyz,NTH-ieoxyz ! it = 2*itb + ieoxyz ! ieoxyz = 0 => itb = 1..NTH => it = 2,4,6,..,NTT ! itb = 0 => it = 0 ! ! ieoxyz = 1 => itb = 0,..NTH-1 => it = 1,3,5,...,NTT-1 ! itb = NTH => it = NTT+1 ! itb = 0 => it = 1 !======================================================== ye%s( 0,iz,iy,ix)%y(:,:) = Z0 ! t = 0 or 1 if (1 == ieoxyz) then ye%s(NTH,iz,iy,ix)%y(:,:) = Z0 ! t = NTT+1 endif enddo enddo enddo return end subroutine
Subroutine : | |
v(COL*SPIN*NT*NZ*NY*NX) : | complex(DP), intent(in) |
y : | type(field_quark_wg), intent(inout) |
Unpack one-dimensional complex array to fermion field. User should set the even/odd-ness on the fermion field.
subroutine unpack_fq_wg(v,y) ! ! \Unpack one-dimensional complex array to fermion field. ! User should set the even/odd-ness on the fermion field. ! implicit none complex(DP), intent(in) :: v(COL*SPIN*NT*NZ*NY*NX) type(field_quark_wg), intent(inout) :: y integer, parameter :: NSITE = COL*SPIN*NT*NZ*NY*NX integer, parameter :: NHSITE = COL*SPIN*NTH*NZ*NY*NX call new(y) call unpack(v(1:NHSITE),y%eo(0)) call unpack(v(NHSITE+1:NSITE),y%eo(1)) return end subroutine
Subroutine : | |||
vec(COL*SPIN*NTH*NZ*NY*NX) : | complex(DP), intent(in) | ||
yeo : | type(field_quark_eo_wg),
intent(inout)
|
Unpack one-dimensional complex array to fermion field. User should set the even/odd-ness on the fermion field.
subroutine unpack_fq_eo_wg(vec,yeo) ! ! \Unpack one-dimensional complex array to fermion field. ! User should set the even/odd-ness on the fermion field. ! implicit none complex(DP), intent(in) :: vec(COL*SPIN*NTH*NZ*NY*NX) type(field_quark_eo_wg), intent(inout) :: yeo ! even/odd fermion field. one shout set even/odd-ness. integer :: ix,iy,iz,itb,ic,is,ieoxyz,isite integer :: ix_base integer :: iy_base integer :: iz_base integer :: itb_base integer :: is_base !$OMP PARALLEL DO & !$OMP PRIVATE(ix,iy,iz,ieoxyz,itb,ic,is, & !$OMP ix_base,iy_base,iz_base,itb_base,is_base,isite) do ix=1,NX ix_base=(ix -1)*COL*SPIN*NTH*NZ*NY do iy=1,NY iy_base=(iy -1)*COL*SPIN*NTH*NZ +ix_base do iz=1,NZ iz_base=(iz -1)*COL*SPIN*NTH +iy_base ieoxyz=mod(ipeo+ix+iy+iz+yeo%ieo,2) #ifdef _SF if (1==ieoxyz) then !=================================== ! itt=1 is SF boundary ! set zero on fermion field !=================================== itb = 0 itb_base=(itb+ieoxyz-1)*COL*SPIN +iz_base do is=1,SPIN is_base=(is -1)*COL +itb_base do ic=1,COL isite= ic+is_base yeo%s(itb,iz,iy,ix)%y(ic,is) = Z0 enddo enddo endif do itb=1,NTH-ieoxyz itb_base=(itb+ieoxyz-1)*COL*SPIN +iz_base do is=1,SPIN is_base=(is -1)*COL +itb_base do ic=1,COL isite= ic+is_base yeo%s(itb,iz,iy,ix)%y(ic,is)=vec(isite) enddo enddo enddo #else do itb=1-ieoxyz,NTH-ieoxyz itb_base=(itb+ieoxyz-1)*COL*SPIN +iz_base do is=1,SPIN is_base=(is -1)*COL +itb_base do ic=1,COL isite= ic+is_base yeo%s(itb,iz,iy,ix)%y(ic,is)=vec(isite) enddo enddo enddo #endif enddo enddo enddo return end subroutine