portelli/Grid
1
0
Fork 0
mirror of https://github.com/paboyle/Grid.git synced 2026-05-21 09:34:17 +01:00
Code Releases Activity

found bug in unprec DWF: was using |\cdot| in comparison for the eigenvalue sorting

Browse Source
This commit is contained in:
Patrick Oare
2025-09-12 13:31:39 -04:00
parent 82d411ca7b
commit 0b92ef990c
1 changed files with 13 additions and 4 deletions
Download Patch File Download Diff File Expand all files Collapse all files
Grid/algorithms/iterative/KrylovSchur.h
+13 -4
View File
@@ -86,14 +86,22 @@ struct ComplexComparator
return std::abs(z1) < std::abs(z2); return std::abs(z1) < std::abs(z2);
case EvalNormLarge: case EvalNormLarge:
return std::abs(z1) > std::abs(z2); return std::abs(z1) > std::abs(z2);
// case EvalReSmall:
// return std::abs(std::real(z1)) < std::abs(std::real(z2)); // DELETE THE ABS HERE!!!
// case EvalReLarge:
// return std::abs(std::real(z1)) > std::abs(std::real(z2));
// case EvalImSmall:
// return std::abs(std::imag(z1)) < std::abs(std::imag(z2));
// case EvalImLarge:
// return std::abs(std::imag(z1)) > std::abs(std::imag(z2));
case EvalReSmall: case EvalReSmall:
return std::abs(std::real(z1)) < std::abs(std::real(z2)); return std::real(z1) < std::real(z2); // DELETE THE ABS HERE!!!
case EvalReLarge: case EvalReLarge:
return std::abs(std::real(z1)) > std::abs(std::real(z2)); return std::real(z1) > std::real(z2);
case EvalImSmall: case EvalImSmall:
return std::abs(std::imag(z1)) < std::abs(std::imag(z2)); return std::imag(z1) < std::imag(z2);
case EvalImLarge: case EvalImLarge:
return std::abs(std::imag(z1)) > std::abs(std::imag(z2)); return std::imag(z1) > std::imag(z2);
default: default:
assert(0); assert(0);
} }
@@ -311,6 +319,7 @@ class KrylovSchur {
ssq = norm2(v0); ssq = norm2(v0);
RealD approxLambdaMax = approxMaxEval(v0); RealD approxLambdaMax = approxMaxEval(v0);
rtol = Tolerance * approxLambdaMax; rtol = Tolerance * approxLambdaMax;
std::cout << GridLogMessage << "Approximate max eigenvalue: " << approxLambdaMax << std::endl;
// rtol = Tolerance; // rtol = Tolerance;
b = Eigen::VectorXcd::Zero(Nm); // start as e_{k+1} b = Eigen::VectorXcd::Zero(Nm); // start as e_{k+1}