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