portelli/Grid
1
0
Fork 0
mirror of https://github.com/paboyle/Grid.git synced 2025-06-19 08:17:05 +01:00
Code Releases Activity

Fixed an obscure but reproducible hang in the RHMC caused by the bounds check being activated by a random number that wasn't synchronized over the nodes

Browse Source 
HMC now also reports the "L-infinity norm" of the impulse, aka the largest site norm
This commit is contained in:
Christopher Kelly
2021-02-09 12:55:46 -05:00
parent c2676853ca
commit daa095c519
3 changed files with 22 additions and 5 deletions
Download Patch File Download Diff File Expand all files Collapse all files

8
Grid/qcd/action/pseudofermion/GeneralEvenOddRationalRatio.h
View File

@ -232,13 +232,19 @@ NAMESPACE_BEGIN(Grid);
multiShiftInverse(Denominator, ApproxNegHalfPowerAction, param.MaxIter, X,Y);
// Randomly apply rational bounds checks.
if ( param.BoundsCheckFreq != 0 && (rand()%param.BoundsCheckFreq)==0 ) {
int rcheck = rand();
CartesianCommunicator::BroadcastWorld(0,(void *)&rcheck,sizeof(int)); //make sure all nodes have the same number or you will sporadically hang and spend days trying to find out why (trust me - CK)
if ( param.BoundsCheckFreq != 0 && (rcheck % param.BoundsCheckFreq)==0 ) {
std::cout<<GridLogMessage << action_name() << " compute action: doing bounds check" << std::endl;
FermionField gauss(NumOp.FermionRedBlackGrid());
gauss = PhiOdd;
SchurDifferentiableOperator<Impl> MdagM(DenOp);
std::cout<<GridLogMessage << action_name() << " compute action: checking high bounds" << std::endl;
HighBoundCheck(MdagM,gauss,param.hi);
std::cout<<GridLogMessage << action_name() << " compute action: full approximation" << std::endl;
InversePowerBoundsCheck(param.inv_pow,param.MaxIter,param.action_tolerance*100,MdagM,gauss,ApproxNegPowerAction);
std::cout<<GridLogMessage << action_name() << " compute action: bounds check complete" << std::endl;
}
// Phidag VdagV^1/(2*inv_pow) MdagM^-1/(2*inv_pow) MdagM^-1/(2*inv_pow) VdagV^1/(2*inv_pow) Phi

10
Grid/qcd/hmc/integrators/Integrator.h
View File

@ -124,8 +124,14 @@ protected:
if (as[level].actions.at(a)->is_smeared) Smearer.smeared_force(force);
force = FieldImplementation::projectForce(force); // Ta for gauge fields
double end_force = usecond();
Real force_abs = std::sqrt(norm2(force)/U.Grid()->gSites());
std::cout << GridLogIntegrator << "["<<level<<"]["<<a<<"] Force average: " << force_abs << " Time step: " << ep << " Impulse average: " << force_abs * ep * HMC_MOMENTUM_DENOMINATOR << std::endl;
Real force_abs = std::sqrt(norm2(force)/U.Grid()->gSites()); //average per-site norm. nb. norm2(latt) = \sum_x norm2(latt[x])
Real impulse_abs = force_abs * ep * HMC_MOMENTUM_DENOMINATOR;
Real max_force_abs = std::sqrt(maxLocalNorm2(force));
Real max_impulse_abs = max_force_abs * ep * HMC_MOMENTUM_DENOMINATOR;
std::cout << GridLogIntegrator << "["<<level<<"]["<<a<<"] Force average: " << force_abs << " Max force: " << max_force_abs << " Time step: " << ep << " Impulse average: " << impulse_abs << " Max impulse: " << max_impulse_abs << std::endl;
Mom -= force * ep* HMC_MOMENTUM_DENOMINATOR;;
double end_full = usecond();
double time_full = (end_full - start_full) / 1e3;