Reference documentation for deal.II version 9.2.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Classes | Public Member Functions | Protected Member Functions | Private Attributes | List of all members
SolverIDR< VectorType > Class Template Reference

#include <deal.II/lac/solver_idr.h>

Inheritance diagram for SolverIDR< VectorType >:
[legend]

Classes

struct  AdditionalData
 

Public Member Functions

 SolverIDR (SolverControl &cn, VectorMemory< VectorType > &mem, const AdditionalData &data=AdditionalData())
 
 SolverIDR (SolverControl &cn, const AdditionalData &data=AdditionalData())
 
virtual ~SolverIDR () override=default
 
template<typename MatrixType , typename PreconditionerType >
void solve (const MatrixType &A, VectorType &x, const VectorType &b, const PreconditionerType &preconditioner)
 
- Public Member Functions inherited from SolverBase< Vector< double > >
 SolverBase (SolverControl &solver_control, VectorMemory< Vector< double > > &vector_memory)
 
 SolverBase (SolverControl &solver_control)
 
boost::signals2::connection connect (const std::function< SolverControl::State(const unsigned int iteration, const double check_value, const Vector< double > &current_iterate)> &slot)
 
- Public Member Functions inherited from Subscriptor
 Subscriptor ()
 
 Subscriptor (const Subscriptor &)
 
 Subscriptor (Subscriptor &&) noexcept
 
virtual ~Subscriptor ()
 
Subscriptoroperator= (const Subscriptor &)
 
Subscriptoroperator= (Subscriptor &&) noexcept
 
void subscribe (std::atomic< bool > *const validity, const std::string &identifier="") const
 
void unsubscribe (std::atomic< bool > *const validity, const std::string &identifier="") const
 
unsigned int n_subscriptions () const
 
template<typename StreamType >
void list_subscribers (StreamType &stream) const
 
void list_subscribers () const
 
template<class Archive >
void serialize (Archive &ar, const unsigned int version)
 

Protected Member Functions

virtual void print_vectors (const unsigned int step, const VectorType &x, const VectorType &r, const VectorType &d) const
 

Private Attributes

AdditionalData additional_data
 

Additional Inherited Members

- Public Types inherited from SolverBase< Vector< double > >
using vector_type = Vector< double >
 
- Static Public Member Functions inherited from Subscriptor
static ::ExceptionBaseExcInUse (int arg1, std::string arg2, std::string arg3)
 
static ::ExceptionBaseExcNoSubscriber (std::string arg1, std::string arg2)
 
- Protected Attributes inherited from SolverBase< Vector< double > >
GrowingVectorMemory< Vector< double > > static_vector_memory
 
VectorMemory< Vector< double > > & memory
 
boost::signals2::signal< SolverControl::State(const unsigned int iteration, const double check_value, const Vector< double > &current_iterate), StateCombiner > iteration_status
 

Detailed Description

template<class VectorType = Vector<double>>
class SolverIDR< VectorType >

This class implements the IDR(s) method used for solving nonsymmetric, indefinite linear systems, developed in IDR(s): A Family of Simple and Fast Algorithms for Solving Large Nonsymmetric Systems of Linear Equations by Martin B. van Gijzen and Peter Sonneveld . The implementation here is the preconditioned version from Algorithm 913: An Elegant IDR(s) Variant that Efficiently Exploits Biorthogonality Properties by Martin B. van Gijzen and Peter Sonneveld. The local structure AdditionalData takes the value for the parameter s which can be any integer greater than or equal to 1. For s=1, this method has similar convergence to BiCGStab.

Note
Each iteration of IDR(s) requires s+1 preconditioning steps and matrix-vector products. In this implementation the residual is updated and convergence is checked after each of these inner steps inside the outer iteration. If the user enables the history data, the residual at each of these steps is stored and therefore there will be multiple values per iteration.
Author
Conrad Clevenger, 2019

Definition at line 118 of file solver_idr.h.

Constructor & Destructor Documentation

◆ SolverIDR() [1/2]

template<class VectorType = Vector<double>>
SolverIDR< VectorType >::SolverIDR ( SolverControl cn,
VectorMemory< VectorType > &  mem,
const AdditionalData data = AdditionalData() 
)

Constructor.

◆ SolverIDR() [2/2]

template<class VectorType = Vector<double>>
SolverIDR< VectorType >::SolverIDR ( SolverControl cn,
const AdditionalData data = AdditionalData() 
)
explicit

Constructor. Use an object of type GrowingVectorMemory as a default to allocate memory.

◆ ~SolverIDR()

template<class VectorType = Vector<double>>
virtual SolverIDR< VectorType >::~SolverIDR ( )
overridevirtualdefault

Virtual destructor.

Member Function Documentation

◆ solve()

template<class VectorType = Vector<double>>
template<typename MatrixType , typename PreconditionerType >
void SolverIDR< VectorType >::solve ( const MatrixType &  A,
VectorType x,
const VectorType b,
const PreconditionerType &  preconditioner 
)

Solve the linear system Ax=b for x.

◆ print_vectors()

template<class VectorType = Vector<double>>
virtual void SolverIDR< VectorType >::print_vectors ( const unsigned int  step,
const VectorType x,
const VectorType r,
const VectorType d 
) const
protectedvirtual

Interface for derived class. This function gets the current iteration vector, the residual and the update vector in each step. It can be used for graphical output of the convergence history.

Member Data Documentation

◆ additional_data

template<class VectorType = Vector<double>>
AdditionalData SolverIDR< VectorType >::additional_data
private

Additional solver parameters.

Definition at line 181 of file solver_idr.h.


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