portelli/Grid
1
0
Fork 0
mirror of https://github.com/paboyle/Grid.git synced 2025-12-16 10:44:40 +00:00
Code Releases Activity

Make IterationsToComplete and TrueResidual externally accessible

Browse Source
This commit is contained in:
Carleton DeTar
2020-02-19 17:43:57 -06:00
parent 5d834486c9
commit 9479bc8486
3 changed files with 45 additions and 7 deletions
Download Patch File Download Diff File Expand all files Collapse all files

19
Grid/algorithms/iterative/ConjugateGradient.h
View File

@@ -49,6 +49,7 @@ public:
RealD Tolerance;
Integer MaxIterations;
Integer IterationsToComplete; //Number of iterations the CG took to finish. Filled in upon completion
RealD TrueResidual;
ConjugateGradient(RealD tol, Integer maxit, bool err_on_no_conv = true)
: Tolerance(tol),
@@ -82,6 +83,14 @@ public:
cp = a;
ssq = norm2(src);
// Handle trivial case of zero src
if (ssq == 0.){
psi = 0;
IterationsToComplete = 1;
TrueResidual = 0.;
return;
}
std::cout << GridLogIterative << std::setprecision(8) << "ConjugateGradient: guess " << guess << std::endl;
std::cout << GridLogIterative << std::setprecision(8) << "ConjugateGradient: src " << ssq << std::endl;
std::cout << GridLogIterative << std::setprecision(8) << "ConjugateGradient: mp " << d << std::endl;
@@ -93,6 +102,7 @@ public:
// Check if guess is really REALLY good :)
if (cp <= rsq) {
TrueResidual = sqrt(a/ssq);
std::cout << GridLogMessage << "ConjugateGradient guess is converged already " << std::endl;
IterationsToComplete = 0;
return;
@@ -142,7 +152,7 @@ public:
LinalgTimer.Stop();
std::cout << GridLogIterative << "ConjugateGradient: Iteration " << k
<< " residual^2 " << sqrt(cp/ssq) << " target " << Tolerance << std::endl;
<< " residual " << sqrt(cp/ssq) << " target " << Tolerance << std::endl;
// Stopping condition
if (cp <= rsq) {
@@ -170,10 +180,17 @@ public:
if (ErrorOnNoConverge) assert(true_residual / Tolerance < 10000.0);
IterationsToComplete = k;
TrueResidual = true_residual;
return;
}
}
// Failed. Calculate true residual before giving up
Linop.HermOpAndNorm(psi, mmp, d, qq);
p = mmp - src;
TrueResidual = sqrt(norm2(p)/ssq);
std::cout << GridLogMessage << "ConjugateGradient did NOT converge "<<k<<" / "<< MaxIterations<< std::endl;
if (ErrorOnNoConverge) assert(0);