Hex  1.0
Hydrogen-electron collision solver
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Member Functions | Static Public Attributes | Protected Attributes
ILUCGPreconditioner Class Reference

ILU-preconditioned CG-based preconditioner. More...

#include <preconditioners.h>

Inheritance diagram for ILUCGPreconditioner:
Inheritance graph
[legend]
Collaboration diagram for ILUCGPreconditioner:
Collaboration graph
[legend]

Public Member Functions

 ILUCGPreconditioner (Parallel const &par, InputFile const &inp, std::vector< std::pair< int, int >> const &ll, Bspline const &bspline, CommandLine const &cmd)
 
virtual RadialIntegrals const & rad () const
 Get radial integrals. More...
 
virtual void multiply (const cArrayView p, cArrayView q) const
 Multiply by the matrix equation. More...
 
virtual void rhs (cArrayView chi, int ienergy, int instate) const
 Return the right-hand side. More...
 
virtual void precondition (const cArrayView r, cArrayView z) const
 Precondition the equation. More...
 
virtual void setup ()
 Initialize the preconditioner. More...
 
virtual void update (double E)
 Update the preconditioner for the next energy. More...
 
virtual void CG_prec (int iblock, const cArrayView r, cArrayView z) const
 
- Public Member Functions inherited from CGPreconditioner
 CGPreconditioner (Parallel const &par, InputFile const &inp, std::vector< std::pair< int, int >> const &ll, Bspline const &bspline, CommandLine const &cmd)
 
virtual void CG_mmul (int iblock, const cArrayView p, cArrayView q) const
 
- Public Member Functions inherited from NoPreconditioner
 NoPreconditioner (Parallel const &par, InputFile const &inp, std::vector< std::pair< int, int >> const &ll, Bspline const &bspline, CommandLine const &cmd)
 

Static Public Attributes

static const std::string name = "ILU"
 
static const std::string description = "Block inversion using ILU-preconditioned conjugate gradients. The drop tolerance can be given as the --droptol parameter."
 
- Static Public Attributes inherited from CGPreconditioner
static const std::string name = "cg"
 
static const std::string description = "Block inversion using plain conjugate gradients."
 
- Static Public Attributes inherited from NoPreconditioner
static const std::string name = "none"
 
static const std::string description = "Preconditioning by the identity matrix."
 

Protected Attributes

double droptol_
 
std::vector< CsrMatrixcsr_blocks_
 
std::vector< CsrMatrix::LUftlu_
 
- Protected Attributes inherited from NoPreconditioner
CommandLine const & cmd_
 
Parallel const & par_
 
InputFile const & inp_
 
std::vector< std::pair< int,
int > > const & 
l1_l2_
 
std::vector< SymDiaMatrixdia_blocks_
 
Bspline s_bspline_
 
RadialIntegrals s_rad_
 
SymDiaMatrix S_kron_S_
 
SymDiaMatrix S_kron_Mm1_tr_
 
SymDiaMatrix S_kron_Mm2_
 
SymDiaMatrix Mm1_tr_kron_S_
 
SymDiaMatrix Mm2_kron_S_
 
SymDiaMatrix half_D_minus_Mm1_tr_
 
SymDiaMatrix half_D_minus_Mm1_tr_kron_S_
 
SymDiaMatrix S_kron_half_D_minus_Mm1_tr_
 

Detailed Description

Enhances CGPreconditioner conjugate gradients solver by incomplete LU factorization preconditioning. This is done by redefining virtual function CG_prec. The factorization is drop tolerance based and is computed by UMFPACK.

Constructor & Destructor Documentation

ILUCGPreconditioner::ILUCGPreconditioner ( Parallel const &  par,
InputFile const &  inp,
std::vector< std::pair< int, int >> const &  ll,
Bspline const &  bspline,
CommandLine const &  cmd 
)
inline

Member Function Documentation

void ILUCGPreconditioner::CG_prec ( int  iblock,
const cArrayView  r,
cArrayView  z 
) const
virtual

Reimplemented from CGPreconditioner.

virtual void ILUCGPreconditioner::multiply ( const cArrayView  p,
cArrayView  q 
) const
inlinevirtual

This function implements matrix multiplication by the matrix of the set of equations that is to be solved.

Reimplemented from CGPreconditioner.

virtual void ILUCGPreconditioner::precondition ( const cArrayView  r,
cArrayView  z 
) const
inlinevirtual

This function preconditions the equation, solving the preconditioner equation

\[ \mathbf{M}\mathbf{z} = \mathbf{r} \ . \]

It may use the MPI environment.

Reimplemented from CGPreconditioner.

virtual RadialIntegrals const& ILUCGPreconditioner::rad ( ) const
inlinevirtual

Reimplemented from CGPreconditioner.

virtual void ILUCGPreconditioner::rhs ( cArrayView  chi,
int  ienergy,
int  instate 
) const
inlinevirtual

Reimplemented from CGPreconditioner.

void ILUCGPreconditioner::setup ( )
virtual

This function contains all computation intensive preparations for the preconditioner, e.g. computation of radial integrals. It may use only SMP environment.

Reimplemented from CGPreconditioner.

void ILUCGPreconditioner::update ( double  E)
virtual

This function updates the preconditioner for another right hand side. It may use the MPI environment.

Reimplemented from CGPreconditioner.

Field Documentation

std::vector<CsrMatrix> ILUCGPreconditioner::csr_blocks_
mutableprotected
const std::string ILUCGPreconditioner::description = "Block inversion using ILU-preconditioned conjugate gradients. The drop tolerance can be given as the --droptol parameter."
static
double ILUCGPreconditioner::droptol_
protected
std::vector<CsrMatrix::LUft> ILUCGPreconditioner::lu_
mutableprotected
const std::string ILUCGPreconditioner::name = "ILU"
static

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