MFEM  v4.6.0
Finite element discretization library
Public Member Functions | List of all members
mfem::ConvergenceStudy Class Reference

Class to compute error and convergence rates. It supports H1, H(curl) (ND elements), H(div) (RT elements) and L2 (DG). More...

#include <convergence.hpp>

Public Member Functions

void Reset ()
 Clear any internal data. More...
 
void AddL2GridFunction (GridFunction *gf, Coefficient *scalar_u, VectorCoefficient *grad=nullptr, Coefficient *ell_coeff=nullptr, JumpScaling jump_scaling={1.0, JumpScaling::ONE_OVER_H})
 
void AddH1GridFunction (GridFunction *gf, Coefficient *scalar_u, VectorCoefficient *grad=nullptr)
 Add H1 GridFunction, the exact solution and possibly its gradient. More...
 
void AddHcurlGridFunction (GridFunction *gf, VectorCoefficient *vector_u, VectorCoefficient *curl=nullptr)
 Add H(curl) GridFunction, the exact solution and possibly its curl. More...
 
void AddHdivGridFunction (GridFunction *gf, VectorCoefficient *vector_u, Coefficient *div=nullptr)
 Add H(div) GridFunction, the exact solution and possibly its div. More...
 
double GetL2Error (int n)
 Get the L2 error at step n. More...
 
void GetL2Errors (Array< double > &L2Errors_)
 Get all L2 errors. More...
 
double GetDError (int n)
 Get the Grad/Curl/Div error at step n. More...
 
void GetDErrors (Array< double > &DErrors_)
 Get all Grad/Curl/Div errors. More...
 
double GetDGFaceJumpsError (int n)
 Get the DGFaceJumps error at step n. More...
 
void GetDGFaceJumpsErrors (Array< double > &DGFaceErrors_)
 Get all DGFaceJumps errors. More...
 
void Print (bool relative=false, std::ostream &out=mfem::out)
 Print rates and errors. More...
 

Detailed Description

Class to compute error and convergence rates. It supports H1, H(curl) (ND elements), H(div) (RT elements) and L2 (DG).

For "smooth enough" solutions the Galerkin error measured in the appropriate norm satisfies || u - u_h || ~ h^k

Here, k is called the asymptotic rate of convergence

For successive h-refinements the rate can be estimated by k = log(||u - u_h|| / ||u - u_{h/2}||)/(1/dim * log(N_{h/2}/N_{h})

Definition at line 35 of file convergence.hpp.

Member Function Documentation

◆ AddH1GridFunction()

void mfem::ConvergenceStudy::AddH1GridFunction ( GridFunction gf,
Coefficient scalar_u,
VectorCoefficient grad = nullptr 
)
inline

Add H1 GridFunction, the exact solution and possibly its gradient.

Definition at line 86 of file convergence.hpp.

◆ AddHcurlGridFunction()

void mfem::ConvergenceStudy::AddHcurlGridFunction ( GridFunction gf,
VectorCoefficient vector_u,
VectorCoefficient curl = nullptr 
)
inline

Add H(curl) GridFunction, the exact solution and possibly its curl.

Definition at line 93 of file convergence.hpp.

◆ AddHdivGridFunction()

void mfem::ConvergenceStudy::AddHdivGridFunction ( GridFunction gf,
VectorCoefficient vector_u,
Coefficient div = nullptr 
)
inline

Add H(div) GridFunction, the exact solution and possibly its div.

Definition at line 100 of file convergence.hpp.

◆ AddL2GridFunction()

void mfem::ConvergenceStudy::AddL2GridFunction ( GridFunction gf,
Coefficient scalar_u,
VectorCoefficient grad = nullptr,
Coefficient ell_coeff = nullptr,
JumpScaling  jump_scaling = {1.0, JumpScaling::ONE_OVER_H} 
)
inline

Add L2 GridFunction, the exact solution and possibly its gradient and/or DG face jumps parameters

Definition at line 77 of file convergence.hpp.

◆ GetDError()

double mfem::ConvergenceStudy::GetDError ( int  n)
inline

Get the Grad/Curl/Div error at step n.

Definition at line 120 of file convergence.hpp.

◆ GetDErrors()

void mfem::ConvergenceStudy::GetDErrors ( Array< double > &  DErrors_)
inline

Get all Grad/Curl/Div errors.

Definition at line 127 of file convergence.hpp.

◆ GetDGFaceJumpsError()

double mfem::ConvergenceStudy::GetDGFaceJumpsError ( int  n)
inline

Get the DGFaceJumps error at step n.

Definition at line 133 of file convergence.hpp.

◆ GetDGFaceJumpsErrors()

void mfem::ConvergenceStudy::GetDGFaceJumpsErrors ( Array< double > &  DGFaceErrors_)
inline

Get all DGFaceJumps errors.

Definition at line 140 of file convergence.hpp.

◆ GetL2Error()

double mfem::ConvergenceStudy::GetL2Error ( int  n)
inline

Get the L2 error at step n.

Definition at line 107 of file convergence.hpp.

◆ GetL2Errors()

void mfem::ConvergenceStudy::GetL2Errors ( Array< double > &  L2Errors_)
inline

Get all L2 errors.

Definition at line 114 of file convergence.hpp.

◆ Print()

void mfem::ConvergenceStudy::Print ( bool  relative = false,
std::ostream &  out = mfem::out 
)

Print rates and errors.

Definition at line 226 of file convergence.cpp.

◆ Reset()

void mfem::ConvergenceStudy::Reset ( )

Clear any internal data.

Definition at line 9 of file convergence.cpp.


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