12#ifndef MFEM_ERROR_ESTIMATORS
13#define MFEM_ERROR_ESTIMATORS
296 MFEM_VERIFY(tcoeff >= 0.0,
"Tichonov coefficient cannot be negative");
564 std::function<
real_t(
Mesh*,
const int,
const bool)>;
567 int current_sequence = -1;
599 bool own_flux_fespace;
602 const bool isParallel;
606 bool MeshIsModified()
609 MFEM_ASSERT(mesh_sequence >= current_sequence,
610 "improper mesh update sequence");
611 return (mesh_sequence > current_sequence);
622 void ComputeEstimates();
654 if (MeshIsModified())
658 return error_estimates;
662 void Reset()
override { current_sequence = -1; };
672 compute_element_coefficient_)
674 compute_element_coefficient = compute_element_coefficient_;
684 compute_face_coefficient_)
686 compute_face_coefficient = compute_face_coefficient_;
Base class for all error estimators.
virtual ~AbstractErrorEstimator()
The AnisotropicErrorEstimator class is the base class for all error estimators that compute one non-n...
virtual const Array< int > & GetAnisotropicFlags()=0
Get an Array<int> with anisotropic flags for all mesh elements.
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
Base class for all element based error estimators.
virtual real_t GetTotalError() const
Return the total error from the last error estimate.
virtual void Reset()=0
Force recomputation of the estimates on the next call to GetLocalErrors.
virtual ~ErrorEstimator()
Destruct the error estimator.
virtual const Vector & GetLocalErrors()=0
Get a Vector with all element errors.
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Mesh * GetMesh() const
Returns the mesh.
Class for grid function - Vector with associated FE space.
FiniteElementSpace * FESpace()
The KellyErrorEstimator class provides a fast error indication strategy for smooth scalar parallel pr...
void ResetCoefficientFunctions()
Change the coefficients back to default as described above.
void Reset() override
Reset the error estimator.
const Vector & GetLocalErrors() override
Get a Vector with all element errors.
std::function< real_t(Mesh *, const int)> ElementCoefficientFunction
Function type to compute the local coefficient hₑ of an element.
void SetElementCoefficientFunction(ElementCoefficientFunction compute_element_coefficient_)
Change the method to compute hₑ on a per-element basis.
KellyErrorEstimator(BilinearFormIntegrator &di_, GridFunction &sol_, FiniteElementSpace &flux_fes_, const Array< int > &attributes_=Array< int >())
Construct a new KellyErrorEstimator object for a scalar field.
real_t GetTotalError() const override
Return the total error from the last error estimate.
std::function< real_t(Mesh *, const int, const bool)> FaceCoefficientFunction
Function type to compute the local coefficient hₖ of a face. The third argument is true for shared fa...
void SetFaceCoefficientFunction(FaceCoefficientFunction compute_face_coefficient_)
Change the method to compute hₖ on a per-element basis.
The L2ZienkiewiczZhuEstimator class implements the Zienkiewicz-Zhu error estimation procedure where t...
L2ZienkiewiczZhuEstimator(BilinearFormIntegrator &integ, ParGridFunction &sol, ParFiniteElementSpace *flux_fes, ParFiniteElementSpace *smooth_flux_fes)
Construct a new L2ZienkiewiczZhuEstimator object.
void ComputeEstimates()
Compute the element error estimates.
void Reset() override
Reset the error estimator.
void SetLocalErrorNormP(int p)
Set the exponent, p, of the Lp norm used for computing the local element errors. Default value is 1.
virtual ~L2ZienkiewiczZhuEstimator()
Destroy a L2ZienkiewiczZhuEstimator object. Destroys, if owned, the FiniteElementSpace,...
ParFiniteElementSpace * flux_space
Ownership based on the flag own_flux_fes. Its Update() method is called automatically by this class w...
const Vector & GetLocalErrors() override
Get a Vector with all element errors.
real_t GetTotalError() const override
Return the total error from the last error estimate.
ParFiniteElementSpace * smooth_flux_space
Ownership based on the flag own_flux_fes. Its Update() method is called automatically by this class w...
L2ZienkiewiczZhuEstimator(BilinearFormIntegrator &integ, ParGridFunction &sol, ParFiniteElementSpace &flux_fes, ParFiniteElementSpace &smooth_flux_fes)
Construct a new L2ZienkiewiczZhuEstimator object.
ParGridFunction & solution
bool own_flux_fes
Ownership flag for flux_space and smooth_flux_space.
bool MeshIsModified()
Check if the mesh of the solution was modified.
int local_norm_p
Local L_p norm to use, default is 1.
BilinearFormIntegrator & integ
The LSZienkiewiczZhuEstimator class implements the Zienkiewicz-Zhu error estimation procedure [1,...
const Vector & GetLocalErrors() override
Get a Vector with all element errors.
virtual ~LSZienkiewiczZhuEstimator()
LSZienkiewiczZhuEstimator(BilinearFormIntegrator &integ, GridFunction &sol)
Construct a new LSZienkiewiczZhuEstimator object.
BilinearFormIntegrator & integ
void Reset() override
Reset the error estimator.
void DisableReconstructionAcrossSubdomains()
Disable reconstructing the flux in patches spanning different subdomains.
bool MeshIsModified()
Check if the mesh of the solution was modified.
void ComputeEstimates()
Compute the element error estimates.
bool subdomain_reconstruction
void SetTichonovRegularization(real_t tcoeff=1.0e-8)
Solve a Tichonov-regularized least-squares problem for the reconstructed fluxes. This is especially h...
void SetWithCoeff(bool w_coeff=true)
Consider the coefficient in BilinearFormIntegrator to calculate the fluxes for the error estimator.
real_t GetTotalError() const override
Return the total error from the last error estimate.
The LpErrorEstimator class compares the solution to a known coefficient.
LpErrorEstimator(int p, Coefficient &coef, GridFunction &sol)
Construct a new LpErrorEstimator object for a scalar field.
LpErrorEstimator(int p, GridFunction &sol)
Construct a new LpErrorEstimator object for a scalar field.
void SetCoef(VectorCoefficient &A)
LpErrorEstimator(int p, VectorCoefficient &coef, GridFunction &sol)
Construct a new LpErrorEstimator object for a vector field.
void ComputeEstimates()
Compute the element error estimates.
void SetLocalErrorNormP(int p)
Set the exponent, p, of the Lp norm used for computing the local element errors.
bool MeshIsModified()
Check if the mesh of the solution was modified.
const Vector & GetLocalErrors() override
Get a Vector with all element errors.
void SetCoef(Coefficient &A)
virtual ~LpErrorEstimator()
Destructor.
VectorCoefficient * vcoef
void Reset() override
Reset the error estimator.
Abstract parallel finite element space.
Class for parallel grid function.
Base class for vector Coefficients that optionally depend on time and space.
The ZienkiewiczZhuEstimator class implements the Zienkiewicz-Zhu error estimation procedure.
virtual ~ZienkiewiczZhuEstimator()
Destroy a ZienkiewiczZhuEstimator object. Destroys, if owned, the FiniteElementSpace,...
BilinearFormIntegrator & integ
bool MeshIsModified()
Check if the mesh of the solution was modified.
ZienkiewiczZhuEstimator(BilinearFormIntegrator &integ, GridFunction &sol, FiniteElementSpace *flux_fes)
Construct a new ZienkiewiczZhuEstimator object.
const Array< int > & GetAnisotropicFlags() override
Get an Array<int> with anisotropic flags for all mesh elements. Return an empty array when anisotropi...
bool own_flux_fes
Ownership flag for flux_space.
void ComputeEstimates()
Compute the element error estimates.
void SetWithCoeff(bool w_coeff=true)
Consider the coefficient in BilinearFormIntegrator to calculate the fluxes for the error estimator.
ZienkiewiczZhuEstimator(BilinearFormIntegrator &integ, GridFunction &sol, FiniteElementSpace &flux_fes)
Construct a new ZienkiewiczZhuEstimator object.
FiniteElementSpace * flux_space
Ownership based on own_flux_fes. Its Update() method is called automatically by this class when neede...
void SetAnisotropic(bool aniso=true)
Enable/disable anisotropic estimates. To enable this option, the BilinearFormIntegrator must support ...
void Reset() override
Reset the error estimator.
void SetFluxAveraging(int fa)
Set the way the flux is averaged (smoothed) across elements.
const Vector & GetLocalErrors() override
Get a Vector with all element errors.
real_t GetTotalError() const override
Return the total error from the last error estimate.
real_t p(const Vector &x, real_t t)