diff --git a/docs/src/dirac.md b/docs/src/dirac.md index c664063..75a330a 100644 --- a/docs/src/dirac.md +++ b/docs/src/dirac.md @@ -28,7 +28,7 @@ Wilson-Dirac operator. The action of the Dirac operator `Dw!` is the following: ```math -D_w\psi (\vec{x} = x_1,x_2,x_3,x_4) = (4 + m_0)\psi(\vec{x}) - +D_w\psi (\vec{x} = x_1,x_2,x_3,x_4) = (4 + m_0 + i \mu \gamma_5)\psi(\vec{x}) - ``` ```math - \frac{1}{2}\sum_{\mu = 1}^4 \theta (\mu) (1-\gamma_\mu) U_\mu(\vec{x}) \psi(\vec{x} + \hat{\mu}) + \theta^* (\mu) (1 + \gamma_\mu) U^{-1}_\mu(\vec{x} - \hat{\mu}) \psi(\vec{x} - \hat{\mu}) @@ -108,4 +108,5 @@ g5Dw! DwdagDw! Csw! pfrandomize! +mtwmdpar ``` diff --git a/src/Dirac/Dirac.jl b/src/Dirac/Dirac.jl index a6ff103..74eb63c 100644 --- a/src/Dirac/Dirac.jl +++ b/src/Dirac/Dirac.jl @@ -22,7 +22,7 @@ using ..Spinors """ struct DiracParam{T,R} -Stores the parameters of the Dirac operator. It can be generated via the constructor `function DiracParam{T}(::Type{R},m0,csw,th,ct)`. The first argument can be ommited and is taken to be `SU3fund`. +Stores the parameters of the Dirac operator. It can be generated via the constructor `function DiracParam{T}(::Type{R},m0,csw,th,tm,ct)`. The first argument can be ommited and is taken to be `SU3fund`. The parameters are: - `m0::T` : Mass of the fermion @@ -537,7 +537,7 @@ function DwdagDw!(so, U, si, dpar::DiracParam, dws::DiracWorkspace, lp::Union{Sp SF_bndfix!(dws.st,lp) @timeit "g5Dw" begin CUDA.@sync begin - CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_g5Dwimpr!(so, U, dws.st, dws.csw, dpar.m0, dpar.tm, dpar.th, dpar.csw, dpar.ct, lp) + CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_g5Dwimpr!(so, U, dws.st, dws.csw, dpar.m0, -dpar.tm, dpar.th, dpar.csw, dpar.ct, lp) end end SF_bndfix!(so,lp) @@ -553,7 +553,7 @@ function DwdagDw!(so, U, si, dpar::DiracParam, dws::DiracWorkspace, lp::Union{Sp SF_bndfix!(dws.st,lp) @timeit "g5Dw" begin CUDA.@sync begin - CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_g5Dw!(so, U, dws.st, dpar.m0, dpar.tm, dpar.th, dpar.ct, lp) + CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_g5Dw!(so, U, dws.st, dpar.m0, -dpar.tm, dpar.th, dpar.ct, lp) end end SF_bndfix!(so,lp) @@ -576,7 +576,7 @@ function DwdagDw!(so, U, si, dpar::DiracParam, dws::DiracWorkspace, lp::SpacePar @timeit "g5Dw" begin CUDA.@sync begin - CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_g5Dwimpr!(so, U, dws.st, dws.csw, dpar.m0, dpar.tm, dpar.th, dpar.csw, lp) + CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_g5Dwimpr!(so, U, dws.st, dws.csw, dpar.m0, -dpar.tm, dpar.th, dpar.csw, lp) end end end @@ -591,7 +591,7 @@ function DwdagDw!(so, U, si, dpar::DiracParam, dws::DiracWorkspace, lp::SpacePar @timeit "g5Dw" begin CUDA.@sync begin - CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_g5Dw!(so, U, dws.st, dpar.m0, dpar.tm, dpar.th, lp) + CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_g5Dw!(so, U, dws.st, dpar.m0, -dpar.tm, dpar.th, lp) end end end @@ -600,6 +600,16 @@ function DwdagDw!(so, U, si, dpar::DiracParam, dws::DiracWorkspace, lp::SpacePar return nothing end +""" + function mtwmdpar(dpar::DiracParam) + +Returns `dpar` with oposite value of the twisted mass. +""" +function mtwmdpar(dpar::DiracParam{P,R}) where {P,R} + return DiracParam{P}(R,dpar.m0,dpar.csw,dpar.th,-dpar.tm,dpar.ct) +end + + """ SF_bndfix!(sp, lp::Union{SpaceParm{4,6,BC_SF_ORBI,D},SpaceParm{4,6,BC_SF_AFWB,D}}) @@ -694,7 +704,7 @@ function krnl_assign_pf_su2!(f::AbstractArray, p , lp::SpaceParm, t::Int64) return nothing end -export Dw!, g5Dw!, DwdagDw!, SF_bndfix!, Csw!, pfrandomize! +export Dw!, g5Dw!, DwdagDw!, SF_bndfix!, Csw!, pfrandomize!, mtwmdpar include("DiracIO.jl") diff --git a/src/LatticeGPU.jl b/src/LatticeGPU.jl index c091ac3..acbcc82 100644 --- a/src/LatticeGPU.jl +++ b/src/LatticeGPU.jl @@ -57,7 +57,7 @@ export pmul, gpmul, gdagpmul, dmul include("Dirac/Dirac.jl") using .Dirac export DiracWorkspace, DiracParam -export Dw!, g5Dw!, DwdagDw!, SF_bndfix!, Csw!, pfrandomize! +export Dw!, g5Dw!, DwdagDw!, SF_bndfix!, Csw!, pfrandomize!, mtwmdpar export read_prop, save_prop, read_dpar include("Solvers/Solvers.jl") diff --git a/src/Solvers/Propagators.jl b/src/Solvers/Propagators.jl index 99cf462..58008b2 100644 --- a/src/Solvers/Propagators.jl +++ b/src/Solvers/Propagators.jl @@ -25,7 +25,7 @@ function propagator!(pro, U, dpar::DiracParam{T}, dws::DiracWorkspace, lp::Space CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnlg5!(dws.sp) end - g5Dw!(pro,U,dws.sp,dpar,dws,lp) + g5Dw!(pro,U,dws.sp,mtwmdpar(dpar),dws,lp) CG!(pro,U,DwdagDw!,dpar,lp,dws,maxiter,tol) return nothing @@ -46,7 +46,7 @@ function propagator!(pro, U, dpar::DiracParam{T}, dws::DiracWorkspace, lp::Space CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnlg5!(dws.sp) end - g5Dw!(pro,U,dws.sp,dpar,dws,lp) + g5Dw!(pro,U,dws.sp,mtwmdpar(dpar),dws,lp) CG!(pro,U,DwdagDw!,dpar,lp,dws,maxiter,tol) return nothing @@ -92,7 +92,7 @@ function bndpropagator!(pro, U, dpar::DiracParam{T}, dws::DiracWorkspace, lp::Sp CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnlg5!(dws.sp) end - g5Dw!(pro,U,dpar.ct*dws.sp,dpar,dws,lp) + g5Dw!(pro,U,dpar.ct*dws.sp,mtwmdpar(dpar),dws,lp) CG!(pro,U,DwdagDw!,dpar,lp,dws,maxiter,tol) return nothing @@ -137,7 +137,7 @@ function Tbndpropagator!(pro, U, dpar::DiracParam{T}, dws::DiracWorkspace, lp::S CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnlg5!(dws.sp) end - g5Dw!(pro,U,dpar.ct*dws.sp,dpar,dws,lp) + g5Dw!(pro,U,dpar.ct*dws.sp,mtwmdpar(dpar),dws,lp) CG!(pro,U,DwdagDw!,dpar,lp,dws,maxiter,tol) return nothing diff --git a/test/dirac/test_solver_rand.jl b/test/dirac/test_solver_rand.jl index aaad348..0714d8f 100644 --- a/test/dirac/test_solver_rand.jl +++ b/test/dirac/test_solver_rand.jl @@ -32,7 +32,8 @@ end CUDA.@sync begin CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnlg5!(rpsi) end -g5Dw!(prop,U,rpsi,dpar,dws,lp) + +g5Dw!(prop,U,rpsi,mtwmdpar(dpar),dws,lp) CG!(prop,U,DwdagDw!,dpar,lp,dws,10000,1.0e-14) Dw!(dws.sp,U,prop,dpar,dws,lp) diff --git a/test/runtests.jl b/test/runtests.jl index fe59102..2f68f70 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,3 +1,6 @@ -#include("SAD/test_sad.jl") +include("SAD/test_sad.jl") include("flow/test_adapt.jl") +include("dirac/test_fp_fa.jl") +include("dirac/test_solver_plw.jl") +include("dirac/test_solver_rand.jl")