portelli/Grid
1
0
Fork 0
mirror of https://github.com/paboyle/Grid.git synced 2026-03-19 10:46:10 +00:00
Code Releases Activity

updated RitzFilter enum and the input to run krylov schur

Browse Source
This commit is contained in:
Patrick Oare
2025-09-09 13:02:11 -04:00
parent 9dcd7ca761
commit c5d02e5799
4 changed files with 79 additions and 43 deletions
Download Patch File Download Diff File Expand all files Collapse all files

2
Grid/algorithms/Algorithms.h
View File

@@ -81,7 +81,7 @@ NAMESPACE_CHECK(PowerMethod);
NAMESPACE_CHECK(multigrid);
#include <Grid/algorithms/FFT.h>
#include <Grid/algorithms/iterative/Arnoldi.h>
#include <Grid/algorithms/iterative/KrylovSchur.h>
#include <Grid/algorithms/iterative/Arnoldi.h>
#endif

33
Grid/algorithms/iterative/Arnoldi.h
View File

@@ -8,8 +8,6 @@ Copyright (C) 2015
Author: Peter Boyle <paboyle@ph.ed.ac.uk>
Author: paboyle <paboyle@ph.ed.ac.uk>
Author: Chulwoo Jung <chulwoo@bnl.gov>
Author: Christoph Lehner <clehner@bnl.gov>
Author: Patrick Oare <poare@bnl.gov>
This program is free software; you can redistribute it and/or modify
@@ -34,37 +32,6 @@ See the full license in the file "LICENSE" in the top level distribution directo
NAMESPACE_BEGIN(Grid);
/**
* Options for which Ritz values to keep in implicit restart.
*/
enum RitzFilter {
EvalNormSmall, // Keep evals with smallest norm
EvalNormLarge, // Keep evals with largest norm
EvalReSmall, // Keep evals with smallest real part
EvalReLarge // Keep evals with largest real part
};
// Select comparison function from RitzFilter
struct ComplexComparator
{
RitzFilter f;
ComplexComparator (RitzFilter _f) : f(_f) {}
bool operator()(std::complex<double> z1, std::complex<double> z2) {
switch (f) {
case EvalNormSmall:
return std::abs(z1) < std::abs(z2);
case EvalNormLarge:
return std::abs(z1) > std::abs(z2);
case EvalReSmall:
return std::abs(std::real(z1)) < std::abs(std::real(z2));
case EvalReLarge:
return std::abs(std::real(z1)) > std::abs(std::real(z2));
default:
assert(0);
}
}
};
/**
* Implementation of the Arnoldi algorithm.
*/

76
Grid/algorithms/iterative/KrylovSchur.h
View File

@@ -8,8 +8,6 @@ Copyright (C) 2015
Author: Peter Boyle <paboyle@ph.ed.ac.uk>
Author: paboyle <paboyle@ph.ed.ac.uk>
Author: Chulwoo Jung <chulwoo@bnl.gov>
Author: Christoph Lehner <clehner@bnl.gov>
Author: Patrick Oare <poare@bnl.gov>
This program is free software; you can redistribute it and/or modify
@@ -32,14 +30,76 @@ See the full license in the file "LICENSE" in the top level distribution directo
#ifndef GRID_KRYLOVSCHUR_H
#define GRID_KRYLOVSCHUR_H
// #include <Grid/Grid.h>
// #include <Grid/parallelIO/IldgIOtypes.h>
// #include <Grid/parallelIO/IldgIO.h>
// #include <Grid/serialisation/emptyUserRecord.h>
NAMESPACE_BEGIN(Grid);
/**
* Options for which Ritz values to keep in implicit restart. TODO move this and utilities into a new file
*/
enum RitzFilter {
EvalNormSmall, // Keep evals with smallest norm
EvalNormLarge, // Keep evals with largest norm
EvalReSmall, // Keep evals with smallest real part
EvalReLarge, // Keep evals with largest real part
EvalImSmall, // Keep evals with smallest imaginary part
EvalImLarge // Keep evals with largest imaginary part
};
/** Selects the RitzFilter corresponding to the input string. */
RitzFilter selectRitzFilter(std::string s) {
if (s == "EvalNormSmall") { return EvalNormSmall; } else
if (s == "EvalNormLarge") { return EvalNormLarge; } else
if (s == "EvalReSmall") { return EvalReSmall; } else
if (s == "EvalReLarge") { return EvalReLarge; } else
if (s == "EvalImSmall") { return EvalImSmall; } else
if (s == "EvalImLarge") { return EvalImLarge; } else
{ assert(0); }
}
/** Returns a string saying which RitzFilter it is. */
std::string rfToString(RitzFilter RF) {
switch (RF) {
case EvalNormSmall:
return "EvalNormSmall";
case EvalNormLarge:
return "EvalNormLarge";
case EvalReSmall:
return "EvalReSmall";
case EvalReLarge:
return "EvalReLarge";
case EvalImSmall:
return "EvalImSmall";
case EvalImLarge:
return "EvalImLarge";
default:
assert(0);
}
}
// Select comparison function from RitzFilter
struct ComplexComparator
{
RitzFilter RF;
ComplexComparator (RitzFilter _rf) : RF(_rf) {}
bool operator()(std::complex<double> z1, std::complex<double> z2) {
switch (RF) {
case EvalNormSmall:
return std::abs(z1) < std::abs(z2);
case EvalNormLarge:
return std::abs(z1) > std::abs(z2);
case EvalReSmall:
return std::abs(std::real(z1)) < std::abs(std::real(z2));
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));
default:
assert(0);
}
}
};
/**
* Computes a complex Schur decomposition of a complex matrix A using Eigen's matrix library. The Schur decomposition,
* A = Q^\dag S Q

11
examples/Example_spec_kryschur.cc
View File

@@ -350,6 +350,15 @@ int main (int argc, char ** argv)
std::string file = argv[5];
std::string outDir = argv[6];
RitzFilter RF;
if (argc == 8) {
std::string rf = argv[7];
RF = selectRitzFilter(rf);
} else {
RF = EvalReSmall;
}
std::cout << "Sorting eigenvalues using " << rfToString(RF) << std::endl;
//const int Ls=16;
const int Ls = 8;
@@ -438,7 +447,7 @@ int main (int argc, char ** argv)
<< ", Nstop = " << Nstop << std::endl;
// Arnoldi Arn(PVdagM, FGrid, 1e-8);
// Arn(src, maxIter, Nm, Nk, Nstop);
KrylovSchur KrySchur (PVdagM, FGrid, 1e-8);
KrylovSchur KrySchur (PVdagM, FGrid, 1e-8, RF);
KrySchur(src, maxIter, Nm, Nk, Nstop);
std::cout<<GridLogMessage << "*******************************************" << std::endl;