Block-Structured AMR Software Framework
Loading...
Searching...
No Matches
amrex::GMRES_MV< T > Class Template Reference

#include <AMReX_GMRES_MV.H>

Public Types

using RT = Real
 
using VEC = AlgVector< T >
 
using MAT = SpMatrix< T >
 
using GM = amrex::GMRES< VEC, GMRES_MV< T > >
 
using PC = std::function< void(VEC &, VEC const &)>
 

Public Member Functions

 GMRES_MV (MAT const *a_mat)
 Bind GMRES to a sparse matrix described by a_mat.
 
void setPrecond (PC a_pc)
 Supply an optional right-preconditioner functor.
 
void solve (VEC &a_sol, VEC const &a_rhs, T a_tol_rel, T a_tol_abs)
 Solve the linear system.
 
void setVerbose (int v)
 Set verbosity level v for the underlying GMRES solver.
 
GMgetGMRES ()
 Access the underlying GMRES object for additional tuning.
 
VEC makeVecRHS () const
 Return an AlgVector that mirrors the matrix partition for RHS storage.
 
VEC makeVecLHS () const
 Return another AlgVector with the same partition for LHS storage.
 
void apply (VEC &lhs, VEC &rhs) const
 Apply the sparse operator to rhs and store the result in lhs.
 
void precond (VEC &lhs, VEC const &rhs) const
 Apply the optional preconditioner (or copy if none is provided).
 

Static Public Member Functions

static T norm2 (VEC const &vec)
 Euclidean norm of vec.
 
static void scale (VEC &vec, T scale_factor)
 Scale vec in place by scale_factor.
 
static T dotProduct (VEC const &vec1, VEC const &vec2)
 Dot product between vec1 and vec2.
 
static void setToZero (VEC &lhs)
 Reset lhs to zero.
 
static void assign (VEC &lhs, VEC const &rhs)
 Copy rhs into lhs.
 
static void increment (VEC &lhs, VEC const &rhs, T a)
 Accumulate $\text{lhs} \mathrel{+}= a\,\text{rhs}$.
 
static void linComb (VEC &lhs, T a, VEC const &rhs_a, T b, VEC const &rhs_b)
 Form the linear combination $\text{lhs} = a\,\text{rhs}_a + b\,\text{rhs}_b$.
 

Member Typedef Documentation

◆ GM

template<typename T >
using amrex::GMRES_MV< T >::GM = amrex::GMRES<VEC,GMRES_MV<T> >

◆ MAT

template<typename T >
using amrex::GMRES_MV< T >::MAT = SpMatrix<T>

◆ PC

template<typename T >
using amrex::GMRES_MV< T >::PC = std::function<void(VEC&,VEC const&)>

◆ RT

template<typename T >
using amrex::GMRES_MV< T >::RT = Real

◆ VEC

template<typename T >
using amrex::GMRES_MV< T >::VEC = AlgVector<T>

Constructor & Destructor Documentation

◆ GMRES_MV()

template<typename T >
amrex::GMRES_MV< T >::GMRES_MV ( MAT const *  a_mat)

Bind GMRES to a sparse matrix described by a_mat.

Parameters
a_matMatrix supplying the SpMV operation and partition.

Member Function Documentation

◆ apply()

template<typename T >
void amrex::GMRES_MV< T >::apply ( VEC lhs,
VEC rhs 
) const

Apply the sparse operator to rhs and store the result in lhs.

Parameters
lhsDestination vector updated in place.
rhsSource vector produced by the preconditioner (passed by non-const reference to satisfy the GMRES interface).

◆ assign()

template<typename T >
void amrex::GMRES_MV< T >::assign ( VEC lhs,
VEC const &  rhs 
)
static

Copy rhs into lhs.

Parameters
lhsDestination vector updated in place.
rhsSource vector.

◆ dotProduct()

template<typename T >
T amrex::GMRES_MV< T >::dotProduct ( VEC const &  vec1,
VEC const &  vec2 
)
static

Dot product between vec1 and vec2.

◆ getGMRES()

template<typename T >
GM & amrex::GMRES_MV< T >::getGMRES ( )
inline

Access the underlying GMRES object for additional tuning.

◆ increment()

template<typename T >
void amrex::GMRES_MV< T >::increment ( VEC lhs,
VEC const &  rhs,
a 
)
static

Accumulate $\text{lhs} \mathrel{+}= a\,\text{rhs}$.

Parameters
lhsAccumulation target.
rhsSource vector.
aScalar factor applied to rhs.

◆ linComb()

template<typename T >
void amrex::GMRES_MV< T >::linComb ( VEC lhs,
a,
VEC const &  rhs_a,
b,
VEC const &  rhs_b 
)
static

Form the linear combination $\text{lhs} = a\,\text{rhs}_a + b\,\text{rhs}_b$.

Parameters
lhsDestination vector.
aScalar applied to rhs_a.
rhs_aFirst source vector.
bScalar applied to rhs_b.
rhs_bSecond source vector.

◆ makeVecLHS()

template<typename T >
auto amrex::GMRES_MV< T >::makeVecLHS ( ) const

Return another AlgVector with the same partition for LHS storage.

◆ makeVecRHS()

template<typename T >
auto amrex::GMRES_MV< T >::makeVecRHS ( ) const

Return an AlgVector that mirrors the matrix partition for RHS storage.

◆ norm2()

template<typename T >
T amrex::GMRES_MV< T >::norm2 ( VEC const &  vec)
static

Euclidean norm of vec.

◆ precond()

template<typename T >
void amrex::GMRES_MV< T >::precond ( VEC lhs,
VEC const &  rhs 
) const

Apply the optional preconditioner (or copy if none is provided).

◆ scale()

template<typename T >
void amrex::GMRES_MV< T >::scale ( VEC vec,
scale_factor 
)
static

Scale vec in place by scale_factor.

◆ setPrecond()

template<typename T >
void amrex::GMRES_MV< T >::setPrecond ( PC  a_pc)
inline

Supply an optional right-preconditioner functor.

The functor must follow the signature void(VEC& lhs, VEC const& rhs).

Parameters
a_pcCallable that computes lhs = P^{-1}(rhs).

◆ setToZero()

template<typename T >
void amrex::GMRES_MV< T >::setToZero ( VEC lhs)
static

Reset lhs to zero.

Parameters
lhsVector whose entries are zeroed asynchronously.

◆ setVerbose()

template<typename T >
void amrex::GMRES_MV< T >::setVerbose ( int  v)
inline

Set verbosity level v for the underlying GMRES solver.

◆ solve()

template<typename T >
void amrex::GMRES_MV< T >::solve ( VEC a_sol,
VEC const &  a_rhs,
a_tol_rel,
a_tol_abs 
)

Solve the linear system.

Parameters
a_solunknowns, i.e., x in A x = b.
a_rhsRHS, i.e., b in A x = b.
a_tol_relrelative tolerance.
a_tol_absabsolute tolerance.

The documentation for this class was generated from the following file: