GMRESSolver.h

Go to the documentation of this file.

#ifndef included_AMP_GMRESSolver
#define included_AMP_GMRESSolver
 
#include "AMP/solvers/SolverStrategy.h"
#include "AMP/solvers/SolverStrategyParameters.h"
#include "AMP/utils/AMP_MPI.h"
#include "AMP/utils/Array.h"
 
namespace AMP::Solver {
 
template<typename T = double>
class GMRESSolver : public SolverStrategy
{
public:
    GMRESSolver();
 
    explicit GMRESSolver( std::shared_ptr<SolverStrategyParameters> params );
 
    static std::unique_ptr<SolverStrategy>
    createSolver( std::shared_ptr<SolverStrategyParameters> params )
    {
        return std::make_unique<GMRESSolver<T>>( params );
    }
 
    virtual ~GMRESSolver() = default;
 
    std::string type() const override { return "GMRESSolver"; }
 
    void apply( AMP::LinearAlgebra::Vector::const_shared_ptr f,
                AMP::LinearAlgebra::Vector::shared_ptr u ) override;
 
    void initialize( std::shared_ptr<const SolverStrategyParameters> params ) override;
 
    void registerOperator( std::shared_ptr<AMP::Operator::Operator> op ) override;
 
    bool restarted() { return d_bRestart; }
 
protected:
    // stored internally. Store the coefficients of the Arnoldi
    // iteration internally in a upper Hessenberg matrix
    virtual void orthogonalize( const int k, std::shared_ptr<AMP::LinearAlgebra::Vector> v );
 
    std::vector<T> basisInnerProducts( const int k, std::shared_ptr<AMP::LinearAlgebra::Vector> v );
 
    // stored internally using classical Gram-Schmidt. Store the coefficients of the Arnoldi
    // iteration internally in a upper Hessenberg matrix
    void cgs( const int k, std::shared_ptr<AMP::LinearAlgebra::Vector> v );
 
    // stored internally using classical Gram-Schmidt with re-orthogonalization. Store the
    // coefficients of the Arnoldi iteration internally in a upper Hessenberg matrix
    void cgs2( const int k, std::shared_ptr<AMP::LinearAlgebra::Vector> v );
 
    void applyGivensRotation( const int i, const int k );
 
    void computeGivensRotation( const int k );
 
    void backwardSolve( const int nr );
 
    void addCorrection( const int nr,
                        std::shared_ptr<AMP::LinearAlgebra::Vector> z,
                        std::shared_ptr<AMP::LinearAlgebra::Vector> z1,
                        std::shared_ptr<AMP::LinearAlgebra::Vector> u );
 
    void computeInitialResidual( bool use_zero_guess,
                                 std::shared_ptr<const AMP::LinearAlgebra::Vector> f,
                                 std::shared_ptr<AMP::LinearAlgebra::Vector> u,
                                 std::shared_ptr<AMP::LinearAlgebra::Vector> tmp,
                                 std::shared_ptr<AMP::LinearAlgebra::Vector> res );
 
    void getFromInput( std::shared_ptr<AMP::Database> db );
 
private:
    void allocateBasis( std::shared_ptr<const AMP::LinearAlgebra::Vector> u = nullptr );
 
    bool d_bRestart = false; 
 
    int d_iMaxKrylovDimension = 50; 
 
    int d_iBasisAllocSize = 4; 
 
    int d_restarts; 
 
    std::string d_sOrthogonalizationMethod = "MGS";
 
    std::string d_preconditioner_side = "right";
 
    bool d_bUsesPreconditioner = false;
 
    bool d_bFlexibleGMRES = false;
 
    AMP::Array<T> d_dHessenberg;
 
    // when Givens is being used
    std::vector<T> d_dcos;
    std::vector<T> d_dsin;
 
    std::vector<T> d_dw;
 
    std::vector<T> d_dy;
 
    std::vector<std::shared_ptr<AMP::LinearAlgebra::Vector>> d_vBasis;
 
    std::vector<std::shared_ptr<AMP::LinearAlgebra::Vector>> d_zBasis;
 
    std::shared_ptr<AMP::LinearAlgebra::Vector> d_z;
    std::shared_ptr<AMP::LinearAlgebra::Vector> d_z1;
};
} // namespace AMP::Solver
 
#endif