diff --git a/lib/qcd/hmc/integrators/Integrator.h b/lib/qcd/hmc/integrators/Integrator.h index bba8b513..b162bcd2 100644 --- a/lib/qcd/hmc/integrators/Integrator.h +++ b/lib/qcd/hmc/integrators/Integrator.h @@ -125,10 +125,10 @@ class Integrator { Mom -= MomDer * ep; // Auxiliary fields - //P.update_auxiliary_momenta(ep*0.5); - //P.AuxiliaryFieldsDerivative(MomDer); - //Mom -= MomDer * ep; - //P.update_auxiliary_momenta(ep*0.5); + P.update_auxiliary_momenta(ep*0.5); + P.AuxiliaryFieldsDerivative(MomDer); + Mom -= MomDer * ep; + P.update_auxiliary_momenta(ep*0.5); for (int a = 0; a < as[level].actions.size(); ++a) { @@ -193,9 +193,9 @@ class Integrator { } // Auxiliary fields - //P.update_auxiliary_momenta(ep*0.5); - //P.AuxiliaryFieldsDerivative(AuxDer); - //Msum += AuxDer; + P.update_auxiliary_momenta(ep*0.5); + P.AuxiliaryFieldsDerivative(AuxDer); + Msum += AuxDer; // Here run recursively @@ -222,7 +222,7 @@ class Integrator { P.Mom = NewMom; // update the auxiliary fields momenta - //P.update_auxiliary_momenta(ep*0.5); + P.update_auxiliary_momenta(ep*0.5); } @@ -265,7 +265,7 @@ class Integrator { P.DerivativeP(Mom1); // first term in the derivative - //P.update_auxiliary_fields(ep*0.5); + P.update_auxiliary_fields(ep*0.5); do { @@ -293,7 +293,7 @@ class Integrator { U = NewU; - //P.update_auxiliary_fields(ep*0.5); + P.update_auxiliary_fields(ep*0.5); } @@ -340,6 +340,7 @@ class Integrator { void reverse_momenta(){ P.Mom *= -1.0; + P.AuxMom *= -1.0; } // to be used by the actionlevel class to iterate diff --git a/lib/qcd/hmc/integrators/Integrator_algorithm.h b/lib/qcd/hmc/integrators/Integrator_algorithm.h index 3e41a266..d61e28ee 100644 --- a/lib/qcd/hmc/integrators/Integrator_algorithm.h +++ b/lib/qcd/hmc/integrators/Integrator_algorithm.h @@ -344,7 +344,8 @@ class ImplicitLeapFrog : public Integratorupdate_P(U, level, eps / 2.0); } else { - this->implicit_update_P(U, level, eps, true);// intermediate step + this->implicit_update_P(U, level, eps, true);// works intermediate step + // this->update_P(U, level, eps); // looks not reversible } } } diff --git a/lib/qcd/utils/Metric.h b/lib/qcd/utils/Metric.h index 8eeb5846..60a9bfc5 100644 --- a/lib/qcd/utils/Metric.h +++ b/lib/qcd/utils/Metric.h @@ -105,7 +105,7 @@ public: // Modify the distribution with the metric M.MSquareRoot(Mom); - if (0) { + if (1) { // Auxiliary momenta // do nothing if trivial, so hide in the metric MomentaField AuxMomTemp(Mom._grid); @@ -132,7 +132,7 @@ public: Hloc += trace(Mom_mu * inv_mu); } - if (0) { + if (1) { // Auxiliary Fields // hide in the metric M.M(AuxMom, inv); @@ -167,7 +167,7 @@ public: void AuxiliaryFieldsDerivative(MomentaField& der){ der = zero; - if (0){ + if (1){ // Auxiliary fields MomentaField der_temp(der._grid); MomentaField X(der._grid); @@ -194,13 +194,13 @@ public: } void update_auxiliary_momenta(RealD ep){ - if(0){ + if(1){ AuxMom -= ep * AuxField; } } void update_auxiliary_fields(RealD ep){ - if (0) { + if (1) { MomentaField tmp(AuxMom._grid); MomentaField tmp2(AuxMom._grid); M.M(AuxMom, tmp);