From b9e61565551715358f4eee4b27bc5d7cee62e62a Mon Sep 17 00:00:00 2001
From: James Osborn <osborn@alcf.anl.gov>
Date: Thu, 19 Dec 2024 12:34:05 -0600
Subject: [PATCH] improve hisqsmear tests

---
 src/gauge/hisqsmear.nim   | 50 +++++++++++++++++++++++++--------------
 src/maths/complexType.nim |  2 +-
 2 files changed, 33 insertions(+), 19 deletions(-)

diff --git a/src/gauge/hisqsmear.nim b/src/gauge/hisqsmear.nim
index 77b7749..a7bd4a5 100644
--- a/src/gauge/hisqsmear.nim
+++ b/src/gauge/hisqsmear.nim
@@ -65,6 +65,8 @@ when isMainModule:
     ll2 = lo.newGauge()
     fc = lo.newGauge()
     lc = lo.newGauge()
+    ch = lo.newGauge()
+    chl = lo.newGauge()
     gc = GaugeActionCoeffs(plaq:1.0)
     eps = floatParam("eps", 1e-6)
     #warm = floatParam("warm", 1e-5)
@@ -76,6 +78,8 @@ when isMainModule:
   #g.warm warm, r
   g.stagPhase
   dg.gaussian r
+  ch.gaussian r
+  chl.gaussian r
   for mu in 0..<g.len:
     for e in g[mu]:
       dg[mu][e] *= eps
@@ -91,30 +95,40 @@ when isMainModule:
       echo &"> ERROR rel error |{r}| > {tol*eps}"
 
   proc checkG =
-    echo "Checking GaugeDeriv"
+    #echo "Checking GaugeDeriv"
+    #for mu in 0..<fd.len:
+    #  fd[mu] := 0
+    #let a = gc.gaugeAction2(g)
+    #let a2 = gc.gaugeAction2(g2)
+    #gc.gaugeDeriv2(g, fd)
+    echo "Checking redot deriv"
+    var a, a2 = 0.0
     for mu in 0..<fd.len:
-      fd[mu] := 0
-    let a = gc.gaugeAction2(g)
-    let a2 = gc.gaugeAction2(g2)
-    gc.gaugeDeriv2(g, fd)
-    check(a2-a, 10)
+      a += redot(ch[mu], g[mu])
+      a2 += redot(ch[mu], g2[mu])
+      fd[mu] := ch[mu]
+    check(a2-a, 1)
   checkG()
 
   proc checkL(name: string, tol: float) =
     echo "Checking ", name
-    for mu in 0..<fd.len:
-      fd[mu] := 0
-      fc[mu] := 0
-      lc[mu] := 0
     resetTimers()
     let f = hisq.smearGetForce(g, fl, ll)
     let f2 = hisq.smearGetForce(g2, fl2, ll2)
-    gc.gaugeDeriv2(fl, fc)
-    gc.gaugeDeriv2(ll, lc)
+    var a, a2 = 0.0
+    for mu in 0..<fd.len:
+      a += redot(ch[mu], fl[mu]) + redot(chl[mu], ll[mu])
+      a2 += redot(ch[mu], fl2[mu]) + redot(chl[mu], ll2[mu])
+      fd[mu] := 0
+      fc[mu] := ch[mu]
+      lc[mu] := chl[mu]
+    #gc.gaugeDeriv2(fl, fc)
+    #gc.gaugeDeriv2(ll, lc)
+    #let a = gc.gaugeAction2(fl) + gc.gaugeAction2(ll)
+    #let a2 = gc.gaugeAction2(fl2) + gc.gaugeAction2(ll2)
     f(fd, fc, lc)
-    let a = gc.gaugeAction2(fl) + gc.gaugeAction2(ll)
-    let a2 = gc.gaugeAction2(fl2) + gc.gaugeAction2(ll2)
-    check(a2-a, tol)
+    f2(fd, fc, lc)  # combine forces for better accuracy
+    check(2.0*(a2-a), tol)  # 2.0 due to combined forces
 
   hisq.fat7first.oneLink = 1.0
   hisq.fat7first.threeStaple = 0.0
@@ -127,15 +141,15 @@ when isMainModule:
   hisq.fat7second.sevenStaple = 0.0
   hisq.fat7second.lepage = 0.0
   hisq.naik = 0.0
-  checkL("oneLink", 10)
+  checkL("oneLink", 1)
 
   hisq.fat7first.setHisqFat7(0.0, 0.0)
   hisq.fat7second.setHisqFat7(0.0, 0.0)
   hisq.naik = 0.0
-  checkL("all", 300)
+  checkL("all", 1)
 
   hisq = newHisq()
-  checkL("all", 700)
+  checkL("all", 1)
 
   echoProf()
   qexFinalize()
diff --git a/src/maths/complexType.nim b/src/maths/complexType.nim
index 22d943c..098d4dd 100644
--- a/src/maths/complexType.nim
+++ b/src/maths/complexType.nim
@@ -1,5 +1,5 @@
 import macros
-import base/basicOps
+import base/[basicOps, wrapperTypes]
 import complexProxy
 export complexProxy