Hypre::HyprePreconditioner< M, X, Y, Comm > Class Template Reference

Wrapper for Hypre's BoomerAMG preconditioner. More...

#include <HyprePreconditioner.hpp>

Inheritance diagram for Hypre::HyprePreconditioner< M, X, Y, Comm >:
Inheritance graph

Public Types

using matrix_type = M
 The matrix type the preconditioner is for. More...
 
using matrix_field_type = typename M::field_type
 The field type of the matrix. More...
 
using domain_type = X
 The domain type of the preconditioner. More...
 
using range_type = Y
 The range type of the preconditioner. More...
 
using vector_field_type = typename X::field_type
 The field type of the vectors. More...
 

Public Member Functions

 HyprePreconditioner (const M &A, const Opm::PropertyTree prm, const Comm &comm)
 Constructor for the HyprePreconditioner class. More...
 
 ~HyprePreconditioner ()
 Destructor for HyprePreconditioner. More...
 
void update () override
 Updates the preconditioner with the current matrix values. More...
 
void pre (X &v, Y &) override
 Pre-processing step before applying the preconditioner. More...
 
void apply (X &v, const Y &d) override
 Applies the preconditioner to a vector. More...
 
void post (X &) override
 Post-processing step after applying the preconditioner. More...
 
Dune::SolverCategory::Category category () const override
 Returns the solver category. More...
 
bool hasPerfectUpdate () const override
 Checks if the preconditioner has a perfect update. More...
 

Detailed Description

template<class M, class X, class Y, class Comm>
class Hypre::HyprePreconditioner< M, X, Y, Comm >

Wrapper for Hypre's BoomerAMG preconditioner.

This class provides an interface to the BoomerAMG preconditioner from the Hypre library. It is designed to work with matrices, update vectors, and defect vectors specified by the template parameters. The HypreInterface class provides a unified interface to Hypre's functionality, allowing for easy switching between CPU and GPU input data types and backend acceleration.

Supports four use cases:

  1. Input type is CPU and backend acceleration is CPU
  2. Input type is CPU and backend acceleration is GPU
  3. Input type is GPU and backend acceleration is GPU
  4. Input type is GPU and backend acceleration is CPU
Template Parameters
MThe matrix type
XThe vector type for the solution
YThe vector type for the right-hand side

Member Typedef Documentation

◆ domain_type

template<class M , class X , class Y , class Comm >
using Hypre::HyprePreconditioner< M, X, Y, Comm >::domain_type = X

The domain type of the preconditioner.

◆ matrix_field_type

template<class M , class X , class Y , class Comm >
using Hypre::HyprePreconditioner< M, X, Y, Comm >::matrix_field_type = typename M::field_type

The field type of the matrix.

◆ matrix_type

template<class M , class X , class Y , class Comm >
using Hypre::HyprePreconditioner< M, X, Y, Comm >::matrix_type = M

The matrix type the preconditioner is for.

◆ range_type

template<class M , class X , class Y , class Comm >
using Hypre::HyprePreconditioner< M, X, Y, Comm >::range_type = Y

The range type of the preconditioner.

◆ vector_field_type

template<class M , class X , class Y , class Comm >
using Hypre::HyprePreconditioner< M, X, Y, Comm >::vector_field_type = typename X::field_type

The field type of the vectors.

Constructor & Destructor Documentation

◆ HyprePreconditioner()

template<class M , class X , class Y , class Comm >
Hypre::HyprePreconditioner< M, X, Y, Comm >::HyprePreconditioner ( const M &  A,
const Opm::PropertyTree  prm,
const Comm comm 
)
inline

Constructor for the HyprePreconditioner class.

Initializes the preconditioner with the given matrix and property tree.

Parameters
AThe matrix for which the preconditioner is constructed.
prmThe property tree containing configuration parameters.
commParallel communicator.

References Opm::gpuistl::SparsityPattern::cols, Opm::gpuistl::HypreDeviceDataArrays::cols_device, Opm::gpuistl::HypreInterface::computeRowIndexes(), Opm::gpuistl::HypreHostDataArrays::continuous_vector_values, Opm::gpuistl::HypreInterface::createAMGSolver(), Opm::gpuistl::HypreInterface::createMatrix(), Opm::gpuistl::HypreInterface::createVector(), Opm::gpuistl::ParallelInfo::dof_offset, Opm::gpuistl::HypreHostDataArrays::indices, Opm::gpuistl::HypreDeviceDataArrays::indices_device, Opm::gpuistl::HypreInterface::initialize(), Opm::gpuistl::ParallelInfo::local_dune_to_local_hypre, Opm::gpuistl::HypreDeviceDataArrays::matrix_buffer_device, Opm::gpuistl::ParallelInfo::N_owned, Opm::gpuistl::SparsityPattern::ncols, Opm::gpuistl::HypreDeviceDataArrays::ncols_device, Opm::gpuistl::SparsityPattern::nnz, Opm::gpuistl::ParallelInfo::owner_first, Opm::gpuistl::HypreHostDataArrays::row_indexes, Opm::gpuistl::HypreDeviceDataArrays::row_indexes_device, Opm::gpuistl::SparsityPattern::rows, Opm::gpuistl::HypreDeviceDataArrays::rows_device, Opm::gpuistl::HypreInterface::setSolverParameters(), Opm::gpuistl::HypreInterface::setupHypreParallelInfo(), Opm::gpuistl::HypreInterface::setupSparsityPattern(), Hypre::HyprePreconditioner< M, X, Y, Comm >::update(), and Opm::gpuistl::HypreDeviceDataArrays::vector_buffer_device.

◆ ~HyprePreconditioner()

Member Function Documentation

◆ apply()

template<class M , class X , class Y , class Comm >
void Hypre::HyprePreconditioner< M, X, Y, Comm >::apply ( X &  v,
const Y &  d 
)
inlineoverride

Applies the preconditioner to a vector.

Performs one AMG V-cycle to solve the system. Involves uploading vectors to Hypre, applying the preconditioner, and transferring the result back to the vector.

Parameters
vThe update vector.
dThe defect vector.

References HYPRE_SAFE_CALL, Opm::gpuistl::HypreInterface::transferVectorFromHypre(), and Opm::gpuistl::HypreInterface::transferVectorToHypre().

◆ category()

template<class M , class X , class Y , class Comm >
Dune::SolverCategory::Category Hypre::HyprePreconditioner< M, X, Y, Comm >::category ( ) const
inlineoverride

Returns the solver category.

Returns
The solver category, which is sequential.

◆ hasPerfectUpdate()

template<class M , class X , class Y , class Comm >
bool Hypre::HyprePreconditioner< M, X, Y, Comm >::hasPerfectUpdate ( ) const
inlineoverridevirtual

Checks if the preconditioner has a perfect update.

Returns
True, indicating that the preconditioner can be perfectly updated.

Implements Dune::PreconditionerWithUpdate< X, Y >.

◆ post()

template<class M , class X , class Y , class Comm >
void Hypre::HyprePreconditioner< M, X, Y, Comm >::post ( X &  )
inlineoverride

Post-processing step after applying the preconditioner.

This method is currently a no-op.

Parameters
vThe update vector.

◆ pre()

template<class M , class X , class Y , class Comm >
void Hypre::HyprePreconditioner< M, X, Y, Comm >::pre ( X &  v,
Y &   
)
inlineoverride

Pre-processing step before applying the preconditioner.

This method is currently a no-op.

Parameters
vThe update vector.
dThe defect vector.

◆ update()

template<class M , class X , class Y , class Comm >
void Hypre::HyprePreconditioner< M, X, Y, Comm >::update ( )
inlineoverridevirtual

Updates the preconditioner with the current matrix values.

This method should be called whenever the matrix values change.

Implements Dune::PreconditionerWithUpdate< X, Y >.

References HYPRE_SAFE_CALL, and Opm::gpuistl::HypreInterface::updateMatrixValues().

Referenced by Hypre::HyprePreconditioner< M, X, Y, Comm >::HyprePreconditioner().


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