Inconsistency when using QUDA_TWISTED_MASS_DSLASH with QUDA_MASS_NORMALIZATION

[m-schroeck]

A simple residual test for QUDA_TWISTED_MASS_DSLASH (with the default setup QUDA_MAT_SOLUTION and QUDA_NORMOP_PC_SOLVE) of invert_test.cpp fails if one uses QUDA_MASS_NORMALIZATION instead of the default QUDA_KAPPA_NORMALIZATION.

Directly after the inversion (line 393 in invert_test.cpp), calculate ||Ax-b||^2:

///////////////////////
// test ||Ax - b||^2 //
///////////////////////
int vol = inv_param.solution_type == QUDA_MAT_SOLUTION ? V : Vh;
MatQuda( spinorCheck, spinorOut, &inv_param); // Ax
mxpy(spinorIn, spinorCheck, vol_spinorSiteSize_inv_param.Ls, inv_param.cpu_prec); // Ax -= b
double Axmb_nrm2 = norm_2(spinorCheck, vol_spinorSiteSize_inv_param.Ls, inv_param.cpu_prec);
printfQuda("\n\n>>> Residual ||Ax - b||^2: %e\n\n", Axmb_nrm2);
///////////////////////

This test will succeed for the default setup. When switching to QUDA_MASS_NORMALIZATION it will fail. But when using QUDA_MASS_NORMALIZATION together with QUDA_NORMOP_SOLVE (instead of the preconditioned QUDA_NORMOP_PC_SOLVE) it will succeed again. Therefore the guess that the bug may lie in Dirac::prepare() and/or Dirac::reconstruct().

Tested in the quda-0.7 branch revision 5255 (last).

[AlexVaq]

Oh yeah, I think this has been around for a while. I remember bugging Alexei about this in my early days in Cyprus.

On 28/08/2014, at 15:42, m-schroeck notifications@github.com wrote:

A simple residual test for QUDA_TWISTED_MASS_DSLASH (with the default setup QUDA_MAT_SOLUTION and QUDA_NORMOP_PC_SOLVE) of invert_test.cpp fails if one uses QUDA_MASS_NORMALIZATION instead of the default QUDA_KAPPA_NORMALIZATION.

Directly after the inversion (line 393 in invert_test.cpp), calculate ||Ax-b||^2:

///////////////////////
// test ||Ax - b||^2 //
///////////////////////
int vol = inv_param.solution_type == QUDA_MAT_SOLUTION ? V : Vh;
MatQuda( spinorCheck, spinorOut, &inv_param); // Ax
mxpy(spinorIn, spinorCheck, vol_spinorSiteSize_inv_param.Ls, inv_param.cpu_prec); // Ax -= b
double Axmb_nrm2 = norm_2(spinorCheck, vol_spinorSiteSize_inv_param.Ls, inv_param.cpu_prec);
printfQuda("\n\n>>> Residual ||Ax - b||^2: %e\n\n", Axmb_nrm2);
///////////////////////

This test will succeed for the default setup. When switching to QUDA_MASS_NORMALIZATION it will fail. But when using QUDA_MASS_NORMALIZATION together with QUDA_NORMOP_SOLVE (instead of the preconditioned QUDA_NORMOP_PC_SOLVE) it will succeed again. Therefore the guess that the bug may lie in Dirac::prepare() and/or Dirac::reconstruct().

Tested in the quda-0.7 branch revision 5255 (last).

—
Reply to this email directly or view it on GitHub.

[maddyscientist]

I just tried to reproduce this issue, and it doesn't occur. Both mass and kappa normalization appear to work fine. Did someone fix this and not close the issue?

[AlexVaq]

I don't know, but I'm happy to hear it's working. :)

Alex

El 28/10/2014, a las 22:41, mikeaclark notifications@github.com escribió:

I just tried to reproduce this issue, and it doesn't occur. Both mass and
kappa normalization appear to work fine. Did someone fix this and not close
the issue?

—
Reply to this email directly or view it on GitHub
#151 (comment).

[alexstrel]

No, it does not work yet, e.g.:
Residuals: (L2 relative) tol 1e-07, QUDA = 6.88514e-08, host = 2.57123;

[AlexVaq]

No, it does not work yet, e.g.:
Residuals: (L2 relative) tol 1e-07, QUDA = 6.88514e-08, host = 2.57123;

Ops… Steps to reproduce? It seems Mike got it right and you got it wrong, so there must be a difference in the methodology. But of course, the bug is still there.

[m-schroeck]

I agree with Alexei, running the exact same test I described in the initial post I obtain:

Residual ||Ax - b||^2: 1.755031e+07

(hash 769e563)

[alexstrel]

ok, this is not a twisted mass bug, e.g., for wilson I got
Residuals: (L2 relative) tol 1e-07, QUDA = 7.64141e-08, host = 2.57293; (heavy-quark) tol 0.001, QUDA = 5.1883e-08

[AlexVaq]

Then it’s pretty serious. People that are not aware of this bug might have wrong results if they are not careful. I’d disable the feature until it’s fixed.

In any case, we should check this asap. It should be just a rescaling...

[m-schroeck]

what's your setup to get it wrong for wilson?

[alexstrel]

I guess same as for the twisted mass (check QUDA_MAT_SOLUTION).

[maddyscientist]

Do the problem also occur if you use SYMMETRIC preconditioning instead of ASYMMETRIC?

I suspect it's related to setting the correct preconditioning (SYMMETRIC / ASYMMETRIC) with the correct normalization (QUDA_MASS_NORMALIZATION / QUDA_ASYMMETRIC_MASS_NORMALIZATION).

---

This email message is for the sole use of the intended recipient(s) and may contain
confidential information. Any unauthorized review, use, disclosure or distribution
is prohibited. If you are not the intended recipient, please contact the sender by

reply email and destroy all copies of the original message.

[alexstrel]

for the twisted mass that was symmetric preconditioning (and wilson does not use asymmetric, anyway)

[maddyscientist]

Ok, I've reproduced the issue, but I can't see it for Wilson. Running on a 4^4 volume:

twisted-mass
kappa normalization: >>> Residual ||Ax - b||^2: 4.196415e-11
mass normalization: >>> Residual ||Ax - b||^2: 1.330323e+04
wilson
kappa normalization: >>> Residual ||Ax - b||^2: 1.041266e-11
mass normalization: >>> Residual ||Ax - b||^2: 1.246011e-11

Continuing to investigate.

[maddyscientist]

I take that back. I see I was solving a different system with Wilson. I agree then, when Wilson is set to solve MAT then it breaks with mass normalization.

[maddyscientist]

Closing this issue. The problem was caused by the fact that Dirac::prepare was being called before massRescale. Dirac::prepare gives a right hand side that is for the preconditioned system, but the mass rescale was applied thinking it was for the unpreconditioned system. Swapping the order of these two trivially fixes this issue.