Reference documentation for deal.II version 8.4.1
Public Member Functions | Public Attributes | List of all members
SolverGMRES< VectorType >::AdditionalData Struct Reference

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

Public Member Functions

 AdditionalData (const unsigned int max_n_tmp_vectors=30, const bool right_preconditioning=false, const bool use_default_residual=true, const bool force_re_orthogonalization=false)
 
 AdditionalData (const unsigned int max_n_tmp_vectors, const bool right_preconditioning, const bool use_default_residual, const bool force_re_orthogonalization, const bool compute_eigenvalues) DEAL_II_DEPRECATED
 

Public Attributes

unsigned int max_n_tmp_vectors
 
bool right_preconditioning
 
bool use_default_residual
 
bool force_re_orthogonalization
 
bool compute_eigenvalues
 

Detailed Description

template<class VectorType = Vector<double>>
struct SolverGMRES< VectorType >::AdditionalData

Standardized data struct to pipe additional data to the solver.

Definition at line 180 of file solver_gmres.h.

Constructor & Destructor Documentation

template<class VectorType = Vector<double>>
SolverGMRES< VectorType >::AdditionalData::AdditionalData ( const unsigned int  max_n_tmp_vectors = 30,
const bool  right_preconditioning = false,
const bool  use_default_residual = true,
const bool  force_re_orthogonalization = false 
)
explicit

Constructor. By default, set the number of temporary vectors to 30, i.e. do a restart every 28 iterations. Also set preconditioning from left, the residual of the stopping criterion to the default residual, and re-orthogonalization only if necessary.

template<class VectorType = Vector<double>>
SolverGMRES< VectorType >::AdditionalData::AdditionalData ( const unsigned int  max_n_tmp_vectors,
const bool  right_preconditioning,
const bool  use_default_residual,
const bool  force_re_orthogonalization,
const bool  compute_eigenvalues 
)

Constructor.

Deprecated:
To obtain the estimated eigenvalues instead use: connect_eigenvalues_slot

Member Data Documentation

template<class VectorType = Vector<double>>
unsigned int SolverGMRES< VectorType >::AdditionalData::max_n_tmp_vectors

Maximum number of temporary vectors. This parameter controls the size of the Arnoldi basis, which for historical reasons is max_n_tmp_vectors-2.

Definition at line 210 of file solver_gmres.h.

template<class VectorType = Vector<double>>
bool SolverGMRES< VectorType >::AdditionalData::right_preconditioning

Flag for right preconditioning.

Note
Change between left and right preconditioning will also change the way residuals are evaluated. See the corresponding section in the SolverGMRES.

Definition at line 219 of file solver_gmres.h.

template<class VectorType = Vector<double>>
bool SolverGMRES< VectorType >::AdditionalData::use_default_residual

Flag for the default residual that is used to measure convergence.

Definition at line 224 of file solver_gmres.h.

template<class VectorType = Vector<double>>
bool SolverGMRES< VectorType >::AdditionalData::force_re_orthogonalization

Flag to force re-orthogonalization of orthonormal basis in every step. If set to false, the solver automatically checks for loss of orthogonality every 5 iterations and enables re-orthogonalization only if necessary.

Definition at line 232 of file solver_gmres.h.

template<class VectorType = Vector<double>>
bool SolverGMRES< VectorType >::AdditionalData::compute_eigenvalues

Compute all eigenvalues of the Hessenberg matrix generated while solving, i.e., the projected system matrix. This gives an approximation of the eigenvalues of the (preconditioned) system matrix. Since the Hessenberg matrix is thrown away at restart, the eigenvalues are printed for every 30 iterations.

Note
Requires LAPACK support.

Definition at line 243 of file solver_gmres.h.


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