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

Base abstract class for time dependent operators. More...

#include <operator.hpp>

Inherits mfem::Operator.

Inherited by ConductionOperator, ConductionOperator, ConductionOperator, ConductionOperator, FE_Evolution, FE_Evolution, FE_Evolution, FE_Evolution, FE_Evolution, FE_Evolution, HyperelasticOperator, HyperelasticOperator, HyperelasticOperator, HyperelasticOperator, HyperelasticOperator, and mfem::electromagnetics::MagneticDiffusionEOperator.

Collaboration diagram for mfem::TimeDependentOperator:
[legend]

Public Types

enum  Type { EXPLICIT, IMPLICIT, HOMOGENEOUS }
 
- Public Types inherited from mfem::Operator
enum  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 Member Functions

 TimeDependentOperator (int n=0, double t_=0.0, Type type_=EXPLICIT)
 Construct a "square" TimeDependentOperator y = f(x,t), where x and y have the same dimension n. More...
 
 TimeDependentOperator (int h, int w, double t_=0.0, Type type_=EXPLICIT)
 Construct a TimeDependentOperator y = f(x,t), where x and y have dimensions w and h, respectively. More...
 
virtual double GetTime () const
 Read the currently set time. More...
 
virtual void SetTime (const double _t)
 Set the current time. More...
 
bool isExplicit () const
 True if type is EXPLICIT. More...
 
bool isImplicit () const
 True if type is IMPLICIT or HOMOGENEOUS. More...
 
bool isHomogeneous () const
 True if type is HOMOGENEOUS. More...
 
virtual void ExplicitMult (const Vector &x, Vector &y) const
 Perform the action of the explicit part of the operator, G: y = G(x, t) where t is the current time. More...
 
virtual void ImplicitMult (const Vector &x, const Vector &k, Vector &y) const
 Perform the action of the implicit part of the operator, F: y = F(x, k, t) where t is the current time. More...
 
virtual void Mult (const Vector &x, Vector &y) const
 Perform the action of the operator: y = k = f(x, t), where k solves the algebraic equation F(x, k, t) = G(x, t) and t is the current time. More...
 
virtual void ImplicitSolve (const double dt, const Vector &x, Vector &k)
 Solve the equation: k = f(x + dt k, t), for the unknown k at the current time t. More...
 
virtual OperatorGetImplicitGradient (const Vector &x, const Vector &k, double shift) const
 Return an Operator representing (dF/dk shift + dF/dx) at the given x, k, and the currently set time. More...
 
virtual OperatorGetExplicitGradient (const Vector &x) const
 Return an Operator representing dG/dx at the given point x and the currently set time. More...
 
virtual ~TimeDependentOperator ()
 
- 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 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...
 

Protected Attributes

double t
 Current time. More...
 
Type type
 Describes the form of the TimeDependentOperator. More...
 
- 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...
 

Detailed Description

Base abstract class for time dependent operators.

Operator of the form: (x,t) -> f(x,t), where k = f(x,t) generally solves the algebraic equation F(x,k,t) = G(x,t). The functions F and G represent the implicit and explicit parts of the operator, respectively. For explicit operators, F(x,k,t) = k, so f(x,t) = G(x,t).

Definition at line 142 of file operator.hpp.

Member Enumeration Documentation

Enumerator
EXPLICIT 

This type assumes F(x,k,t) = k, i.e. k = f(x,t) = G(x,t).

IMPLICIT 

This is the most general type, no assumptions on F and G.

HOMOGENEOUS 

This type assumes that G(x,t) = 0.

Definition at line 145 of file operator.hpp.

Constructor & Destructor Documentation

mfem::TimeDependentOperator::TimeDependentOperator ( int  n = 0,
double  t_ = 0.0,
Type  type_ = EXPLICIT 
)
inlineexplicit

Construct a "square" TimeDependentOperator y = f(x,t), where x and y have the same dimension n.

Definition at line 159 of file operator.hpp.

mfem::TimeDependentOperator::TimeDependentOperator ( int  h,
int  w,
double  t_ = 0.0,
Type  type_ = EXPLICIT 
)
inline

Construct a TimeDependentOperator y = f(x,t), where x and y have dimensions w and h, respectively.

Definition at line 165 of file operator.hpp.

virtual mfem::TimeDependentOperator::~TimeDependentOperator ( )
inlinevirtual

Definition at line 255 of file operator.hpp.

Member Function Documentation

virtual void mfem::TimeDependentOperator::ExplicitMult ( const Vector x,
Vector y 
) const
inlinevirtual

Perform the action of the explicit part of the operator, G: y = G(x, t) where t is the current time.

Presently, this method is used by some PETSc ODE solvers, for more details, see the PETSc Manual.

Reimplemented in FE_Evolution.

Definition at line 186 of file operator.hpp.

virtual Operator& mfem::TimeDependentOperator::GetExplicitGradient ( const Vector x) const
inlinevirtual

Return an Operator representing dG/dx at the given point x and the currently set time.

Presently, this method is used by some PETSc ODE solvers, for more details, see the PETSc Manual.

Reimplemented in FE_Evolution.

Definition at line 248 of file operator.hpp.

virtual Operator& mfem::TimeDependentOperator::GetImplicitGradient ( const Vector x,
const Vector k,
double  shift 
) const
inlinevirtual

Return an Operator representing (dF/dk shift + dF/dx) at the given x, k, and the currently set time.

Presently, this method is used by some PETSc ODE solvers, for more details, see the PETSc Manual.

Reimplemented in FE_Evolution.

Definition at line 235 of file operator.hpp.

virtual double mfem::TimeDependentOperator::GetTime ( ) const
inlinevirtual

Read the currently set time.

Definition at line 169 of file operator.hpp.

virtual void mfem::TimeDependentOperator::ImplicitMult ( const Vector x,
const Vector k,
Vector y 
) const
inlinevirtual

Perform the action of the implicit part of the operator, F: y = F(x, k, t) where t is the current time.

Presently, this method is used by some PETSc ODE solvers, for more details, see the PETSc Manual.

Reimplemented in FE_Evolution.

Definition at line 196 of file operator.hpp.

virtual void mfem::TimeDependentOperator::ImplicitSolve ( const double  dt,
const Vector x,
Vector k 
)
inlinevirtual

Solve the equation: k = f(x + dt k, t), for the unknown k at the current time t.

For general F and G, the equation for k becomes: F(x + dt k, k, t) = G(x + dt k, t).

The input vector x corresponds to time index (or cycle) n, while the currently set time, t, and the result vector k correspond to time index n+1. The time step dt corresponds to the time interval between cycles n and n+1.

This method allows for the abstract implementation of some time integration methods, including diagonal implicit Runge-Kutta (DIRK) methods and the backward Euler method in particular.

If not re-implemented, this method simply generates an error.

Reimplemented in mfem::electromagnetics::MagneticDiffusionEOperator.

Definition at line 225 of file operator.hpp.

bool mfem::TimeDependentOperator::isExplicit ( ) const
inline

True if type is EXPLICIT.

Definition at line 175 of file operator.hpp.

bool mfem::TimeDependentOperator::isHomogeneous ( ) const
inline

True if type is HOMOGENEOUS.

Definition at line 179 of file operator.hpp.

bool mfem::TimeDependentOperator::isImplicit ( ) const
inline

True if type is IMPLICIT or HOMOGENEOUS.

Definition at line 177 of file operator.hpp.

virtual void mfem::TimeDependentOperator::Mult ( const Vector x,
Vector y 
) const
inlinevirtual

Perform the action of the operator: y = k = f(x, t), where k solves the algebraic equation F(x, k, t) = G(x, t) and t is the current time.

Implements mfem::Operator.

Reimplemented in mfem::electromagnetics::MagneticDiffusionEOperator, FE_Evolution, FE_Evolution, FE_Evolution, FE_Evolution, FE_Evolution, and FE_Evolution.

Definition at line 204 of file operator.hpp.

virtual void mfem::TimeDependentOperator::SetTime ( const double  _t)
inlinevirtual

Set the current time.

Reimplemented in mfem::electromagnetics::MagneticDiffusionEOperator.

Definition at line 172 of file operator.hpp.

Member Data Documentation

double mfem::TimeDependentOperator::t
protected

Current time.

Definition at line 153 of file operator.hpp.

Type mfem::TimeDependentOperator::type
protected

Describes the form of the TimeDependentOperator.

Definition at line 154 of file operator.hpp.


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