Feel++  0.91.4
Protected Attributes
Feel::SolverNonLinear< T > Class Template Reference

Non linear solver base interface. More...

#include <solvernonlinear.hpp>

List of all members.

Public Types

Typedefs
typedef SolverNonLinear< T > self_type
typedef boost::shared_ptr
< SolverNonLinear< T > > 
self_ptrtype
typedef self_type solvernonlinear_type
typedef boost::shared_ptr
< self_type
solvernonlinear_ptrtype
typedef T value_type
typedef type_traits< T >::real_type real_type
typedef boost::shared_ptr
< Preconditioner< T > > 
preconditioner_ptrtype
typedef boost::shared_ptr
< Vector< value_type > > 
vector_ptrtype
typedef boost::shared_ptr
< MatrixSparse< value_type > > 
sparse_matrix_ptrtype
typedef ublas::matrix< value_type > dense_matrix_type
typedef ublas::vector< value_type > dense_vector_type
typedef boost::function< void(const
vector_ptrtype &X,
vector_ptrtype &R)> 
residual_function_type
typedef boost::function< void(const
vector_ptrtype &X,
sparse_matrix_ptrtype &J)> 
jacobian_function_type
typedef boost::function< void(const
vector_ptrtype &X,
vector_ptrtype &R,
sparse_matrix_ptrtype &J)> 
matvec_function_type
typedef boost::function< void(dense_vector_type
const &X, dense_vector_type &R)> 
dense_residual_function_type
typedef boost::function< void(dense_vector_type
const &X, dense_matrix_type &J)> 
dense_jacobian_function_type
typedef boost::function< void(dense_vector_type
const &X, dense_vector_type &R,
dense_matrix_type &J)> 
dense_matvec_function_type

Public Member Functions

Accessors
mpi::communicator const & comm () const
bool initialized () const
virtual void clear ()
virtual MatrixStructure precMatrixStructure () const
SolverNonLinearType getType () const
double getAbsoluteResidualTol () const
double getRelativeResidualTol () const
double getAbsoluteSolutionTol () const
uint getNbItMax () const
int reuseJacobian () const
int reusePreconditioner () const
Mutators
virtual void setPrecMatrixStructure (MatrixStructure mstruct)
void setType (SolverNonLinearType snl_type)
SolverNonLinearType nlSolverType () const
void setKspSolverType (const SolverType st)
SolverType kspSolverType () const
PreconditionerType preconditionerType () const
void setPreconditionerType (const PreconditionerType pct)
void attachPreconditioner (preconditioner_ptrtype preconditioner)
void setMatSolverPackageType (const MatSolverPackageType mspackt)
MatSolverPackageType matSolverPackageType () const
virtual void setReuse (int jac=1, int prec=1)
void setRelativeResidualTol (double tol)
void setAbsoluteResidualTol (double tol)
void setAbsoluteSolutionTol (double tol)
void setNbItMax (uint n)

Protected Attributes

mpi::communicator M_comm
bool M_is_initialized
MatrixStructure M_prec_matrix_structure
SolverNonLinearType M_snl_type
SolverType M_kspSolver_type
PreconditionerType M_preconditioner_type
preconditioner_ptrtype M_preconditioner
MatSolverPackageType M_matSolverPackage_type
double M_relativeResidualTol
double M_absoluteResidualTol
double M_absoluteSolutionTol
uint M_nbItMax
int M_reuse_jac
int M_reuse_prec

Constructors, destructor

 SolverNonLinear ()
 SolverNonLinear (SolverNonLinear const &)
virtual ~SolverNonLinear ()
virtual void init ()=0
static solvernonlinear_ptrtype build (po::variables_map const &vm, std::string const &prefix="")
static solvernonlinear_ptrtype build (SolverPackage solver_package)

Methods

residual_function_type residual
jacobian_function_type jacobian
matvec_function_type matvec
dense_residual_function_type dense_residual
dense_jacobian_function_type dense_jacobian
dense_matvec_function_type dense_matvec
virtual std::pair< int, real_type > solve (sparse_matrix_ptrtype &, vector_ptrtype &, vector_ptrtype &, const double, const unsigned int)=0
virtual std::pair< unsigned
int, real_type > 
solve (dense_matrix_type &, dense_vector_type &, dense_vector_type &, const double, const unsigned int)=0

Detailed Description

template<typename T>
class Feel::SolverNonLinear< T >

Non linear solver base interface.

This class provides a uniform interface for nonlinear solvers. This base class is overloaded to provide nonlinear solvers from different packages like PETSC.

Author:
Christophe Prud'homme

Constructor & Destructor Documentation

template<typename T >
Feel::SolverNonLinear< T >::SolverNonLinear ( ) [inline]

Constructor. Initializes Solver data structures

template<typename T >
Feel::SolverNonLinear< T >::SolverNonLinear ( SolverNonLinear< T > const &  snl) [inline]

copy constructor

template<typename T >
Feel::SolverNonLinear< T >::~SolverNonLinear ( ) [inline, virtual]

Destructor.


Member Function Documentation

template<typename T >
void Feel::SolverNonLinear< T >::attachPreconditioner ( preconditioner_ptrtype  preconditioner) [inline]
template<typename T >
boost::shared_ptr< SolverNonLinear< T > > Feel::SolverNonLinear< T >::build ( po::variables_map const &  vm,
std::string const &  prefix = "" 
) [static]

Builds a NonlinearSolver using the nonlinear solver package specified by the variables_map vm and prefix

template<typename T >
boost::shared_ptr< SolverNonLinear< T > > Feel::SolverNonLinear< T >::build ( SolverPackage  solver_package) [static]

Builds a NonlinearSolver using the nonlinear solver package specified by solver_package

template<typename T >
virtual void Feel::SolverNonLinear< T >::clear ( ) [inline, virtual]

Release all memory and clear data structures.

template<typename T >
mpi::communicator const& Feel::SolverNonLinear< T >::comm ( ) const [inline]
Returns:
the communicator
template<typename T >
virtual void Feel::SolverNonLinear< T >::init ( ) [pure virtual]

Initialize data structures if not done so already.

template<typename T >
bool Feel::SolverNonLinear< T >::initialized ( ) const [inline]
Returns:
true if the data structures are initialized, false otherwise.

References Feel::SolverNonLinear< T >::M_is_initialized.

template<typename T >
SolverType Feel::SolverNonLinear< T >::kspSolverType ( ) const [inline]

Returns the type of solver to use.

References Feel::SolverNonLinear< T >::M_kspSolver_type.

template<typename T >
MatSolverPackageType Feel::SolverNonLinear< T >::matSolverPackageType ( ) const [inline]

Returns the type of preconditioner to use.

References Feel::SolverNonLinear< T >::M_matSolverPackage_type.

template<typename T >
SolverNonLinearType Feel::SolverNonLinear< T >::nlSolverType ( ) const [inline]

Returns the type of solver to use.

References Feel::SolverNonLinear< T >::M_snl_type.

template<typename T >
virtual MatrixStructure Feel::SolverNonLinear< T >::precMatrixStructure ( ) const [inline, virtual]
Returns:
the preconditioner matrix structure it may not be relevant to all non linear solvers
template<typename T >
PreconditionerType Feel::SolverNonLinear< T >::preconditionerType ( ) const [inline]

Returns the type of preconditioner to use.

References Feel::SolverNonLinear< T >::M_preconditioner, and Feel::SolverNonLinear< T >::M_preconditioner_type.

template<typename T >
void Feel::SolverNonLinear< T >::setKspSolverType ( const SolverType  st) [inline]

Sets the type of solver to use.

References Feel::SolverNonLinear< T >::M_kspSolver_type.

template<typename T >
void Feel::SolverNonLinear< T >::setMatSolverPackageType ( const MatSolverPackageType  mspackt) [inline]

Sets the type of preconditioner to use.

References Feel::SolverNonLinear< T >::M_matSolverPackage_type.

template<typename T >
void Feel::SolverNonLinear< T >::setNbItMax ( uint  n) [inline]

Define the number max of iteration

References Feel::SolverNonLinear< T >::M_nbItMax.

template<typename T >
virtual void Feel::SolverNonLinear< T >::setPrecMatrixStructure ( MatrixStructure  mstruct) [inline, virtual]
Returns:
the preconditioner matrix structure it may not be relevant to all non linear solvers
template<typename T >
void Feel::SolverNonLinear< T >::setPreconditionerType ( const PreconditionerType  pct) [inline]
template<typename T >
void Feel::SolverNonLinear< T >::setRelativeResidualTol ( double  tol) [inline]

Define values of tolerance for the non linear solver

References Feel::SolverNonLinear< T >::M_relativeResidualTol.

template<typename T >
virtual void Feel::SolverNonLinear< T >::setReuse ( int  jac = 1,
int  prec = 1 
) [inline, virtual]

set reuse jacobian and/or preconditioner

  • jac=-1: means never rebuilt (preconditioner is not rebuit either)
  • jac=-2: build at each new nonlinear iterations (preconditioner is not rebuit either)
  • jac= 1: build at each new nonlinear iterations (rebuild jacobian every single nonlinear iteration)
  • jac= <n> (n>1): build every n iterations when jac >= n, prec=-1 (rebuilt once at each new nonlinear iteration) when jac >= n, prec=-1 (rebuilt once at each new nonlinear iteration)

References Feel::SolverNonLinear< T >::M_reuse_jac, and Feel::SolverNonLinear< T >::M_reuse_prec.

template<typename T >
void Feel::SolverNonLinear< T >::setType ( SolverNonLinearType  snl_type) [inline]

Select type of non linear solver : LINEAR_SEARCH, TRUST_REGION, ...

References Feel::SolverNonLinear< T >::M_snl_type.

template<typename T >
virtual std::pair<int, real_type> Feel::SolverNonLinear< T >::solve ( sparse_matrix_ptrtype &  ,
vector_ptrtype &  ,
vector_ptrtype &  ,
const double  ,
const unsigned  int 
) [pure virtual]

Solves a sparse nonlinear system.

template<typename T >
virtual std::pair<unsigned int, real_type> Feel::SolverNonLinear< T >::solve ( dense_matrix_type &  ,
dense_vector_type &  ,
dense_vector_type &  ,
const double  ,
const unsigned  int 
) [pure virtual]

Solves a sparse nonlinear system.


Member Data Documentation

template<typename T >
dense_jacobian_function_type Feel::SolverNonLinear< T >::dense_jacobian

Function that computes the Jacobian J(X) of the nonlinear system at the input iterate X.

template<typename T >
dense_matvec_function_type Feel::SolverNonLinear< T >::dense_matvec

Function that computes either the residual $ R(X) $ or the Jacobian $ J(X) $ of the nonlinear system at the input iterate $ X $. Note that either R or J could be XSNULL.

template<typename T >
dense_residual_function_type Feel::SolverNonLinear< T >::dense_residual

Function that computes the residual R(X) of the nonlinear system at the input iterate X.

template<typename T >
jacobian_function_type Feel::SolverNonLinear< T >::jacobian

Function that computes the Jacobian J(X) of the nonlinear system at the input iterate X.

template<typename T >
double Feel::SolverNonLinear< T >::M_absoluteSolutionTol [protected]

Absolute tolerances between successive iteration

template<typename T >
bool Feel::SolverNonLinear< T >::M_is_initialized [protected]

Flag indicating if the data structures have been initialized.

Referenced by Feel::SolverNonLinear< T >::attachPreconditioner(), and Feel::SolverNonLinear< T >::initialized().

template<typename T >
SolverType Feel::SolverNonLinear< T >::M_kspSolver_type [protected]

Enum stating which type of iterative linear solver to use.

Referenced by Feel::SolverNonLinear< T >::kspSolverType(), and Feel::SolverNonLinear< T >::setKspSolverType().

template<typename T >
MatSolverPackageType Feel::SolverNonLinear< T >::M_matSolverPackage_type [protected]

Enum the software that is used to perform the factorization

Referenced by Feel::SolverNonLinear< T >::matSolverPackageType(), and Feel::SolverNonLinear< T >::setMatSolverPackageType().

template<typename T >
uint Feel::SolverNonLinear< T >::M_nbItMax [protected]

The maximum numbers of allowable nonlinear iterations

Referenced by Feel::SolverNonLinear< T >::setNbItMax().

template<typename T >
preconditioner_ptrtype Feel::SolverNonLinear< T >::M_preconditioner [protected]
template<typename T >
PreconditionerType Feel::SolverNonLinear< T >::M_preconditioner_type [protected]
template<typename T >
double Feel::SolverNonLinear< T >::M_relativeResidualTol [protected]

Two differents tolerances on the residual for the resolution of non linear system

Referenced by Feel::SolverNonLinear< T >::setRelativeResidualTol().

template<typename T >
int Feel::SolverNonLinear< T >::M_reuse_jac [protected]

reuse jac level

Referenced by Feel::SolverNonLinear< T >::setReuse().

template<typename T >
int Feel::SolverNonLinear< T >::M_reuse_prec [protected]

reuse preconditioner level

Referenced by Feel::SolverNonLinear< T >::setReuse().

template<typename T >
SolverNonLinearType Feel::SolverNonLinear< T >::M_snl_type [protected]

Define the type of non linear solver

Referenced by Feel::SolverNonLinear< T >::nlSolverType(), and Feel::SolverNonLinear< T >::setType().

template<typename T >
matvec_function_type Feel::SolverNonLinear< T >::matvec

Function that computes either the residual $ R(X) $ or the Jacobian $ J(X) $ of the nonlinear system at the input iterate $ X $. Note that either R or J could be XSNULL.

template<typename T >
residual_function_type Feel::SolverNonLinear< T >::residual

Function that computes the residual R(X) of the nonlinear system at the input iterate X.