libMesh
jump_error_estimator.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2024 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #ifndef LIBMESH_JUMP_ERROR_ESTIMATOR_H
21 #define LIBMESH_JUMP_ERROR_ESTIMATOR_H
22 
23 // Local Includes
24 #include "libmesh/dense_vector.h"
25 #include "libmesh/error_estimator.h"
26 #include "libmesh/fem_context.h"
27 
28 // C++ includes
29 #include <cstddef>
30 #include <string>
31 #include <vector>
32 #include <memory>
33 
34 namespace libMesh
35 {
36 
37 // Forward Declarations
38 class Point;
39 class Elem;
40 
49 {
50 public:
51 
56  : ErrorEstimator(),
57  scale_by_n_flux_faces(false),
60  fine_context(),
62  fine_error(0),
63  coarse_error(0) {}
64 
70  JumpErrorEstimator (const JumpErrorEstimator &) = delete;
72 
78  virtual ~JumpErrorEstimator() = default;
79 
86  virtual void estimate_error (const System & system,
87  ErrorVector & error_per_cell,
88  const NumericVector<Number> * solution_vector = nullptr,
89  bool estimate_parent_error = false) override;
90 
100 
114 
115 protected:
120  void reinit_sides();
121 
126 
131  virtual void init_context(FEMContext & c);
132 
137  virtual void internal_side_integration() = 0;
138 
146  virtual bool boundary_side_integration() { return false; }
147 
153 
158  std::unique_ptr<FEMContext> fine_context, coarse_context;
159 
164 
168  unsigned int var;
169 };
170 
171 
172 } // namespace libMesh
173 
174 #endif // LIBMESH_JUMP_ERROR_ESTIMATOR_H
JumpErrorEstimator & operator=(const JumpErrorEstimator &)=delete
std::unique_ptr< FEMContext > fine_context
Context objects for integrating on the fine and coarse elements sharing a face.
The ErrorVector is a specialization of the StatisticsVector for error data computed on a finite eleme...
Definition: error_vector.h:50
This abstract base class implements utility functions for error estimators which are based on integra...
virtual void internal_side_integration()=0
The function, to be implemented by derived classes, which calculates an error term on an internal sid...
The libMesh namespace provides an interface to certain functionality in the library.
void reinit_sides()
A utility function to reinit the finite element data on elements sharing a side.
virtual bool boundary_side_integration()
The function, to be implemented by derived classes, which calculates an error term on a boundary side...
Manages consistently variables, degrees of freedom, and coefficient vectors.
Definition: system.h:96
This class holds functions that will estimate the error in a finite element solution on a given mesh...
This class provides all data required for a physics package (e.g.
Definition: fem_context.h:62
std::unique_ptr< FEMContext > coarse_context
Real fine_error
The fine and coarse error values to be set by each side_integration();.
virtual ~JumpErrorEstimator()=default
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
float coarse_n_flux_faces_increment()
A utility function to correctly increase n_flux_faces for the coarse element.
bool use_unweighted_quadrature_rules
This boolean flag allows you to use "unweighted" quadrature rules (sized to exactly integrate unweigh...
unsigned int var
The variable number currently being evaluated.
bool scale_by_n_flux_faces
This boolean flag allows you to scale the error indicator result for each element by the number of "f...
virtual void estimate_error(const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=nullptr, bool estimate_parent_error=false) override
This function uses the derived class&#39;s jump error estimate formula to estimate the error on each cell...
bool integrate_boundary_sides
A boolean flag, by default false, to be set to true if integrations with boundary_side_integration() ...
virtual void init_context(FEMContext &c)
An initialization function, to give derived classes a chance to request specific data from the FE obj...