MFEM  v4.0
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
Public Member Functions | Protected Attributes | List of all members
mfem::KinSolver Class Reference

Wrapper for SUNDIALS' KINSOL library – Nonlinear solvers. More...

#include <sundials.hpp>

Inheritance diagram for mfem::KinSolver:
[legend]
Collaboration diagram for mfem::KinSolver:
[legend]

Public Member Functions

 KinSolver (int strategy, bool oper_grad=true)
 Construct a serial KinSolver, a wrapper for SUNDIALS' KINSOL solver. More...
 
 KinSolver (MPI_Comm comm, int strategy, bool oper_grad=true)
 Construct a parallel KinSolver, a wrapper for SUNDIALS' KINSOL solver. More...
 
virtual ~KinSolver ()
 Destroy the associated KINSOL memory. More...
 
virtual void SetOperator (const Operator &op)
 Set the nonlinear Operator of the system. This method calls KINInit(). More...
 
virtual void SetSolver (Solver &solver)
 Set the linear solver for inverting the Jacobian. More...
 
virtual void SetPreconditioner (Solver &solver)
 Equivalent to SetSolver(Solver). More...
 
void SetScaledStepTol (double sstol)
 Set KINSOL's scaled step tolerance. More...
 
void SetFuncNormTol (double ftol)
 Set KINSOL's functional norm tolerance. More...
 
void SetMaxSetupCalls (int max_calls)
 Set maximum number of nonlinear iterations without a Jacobian update. More...
 
virtual void Mult (const Vector &b, Vector &x) const
 Solve the nonlinear system F(x) = 0. More...
 
void Mult (Vector &x, const Vector &x_scale, const Vector &fx_scale) const
 Solve the nonlinear system F(x) = 0. More...
 
- Public Member Functions inherited from mfem::NewtonSolver
 NewtonSolver ()
 
 NewtonSolver (MPI_Comm _comm)
 
virtual double ComputeScalingFactor (const Vector &x, const Vector &b) const
 This method can be overloaded in derived classes to implement line search algorithms. More...
 
- Public Member Functions inherited from mfem::IterativeSolver
 IterativeSolver ()
 
 IterativeSolver (MPI_Comm _comm)
 
void SetRelTol (double rtol)
 
void SetAbsTol (double atol)
 
void SetMaxIter (int max_it)
 
void SetPrintLevel (int print_lvl)
 
int GetNumIterations () const
 
int GetConverged () const
 
double GetFinalNorm () const
 
- Public Member Functions inherited from mfem::Solver
 Solver (int s=0, bool iter_mode=false)
 Initialize a square Solver with size s. More...
 
 Solver (int h, int w, bool iter_mode=false)
 Initialize a Solver with height h and width w. More...
 
- Public Member Functions inherited from mfem::Operator
 Operator (int s=0)
 Construct a square Operator with given size s (default 0). More...
 
 Operator (int h, int w)
 Construct an Operator with the given height (output size) and width (input size). More...
 
int Height () const
 Get the height (size of output) of the Operator. Synonym with NumRows(). More...
 
int NumRows () const
 Get the number of rows (size of output) of the Operator. Synonym with Height(). More...
 
int Width () const
 Get the width (size of input) of the Operator. Synonym with NumCols(). More...
 
int NumCols () const
 Get the number of columns (size of input) of the Operator. Synonym with Width(). More...
 
virtual MemoryClass GetMemoryClass () const
 Return the MemoryClass preferred by the Operator. More...
 
virtual void MultTranspose (const Vector &x, Vector &y) const
 Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an error. More...
 
virtual OperatorGetGradient (const Vector &x) const
 Evaluate the gradient operator at the point x. The default behavior in class Operator is to generate an error. More...
 
virtual const OperatorGetProlongation () const
 Prolongation operator from linear algebra (linear system) vectors, to input vectors for the operator. NULL means identity. More...
 
virtual const OperatorGetRestriction () const
 Restriction operator from input vectors for the operator to linear algebra (linear system) vectors. NULL means identity. More...
 
void FormLinearSystem (const Array< int > &ess_tdof_list, Vector &x, Vector &b, Operator *&A, Vector &X, Vector &B, int copy_interior=0)
 Form a constrained linear system using a matrix-free approach. More...
 
virtual void RecoverFEMSolution (const Vector &X, const Vector &b, Vector &x)
 Reconstruct a solution vector x (e.g. a GridFunction) from the solution X of a constrained linear system obtained from Operator::FormLinearSystem(). More...
 
void PrintMatlab (std::ostream &out, int n=0, int m=0) const
 Prints operator with input size n and output size m in Matlab format. More...
 
virtual ~Operator ()
 Virtual destructor. More...
 
Type GetType () const
 Return the type ID of the Operator class. More...
 
- Public Member Functions inherited from mfem::SundialsSolver
void * SundialsMem () const
 Access the underlying SUNDIALS object. More...
 
int GetFlag () const
 Return the flag returned by the last call to a SUNDIALS function. More...
 

Static Protected Member Functions

Auxiliary callback functions.
static int Mult (const N_Vector u, N_Vector fu, void *user_data)
 
static int GradientMult (N_Vector v, N_Vector Jv, N_Vector u, booleantype *new_u, void *user_data)
 
static int LinSysSetup (KINMemRec *kin_mem)
 
static int LinSysSolve (KINMemRec *kin_mem, N_Vector x, N_Vector b, realtype *sJpnorm, realtype *sFdotJp)
 
- Static Protected Member Functions inherited from mfem::SundialsSolver
static int ODEMult (realtype t, const N_Vector y, N_Vector ydot, void *td_oper)
 Callback function used in CVODESolver and ARKODESolver. More...
 

Protected Attributes

bool use_oper_grad
 
N_Vector y_scale
 
N_Vector f_scale
 
const Operatorjacobian
 
- Protected Attributes inherited from mfem::NewtonSolver
Vector r
 
Vector c
 
- Protected Attributes inherited from mfem::IterativeSolver
const Operatoroper
 
Solverprec
 
int max_iter
 
int print_level
 
double rel_tol
 
double abs_tol
 
int final_iter
 
int converged
 
double final_norm
 
- Protected Attributes inherited from mfem::Operator
int height
 Dimension of the output / number of rows in the matrix. More...
 
int width
 Dimension of the input / number of columns in the matrix. More...
 
- Protected Attributes inherited from mfem::SundialsSolver
void * sundials_mem
 Pointer to the SUNDIALS mem object. More...
 
int flag
 Flag returned by the last call to SUNDIALS. More...
 
N_Vector y
 Auxiliary N_Vector. More...
 

Additional Inherited Members

- Public Types inherited from mfem::Operator
enum  Type {
  ANY_TYPE, MFEM_SPARSEMAT, Hypre_ParCSR, PETSC_MATAIJ,
  PETSC_MATIS, PETSC_MATSHELL, PETSC_MATNEST, PETSC_MATHYPRE,
  PETSC_MATGENERIC
}
 Enumeration defining IDs for some classes derived from Operator. More...
 
- Public Attributes inherited from mfem::Solver
bool iterative_mode
 If true, use the second argument of Mult() as an initial guess. More...
 
- Protected Member Functions inherited from mfem::IterativeSolver
double Dot (const Vector &x, const Vector &y) const
 
double Norm (const Vector &x) const
 
- Protected Member Functions inherited from mfem::SundialsSolver
bool Parallel () const
 
bool Parallel () const
 
 SundialsSolver ()
 
 SundialsSolver (void *mem)
 
- Static Protected Attributes inherited from mfem::SundialsSolver
static const double default_rel_tol = 1e-4
 
static const double default_abs_tol = 1e-9
 

Detailed Description

Wrapper for SUNDIALS' KINSOL library – Nonlinear solvers.

Note
To minimize uncertainty, we advise the user to adhere to the given interface, instead of making similar calls by the KINSOL's internal KINMem object.

Definition at line 312 of file sundials.hpp.

Constructor & Destructor Documentation

mfem::KinSolver::KinSolver ( int  strategy,
bool  oper_grad = true 
)

Construct a serial KinSolver, a wrapper for SUNDIALS' KINSOL solver.

Parameters
[in]strategySpecifies the nonlinear solver strategy: KIN_NONE / KIN_LINESEARCH / KIN_PICARD / KIN_FP.
[in]oper_gradSpecifies whether the solver should use its Operator's GetGradient() method to compute the Jacobian of the system.

Definition at line 788 of file sundials.cpp.

mfem::KinSolver::KinSolver ( MPI_Comm  comm,
int  strategy,
bool  oper_grad = true 
)

Construct a parallel KinSolver, a wrapper for SUNDIALS' KINSOL solver.

Parameters
[in]commThe MPI communicator used to partition the system.
[in]strategySpecifies the nonlinear solver strategy: KIN_NONE / KIN_LINESEARCH / KIN_PICARD / KIN_FP.
[in]oper_gradSpecifies whether the solver should use its Operator's GetGradient() method to compute the Jacobian of the system.

Definition at line 810 of file sundials.cpp.

mfem::KinSolver::~KinSolver ( )
virtual

Destroy the associated KINSOL memory.

Definition at line 1074 of file sundials.cpp.

Member Function Documentation

int mfem::KinSolver::GradientMult ( N_Vector  v,
N_Vector  Jv,
N_Vector  u,
booleantype *  new_u,
void *  user_data 
)
staticprotected

Definition at line 731 of file sundials.cpp.

int mfem::KinSolver::LinSysSetup ( KINMemRec *  kin_mem)
staticprotected

Definition at line 748 of file sundials.cpp.

int mfem::KinSolver::LinSysSolve ( KINMemRec *  kin_mem,
N_Vector  x,
N_Vector  b,
realtype *  sJpnorm,
realtype *  sFdotJp 
)
staticprotected

Definition at line 761 of file sundials.cpp.

int mfem::KinSolver::Mult ( const N_Vector  u,
N_Vector  fu,
void *  user_data 
)
staticprotected

Definition at line 720 of file sundials.cpp.

void mfem::KinSolver::Mult ( const Vector b,
Vector x 
) const
virtual

Solve the nonlinear system F(x) = 0.

Calls the other Mult(Vector&, Vector&, Vector&) const method with x_scale = 1. The values of 'fx_scale' are determined by comparing the chosen relative and functional norm (i.e. absolute) tolerances.

Parameters
[in]bNot used, KINSol always assumes zero RHS.
[in,out]xOn input, initial guess, if iterative_mode = true, otherwise the initial guess is zero; on output, the solution.

Reimplemented from mfem::NewtonSolver.

Definition at line 985 of file sundials.cpp.

void mfem::KinSolver::Mult ( Vector x,
const Vector x_scale,
const Vector fx_scale 
) const

Solve the nonlinear system F(x) = 0.

Calls KINSol() to solve the nonlinear system. Before calling KINSol(), this functions uses the data members inherited from class IterativeSolver to set corresponding KINSOL options.

Parameters
[in,out]xOn input, initial guess, if iterative_mode = true, otherwise the initial guess is zero; on output, the solution.
[in]x_scaleElements of a diagonal scaling matrix D, s.t. D*x has all elements roughly the same when x is close to a solution.
[in]fx_scaleElements of a diagonal scaling matrix E, s.t. D*F(x) has all elements roughly the same when x is not too close to a solution.

Definition at line 1031 of file sundials.cpp.

void mfem::KinSolver::SetFuncNormTol ( double  ftol)
inline

Set KINSOL's functional norm tolerance.

The default tolerance is U^(1/3), where U = machine unit roundoff.

Note
This function is equivalent to SetAbsTol(double).

Definition at line 375 of file sundials.hpp.

void mfem::KinSolver::SetMaxSetupCalls ( int  max_calls)

Set maximum number of nonlinear iterations without a Jacobian update.

The default is 10.

Definition at line 980 of file sundials.cpp.

void mfem::KinSolver::SetOperator ( const Operator op)
virtual

Set the nonlinear Operator of the system. This method calls KINInit().

Reimplemented from mfem::NewtonSolver.

Definition at line 863 of file sundials.cpp.

virtual void mfem::KinSolver::SetPreconditioner ( Solver solver)
inlinevirtual

Equivalent to SetSolver(Solver).

Reimplemented from mfem::IterativeSolver.

Definition at line 367 of file sundials.hpp.

void mfem::KinSolver::SetScaledStepTol ( double  sstol)

Set KINSOL's scaled step tolerance.

The default tolerance is U^(2/3), where U = machine unit roundoff.

Definition at line 975 of file sundials.cpp.

void mfem::KinSolver::SetSolver ( Solver solver)
virtual

Set the linear solver for inverting the Jacobian.

Note
This function assumes that Operator::GetGradient(const Vector &) is implemented by the Operator specified by SetOperator(const Operator &).

Reimplemented from mfem::NewtonSolver.

Definition at line 958 of file sundials.cpp.

Member Data Documentation

N_Vector mfem::KinSolver::f_scale
mutableprotected

Definition at line 316 of file sundials.hpp.

const Operator* mfem::KinSolver::jacobian
protected

Definition at line 317 of file sundials.hpp.

bool mfem::KinSolver::use_oper_grad
protected

Definition at line 315 of file sundials.hpp.

N_Vector mfem::KinSolver::y_scale
mutableprotected

Definition at line 316 of file sundials.hpp.


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