Class | quark_clover_class |
In: |
QuarkCloverClass/quark_clover_class.F90
|
$Id: quark_clover_class.F90,v 1.25 2011/09/25 06:13:11 ishikawa Exp $
Subroutine : | |||
ye : | type(field_quark_eo_wg), intent(inout)
| ||
this : | type(quark_clover),
intent(in)
|
Multiply inverse clover term matrix on fermion field (overwritten) (even/odd sites only)
ye <= Fee \ ye
subroutine accum_mult_clover_term_eo(ye,this) ! ! Multiply inverse clover term matrix on fermion field (overwritten) (even/odd sites only) ! ! ye <= Fee \ ye ! implicit none type(field_quark_eo_wg), intent(inout) :: ye ! even/odd site fermion field (overwritten) type(quark_clover), intent(in) :: this ! clover quark paramters complex(DP):: a1,a2,a3,a4,a5,a6 complex(DP):: b1,b2,b3,b4,b5,b6 complex(DP):: p1,p2,p3,p4,p5,p6 complex(DP):: m1,m2,m3,m4,m5,m6 complex(DP):: y1,y2,y3,y4,y5,y6,y7,y8,y9,ya,yb,yc integer :: ix,iy,iz,itb,ieoxyz,itb0 integer :: ieo ieo = ye%ieo !$OMP PARALLEL DO & !$OMP PRIVATE(ix,iy,iz,ieoxyz,itb,itb0, & !$OMP y1,y2,y3,y4,y5,y6,y7,y8,y9,ya,yb,yc, & !$OMP p1,p2,p3,p4,p5,p6,m1,m2,m3,m4,m5,m6, & !$OMP a1,a2,a3,a4,a5,a6,b1,b2,b3,b4,b5,b6) do ix=1,NX do iy=1,NY do iz=1,NZ ieoxyz=mod(ipeo+ix+iy+iz+ieo,2) #ifdef _SF !=============================== ! Set SF boundary condition ! itt = 0,1 and NTT+1 is zero !=============================== ye%s(0,iz,iy,ix)%y(:,:) = Z0 do itb=1,NTH-ieoxyz #else do itb=1-ieoxyz,NTH-ieoxyz #endif itb0 = itb+ieoxyz #define _YE_ ye%s(itb,iz,iy,ix) #define _YDE_ ye%s(itb,iz,iy,ix) #define _FCLV_ this%fclinv #include "mult_clover_term_eo_kernel.h90" enddo #ifdef _SF if (ieoxyz==1) then ye%s(NTH,iz,iy,ix)%y(:,:) = Z0 endif #endif enddo enddo enddo ! end of do ix,iy,iz return end subroutine
Subroutine : | |
qp1 : | type(quark_clover), intent(inout) |
qp2 : | type(quark_clover), intent(in) |
subroutine assign_quark_clover(qp1,qp2) implicit none type(quark_clover), intent(inout) :: qp1 type(quark_clover), intent(in) :: qp2 call assign(qp1%quark_wilson,qp2%quark_wilson) qp1%csw = qp2%csw qp1%fclinv => qp2%fclinv qp1%clover_term_log = qp2%clover_term_log qp1%clover_term_log_fname = qp2%clover_term_log_fname return end subroutine
Subroutine : | |||
yde : | type(field_quark_eo_wg), intent(inout)
| ||
this : | type(quark_clover),
intent(in)
| ||
ye : | type(field_quark_eo_wg), intent(in)
|
Multiply inverse clover term matrix on fermion field (even/odd sites only)
yde <= Fee \ ye
subroutine assign_mult_clover_term_eo(yde,this,ye) ! ! Multiply inverse clover term matrix on fermion field (even/odd sites only) ! ! yde <= Fee \ ye ! implicit none type(field_quark_eo_wg), intent(inout) :: yde ! even/odd site fermion field type(quark_clover), intent(in) :: this ! clover quark parameters type(field_quark_eo_wg), intent(in) :: ye ! even/odd site fermion field complex(DP):: a1,a2,a3,a4,a5,a6 complex(DP):: b1,b2,b3,b4,b5,b6 complex(DP):: p1,p2,p3,p4,p5,p6 complex(DP):: m1,m2,m3,m4,m5,m6 complex(DP):: y1,y2,y3,y4,y5,y6,y7,y8,y9,ya,yb,yc integer :: ix,iy,iz,itb,ieoxyz,ieo,itb0 ieo = ye%ieo yde%ieo = ye%ieo !*********************************************************** !* Calc: !* !* yde = F^-1ee ye !* !*********************************************************** !$OMP PARALLEL DO & !$OMP PRIVATE(ix,iy,iz,ieoxyz,itb,itb0, & !$OMP a1,a2,a3,a4,a5,a6,b1,b2,b3,b4,b5,b6, & !$OMP p1,p2,p3,p4,p5,p6,m1,m2,m3,m4,m5,m6, & !$OMP y1,y2,y3,y4,y5,y6,y7,y8,y9,ya,yb,yc) do ix=1,NX do iy=1,NY do iz=1,NZ ieoxyz=mod(ipeo+ix+iy+iz+yde%ieo,2) #ifdef _SF !=============================== ! Set SF boundary condition ! itt = 0,1 and NTT+1 is zero !=============================== yde%s(0,iz,iy,ix)%y(:,:) = Z0 do itb=1,NTH-ieoxyz #else do itb=1-ieoxyz,NTH-ieoxyz #endif itb0 = itb+ieoxyz #define _YE_ ye%s(itb,iz,iy,ix) #define _YDE_ yde%s(itb,iz,iy,ix) #define _FCLV_ this%fclinv #include "mult_clover_term_eo_kernel.h90" enddo #ifdef _SF if (ieoxyz==1) then yde%s(NTH,iz,iy,ix)%y(:,:) = Z0 endif #endif enddo enddo enddo ! end of do ix,iy,iz return end subroutine
Subroutine : | |||
this : | type(quark_clover),
intent(in)
| ||
yde : | type(field_quark_eo_wg), intent(inout)
| ||
ye : | type(field_quark_eo_wg), intent(inout)
| ||
u : | type(vfield_gluon_wg), intent(in)
|
Multiply even/odd-site and clover term preconditioned Wilson-Clover-Dirac operator(even/odd sites only):
v yde <= Dee ye, where v Dee = 1ee - kappa^2 Fee \ Meo Foo \ Moe
subroutine assign_mult_eoprec_cd_eo(this,yde,ye,u) ! ! Multiply even/odd-site and clover term preconditioned Wilson-Clover-Dirac operator(even/odd sites only): ! ! v ! yde <= Dee ye, ! ! where v ! Dee = 1ee - kappa^2 Fee \ Meo Foo \ Moe ! implicit none type(quark_clover), intent(in) :: this ! clover quark parameters and inverse clover term type(field_quark_eo_wg), intent(inout) :: yde ! even/odd site fermion field type(field_quark_eo_wg), intent(inout) :: ye ! even/odd site fermion field type(vfield_gluon_wg), intent(in) :: u ! gauge field type(field_quark_eo_wg), allocatable :: Myo real(DP) :: kappa2,kappa allocate(Myo) kappa = get_kappa(this%quark_wilson) kappa2=-(kappa**2) !*********************************** !* Myo = F^-1oo Moe ye !* yde = F^-1ee Meo Myo = F^-1ee Meo F^-1oo Moe ye !*********************************** call assign_mult_hop(this,Myo, ye,u) call assign_mult_hop(this,yde,Myo,u) !*********************************** !* yde = ye - kappa^2 yde !* = ye - kappa^2 Meo Moe ye !*********************************** call accum_mult_add(yde,kappa2,ye) deallocate(Myo) return end subroutine
Subroutine : | |||
this : | type(quark_clover),
intent(in)
| ||
FMye : | type(field_quark_eo_wg), intent(out)
| ||
yo : | type(field_quark_eo_wg), intent(inout)
| ||
u : | type(vfield_gluon_wg), intent(in)
|
Multiply clover term preconditioned Hopping operator (odd->even or even->odd hopping only)
FMye <= Fee \ Meo yo
subroutine assign_mult_clover_hopping_eo(this,FMye,yo,u) ! ! Multiply clover term preconditioned Hopping operator (odd->even or even->odd hopping only) ! ! FMye <= Fee \ Meo yo ! use counter_class implicit none type(quark_clover), intent(in) :: this ! clover quark parameter and inverse clover term type(field_quark_eo_wg), intent(out) :: FMye ! even/odd site fermion field type(field_quark_eo_wg), intent(inout) :: yo ! odd/even site fermion field type(vfield_gluon_wg), intent(in) :: u ! gauge field call copy_boundary(yo) call mult_clv_hopping_tzyx_eo(FMye,yo,u,this%fclinv) call inc(mult_iter) ! call assign_mult_hop(this%quark_wilson,FMye,yo,u) ! call accum_mult_clover_term(FMye,this) return end subroutine
Subroutine : | |
this : | type(quark_clover), intent(inout) |
subroutine delete_quark_clover(this) implicit none type(quark_clover), intent(inout) :: this call delete(this%quark_wilson) call delete(this%clover_term_log) return end subroutine
Subroutine : | |
this : | type(quark_clover), intent(inout) |
delete inverse clover term.
subroutine delete_inverse_clover_term(this) ! ! delete inverse clover term. ! implicit none type(quark_clover), intent(inout) :: this if (associated(this%fclinv)) then deallocate(this%fclinv) endif return end subroutine
Subroutine : | |||
BB : | type(vfield_gluon_wg), intent(inout)
| ||
XX : | type(tfield_gluon_wg), intent(inout)
| ||
u : | type(vfield_gluon_wg), intent(in)
| ||
wl34 : | type(tfield_gluon_wg), intent(in)
| ||
wl98 : | type(tfield_gluon_wg), intent(in)
|
Compute Force for clover term contributions.
BB_{mu}(n) = Sum[vcl_{mu,nu}(n),nu=!mu] X4---<---3X Xn 2X | | | | vcl = v ^ - ^ v | | | | Xn 2X X6---<---5X X : X_{mu,nu}(n) - X^+_{mu,nu}(n), where X_{mu,nu}(n) = \sigma_{mu,nu} (force contributions at site n)
subroutine force_clover_terms(BB,XX,u,wl34,wl98) ! ! Compute Force for clover term contributions. ! ! BB_{mu}(n) = Sum[vcl_{mu,nu}(n),nu=!mu] ! ! X4---<---3X Xn 2X ! | | | | ! vcl = v ^ - ^ v ! | | | | ! Xn 2X X6---<---5X ! ! X : X_{mu,nu}(n) - X^+_{mu,nu}(n), ! where ! X_{mu,nu}(n) = \sigma_{mu,nu} (force contributions at site n) ! implicit none type(vfield_gluon_wg), intent(inout) :: BB ! pre-force contirubution (\dot{u}) accumulated. type(tfield_gluon_wg), intent(inout) :: XX ! pre-force contirubution (\dot{\sigma_{mu,nu}F_{mu,nu}) accumulated. type(vfield_gluon_wg), intent(in) :: u ! gauge field type(tfield_gluon_wg), intent(in) :: wl34,wl98 ! 2-link fields prepared by make_two_links integer :: ipl do ipl=1,NDIM*(NDIM-1)/2 call copy_boundary(XX%eo(0)%munu(ipl)) call copy_boundary(XX%eo(1)%munu(ipl)) enddo call force_clover_terms_eo(BB%eo(0),XX%eo(0),XX%eo(1),u%eo(0),u%eo(1),wl34%eo(0),wl98%eo(1)) call force_clover_terms_eo(BB%eo(1),XX%eo(1),XX%eo(0),u%eo(1),u%eo(0),wl34%eo(1),wl98%eo(0)) return end subroutine
Subroutine : | |||
this : | type(quark_clover),
intent(in)
| ||
XX : | type(tfield_gluon_wg), intent(inout)
|
Calc force contiribution from log determinant of clover term derived from even/odd site preconditioning.
Action:
SF = - Nf log[det[F]]
subroutine force_clover_trlog(this,XX) ! ! Calc force contiribution from log determinant of clover term derived ! from even/odd site preconditioning. ! ! Action: ! ! SF = - Nf log[det[F]] ! implicit none type(quark_clover), intent(in) :: this ! quark parameters type(tfield_gluon_wg), intent(inout) :: XX ! pre-force contribution (\dot{\sigma_{mu,nu}F_{mu,nu}}) call force_clover_trlog_eo(this,this%fclinv%eo(0),XX%eo(0)) call force_clover_trlog_eo(this,this%fclinv%eo(1),XX%eo(1)) return end subroutine
Subroutine : | |||
XX : | type(tfield_gluon_wg), intent(inout)
| ||
zickd8 : | complex(DP), intent(in)
| ||
fx : | type(field_quark_wg), intent(in)
| ||
fz : | type(field_quark_wg), intent(in)
|
Calc force from inverse clover term in even/odd preconditioned Wilson-Dirac(Clover) hopping operator
XX(n) = XX(n) + zickd8 * tr[ sigma_{mu,nu} fz(n) fx(n)^{+}]
where, sigma_{mu,nu} : Dirac rep.
subroutine force_hmc_clover(XX,zickd8,fx,fz) ! ! Calc force from inverse clover term in even/odd preconditioned Wilson-Dirac(Clover) hopping operator ! ! XX(n) = XX(n) + zickd8 * tr[ sigma_{mu,nu} fz(n) fx(n)^{+}] ! ! where, sigma_{mu,nu} : Dirac rep. ! implicit none type(tfield_gluon_wg), intent(inout) :: XX ! pre-force contribution for clvoer term (\dot{T}) complex(DP), intent(in) :: zickd8 ! force coefficient type(field_quark_wg), intent(in) :: fx ! pseudo fermion type(field_quark_wg), intent(in) :: fz ! pseudo fermion call force_hmc_clover_eo(XX%eo(0),zickd8,fx%eo(0),fz%eo(0)) call force_hmc_clover_eo(XX%eo(1),zickd8,fx%eo(1),fz%eo(1)) return end subroutine
Function : | |
csw : | real(DP) |
this : | type(quark_clover), intent(in) |
function get_csw(this) result(csw) implicit none type(quark_clover), intent(in) :: this real(DP) :: csw csw = this%csw return end function
Function : | |
id : | integer |
this : | class(action_parameters), intent(in) |
Original external subprogram is quark_wilson_class#get_id
Function : | |
kappa : | real(DP) |
this : | type(quark_clover), intent(in) |
function get_kappa_quark_clover(this) result(kappa) implicit none type(quark_clover), intent(in) :: this real(DP) :: kappa kappa = get_kappa(this%quark_wilson) return end function
Function : | |
nf : | integer |
this : | type(quark_clover), intent(in) |
function get_nflavor_quark_clover(this) result(nf) implicit none type(quark_clover), intent(in) :: this integer :: nf nf = get_nflavor(this%quark_wilson) return end function
Function : | |
ptr : | type(field_clover_term), pointer |
this : | type(quark_clover), intent(in) |
function get_ptr_inverse_clover_term(this) result(ptr) implicit none type(quark_clover), intent(in) :: this type(field_clover_term), pointer :: ptr ptr => this%fclinv return end function
Subroutine : | |||
this : | type(quark_clover),
intent(inout)
| ||
action : | real(DP), intent(out)
|
Calculate Hamiltonian of log determinant clover term. This also computes inverse clover term by side effect.
Action : SF = - Nf log[det[F]]
subroutine hamil_clover_trlog(this,action) ! ! Calculate Hamiltonian of log determinant clover term. ! This also computes inverse clover term by side effect. ! ! Action : ! ! SF = - Nf log[det[F]] ! implicit none type(quark_clover), intent(inout) :: this ! quark parameters, and contains ! inverse clover term (in chiral rep., linear form) real(DP), intent(out) :: action ! action value real(DP) :: logdetfcl !************************************************** !* calc. inverse clover term and log[det[Fcl]] !************************************************** call make_inverse_clover_term(this,1,logdetfcl) action=-logdetfcl*get_nflavor(this%quark_wilson) return end subroutine
Subroutine : | |||
u : | type(vfield_gluon_wg), intent(in)
| ||
wl34 : | type(tfield_gluon_wg), intent(in)
| ||
wl98 : | type(tfield_gluon_wg), intent(in)
|
Calculate clover leaf (field strength) for Clover fermion
clover field strength is class/module variable and is calculated using the two-link fields.
subroutine make_clover_leaf(u,wl34,wl98) ! ! Calculate clover leaf (field strength) for Clover fermion ! ! clover field strength is class/module variable and is calculated using the two-link fields. ! implicit none type(vfield_gluon_wg), intent(in) :: u ! gauge link type(tfield_gluon_wg), intent(in) :: wl34,wl98 ! 2-links prepared by make_two_links integer :: iflag iflag = 0 if (.not.associated(m_ucl)) then allocate(m_ucl) endif call make_clover_leaf_eo(m_ucl%eo(0),u%eo(0),u%eo(1),wl34%eo(0),wl34%eo(1),wl98%eo(0),wl98%eo(1)) call make_clover_leaf_eo(m_ucl%eo(1),u%eo(1),u%eo(0),wl34%eo(1),wl34%eo(0),wl98%eo(1),wl98%eo(0)) return end subroutine
Subroutine : | |||
this : | type(quark_clover),
intent(inout)
| ||
idet : | integer, intent(in)
| ||
logdetfcl : | real(DP), intent(out)
|
Calculate inverse clover term matrix and log-determinant of clover term.
fcl = [ 1 - csw kappa/2 sigma_{mu,nu} F_{mu,nu}] fclinv = fcl^-1 log[det[fcl]]
This uses module variable, m_ucl : clover leaf field.
[ 1 - csw kappa/2 sigma_{mu,nu} F_{mu,nu}]^-1 = diag(fclinv(1),fclinv(2)), where | 1 2 3 4 5 6 | | 7 8 9 10 11 | | 12 13 14 15 | | 16 17 18 | * 0.5 | 19 20 | | 21 | are stored in fclinv(i). Note: 0.5 is multiplied to simplify conversion from chiral to Dirac rep.
subroutine make_inverse_clover_term(this,idet,logdetfcl) ! ! Calculate inverse clover term matrix and log-determinant of clover term. ! ! fcl = [ 1 - csw kappa/2 sigma_{mu,nu} F_{mu,nu}] ! fclinv = fcl^-1 ! log[det[fcl]] ! ! This uses module variable, m_ucl : clover leaf field. ! ! [ 1 - csw kappa/2 sigma_{mu,nu} F_{mu,nu}]^-1 = diag(fclinv(1),fclinv(2)), ! ! where ! ! | 1 2 3 4 5 6 | ! | 7 8 9 10 11 | ! | 12 13 14 15 | ! | 16 17 18 | * 0.5 ! | 19 20 | ! | 21 | ! ! are stored in fclinv(i). ! ! Note: 0.5 is multiplied to simplify conversion from chiral to Dirac rep. ! use comlib implicit none type(quark_clover), intent(inout) :: this ! quark parameter contains ! inverse clover term (in chiral rep. linear storage) integer, intent(in) :: idet ! if idet == 1 then logdetfcl is computed. real(DP), intent(out) :: logdetfcl ! log determinat of clover term. type(f_c_m_eo), allocatable :: f(:),finv(:) real(DP) :: logdetf1cl(0:1) real(DP) :: logdetf2cl(0:1) integer :: ieo,iout real(DP) :: kappa,csw iout = get_file_unit(this%clover_term_log) kappa = get_kappa(this%quark_wilson) csw = this%csw logdetfcl=0.0_DP if (.not.associated(m_ucl)) then call error_stop("Clover leaf field is not created, Stop!") endif if (.not.associated(this%fclinv)) then allocate(this%fclinv) endif call new(this%fclinv) allocate(f(2),finv(2)) do ieo=0,1 !***************************** !* calc clover term matrix in chiral rep. !* !* fcl = diag[f1cl,f2cl] = 1 - csw kappa/2 sigma_{mu,nu}F_{mu,nu} !* !* solve inverse matrix (even part) !* !***************************** call clover_f1f2(kappa,csw,f(1),f(2),m_ucl%eo(ieo)) call _CLV_SOLVER(f(1)%m,iout,idet,finv(1)%m,logdetf1cl(ieo)) call _CLV_SOLVER(f(2)%m,iout,idet,finv(2)%m,logdetf2cl(ieo)) finv(1)%ieo=f(1)%ieo finv(2)%ieo=f(2)%ieo !***************************** !* chiral full matrix -> chiral linear vector !***************************** call conv_matrix_to_linear_eo(finv(1),finv(2),this%fclinv%eo(ieo)) enddo if (idet == 1) then logdetfcl=logdetf1cl(0)+logdetf2cl(0) +logdetf1cl(1)+logdetf2cl(1) #ifndef _singlePU call comlib_sumcast(logdetfcl) #endif endif deallocate(f,finv) return end subroutine
Subroutine : | |
this : | type(quark_clover), intent(inout) |
id : | integer, intent(in) |
subroutine new_quark_clover(this,id) implicit none type(quark_clover), intent(inout) :: this integer, intent(in) :: id character(len=CHARLEN) :: str call new(this%quark_wilson,id) this%action_name = ACTION_NAME write(str,'(I2.2)')id this%clover_term_log_fname = TRIM(this%clover_term_log_fname)//TRIM(str) call new(this%clover_term_log,this%clover_term_log_fname) return end subroutine
Subroutine : | |
this : | type(quark_clover), intent(inout) |
subroutine print_quark_clover(this) implicit none type(quark_clover), intent(inout) :: this call print(this%quark_wilson) write(*,'(14X," Csw :",F12.8)')this%csw return end subroutine
Derived Type : |
Wilson quark type extended from the above type (extends from action_parameters)
Original external subprogram is quark_wilson_class#quark_wilson
Subroutine : | |
this : | type(quark_clover), intent(inout) |
iout : | integer, intent(in) |
subroutine read_quark_clover(this,iout) use comlib implicit none type(quark_clover), intent(inout) :: this integer, intent(in) :: iout call read(this%quark_wilson,iout) if (0 == nodeid) then read(iout,*)this%csw endif #ifndef _singlePU call comlib_bcast(this%csw,0) #endif return end subroutine
Subroutine : | |
this : | type(quark_clover), intent(in) |
iout : | integer, intent(in) |
subroutine save_config_quark_clover(this,iout) implicit none type(quark_clover), intent(in) :: this integer, intent(in) :: iout call save_config(this%quark_wilson,iout) write(iout)this%csw return end subroutine
Subroutine : | |
this : | type(quark_clover), intent(inout) |
csw : | real(DP), intent(in) |
subroutine set_csw(this,csw) implicit none type(quark_clover), intent(inout) :: this real(DP), intent(in) :: csw this%csw = csw 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)
Original external subprogram is field_gauge_class#set_gaussian_noise
Subroutine : | |
y : | type(field_quark_wg), intent(inout) |
set Gaussian noise on y
Original external subprogram is field_fermion_class#set_gaussian_noise
Subroutine : | |
ye : | type(field_quark_eo_wg), intent(inout) |
set Gaussian noise on even/odd sites only
Original external subprogram is field_fermion_class#set_gaussian_noise
Subroutine : | |
this : | class(action_parameters), intent(inout) |
id : | integer, intent(in) |
Original external subprogram is quark_wilson_class#set_id
Subroutine : | |
this : | type(quark_clover), intent(inout) |
kappa : | real(DP), intent(in) |
subroutine set_kappa_quark_clover(this,kappa) implicit none type(quark_clover), intent(inout) :: this real(DP), intent(in) :: kappa call set_kappa(this%quark_wilson,kappa) return end subroutine
Subroutine : | |
this : | type(quark_clover), intent(inout) |
nf : | integer, intent(in) |
subroutine set_nflavor_quark_clover(this,nf) implicit none type(quark_clover), intent(inout) :: this integer, intent(in) :: nf call set_nflavor(this%quark_wilson,nf) return end subroutine