libMesh
Public Types | Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | List of all members
libMesh::KellyErrorEstimator Class Reference


This class implements the Kelly error indicator which is based on the flux jumps between elements. More...

#include <kelly_error_estimator.h>

Inheritance diagram for libMesh::KellyErrorEstimator:
[legend]

Public Types

typedef std::map< std::pair< const System *, unsigned int >, std::unique_ptr< ErrorVector > > ErrorMap
 When calculating many error vectors at once, we need a data structure to hold them all. More...
 

Public Member Functions

 KellyErrorEstimator ()
 Constructor. More...
 
 KellyErrorEstimator (const KellyErrorEstimator &)=delete
 This class cannot be (default) copy constructed/assigned because its base class has unique_ptr members. More...
 
KellyErrorEstimatoroperator= (const KellyErrorEstimator &)=delete
 
 KellyErrorEstimator (KellyErrorEstimator &&)=default
 Defaulted move ctor, move assignment operator, and destructor. More...
 
KellyErrorEstimatoroperator= (KellyErrorEstimator &&)=default
 
virtual ~KellyErrorEstimator ()=default
 
void attach_flux_bc_function (std::pair< bool, Real > fptr(const System &system, const Point &p, const std::string &var_name))
 Register a user function to use in computing the flux BCs. More...
 
virtual ErrorEstimatorType type () const override
 
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's jump error estimate formula to estimate the error on each cell. More...
 
virtual void estimate_errors (const EquationSystems &equation_systems, ErrorVector &error_per_cell, const std::map< const System *, SystemNorm > &error_norms, const std::map< const System *, const NumericVector< Number > *> *solution_vectors=nullptr, bool estimate_parent_error=false)
 This virtual function can be redefined in derived classes, but by default computes the sum of the error_per_cell for each system in the equation_systems. More...
 
virtual void estimate_errors (const EquationSystems &equation_systems, ErrorMap &errors_per_cell, const std::map< const System *, const NumericVector< Number > *> *solution_vectors=nullptr, bool estimate_parent_error=false)
 This virtual function can be redefined in derived classes, but by default it calls estimate_error repeatedly to calculate the requested error vectors. More...
 

Public Attributes

bool scale_by_n_flux_faces
 This boolean flag allows you to scale the error indicator result for each element by the number of "flux faces" the element actually has. More...
 
bool use_unweighted_quadrature_rules
 This boolean flag allows you to use "unweighted" quadrature rules (sized to exactly integrate unweighted shape functions in master element space) rather than "default" quadrature rules (sized to exactly integrate polynomials of one higher degree than mass matrix terms). More...
 
SystemNorm error_norm
 When estimating the error in a single system, the error_norm is used to control the scaling and norm choice for each variable. More...
 

Protected Member Functions

virtual void init_context (FEMContext &c) override
 An initialization function, for requesting specific data from the FE objects. More...
 
virtual void internal_side_integration () override
 The function which calculates a normal derivative jump based error term on an internal side. More...
 
virtual bool boundary_side_integration () override
 The function which calculates a normal derivative jump based error term on a boundary side. More...
 
void reinit_sides ()
 A utility function to reinit the finite element data on elements sharing a side. More...
 
float coarse_n_flux_faces_increment ()
 A utility function to correctly increase n_flux_faces for the coarse element. More...
 
void reduce_error (std::vector< ErrorVectorReal > &error_per_cell, const Parallel::Communicator &comm) const
 This method takes the local error contributions in error_per_cell from each processor and combines them to get the global error vector. More...
 

Protected Attributes

std::pair< bool, Real >(* _bc_function )(const System &system, const Point &p, const std::string &var_name)
 Pointer to function that provides BC information. More...
 
bool integrate_boundary_sides
 A boolean flag, by default false, to be set to true if integrations with boundary_side_integration() should be performed. More...
 
std::unique_ptr< FEMContextfine_context
 Context objects for integrating on the fine and coarse elements sharing a face. More...
 
std::unique_ptr< FEMContextcoarse_context
 
Real fine_error
 The fine and coarse error values to be set by each side_integration();. More...
 
Real coarse_error
 
unsigned int var
 The variable number currently being evaluated. More...
 

Detailed Description


This class implements the Kelly error indicator which is based on the flux jumps between elements.

See the JumpErrorEstimator class for most user APIs

Full BibteX reference:

* @Article{Kelly83error,
* author = {D.~W.~Kelly and J.~P.~Gago and O.~C.~Zienkiewicz and I.~Babuska},
* title  = {{A posteriori error analysis and adaptive
*            processes in the finite element method: Part I Error analysis}},
* journal = {Int. J. Num. Meth. Engng.},
* volume  = {19},
* pages   = {1593--1619},
* year    = {1983}
* }
* 
Author
Benjamin S. Kirk
Date
2003

Definition at line 59 of file kelly_error_estimator.h.

Member Typedef Documentation

◆ ErrorMap

typedef std::map<std::pair<const System *, unsigned int>, std::unique_ptr<ErrorVector> > libMesh::ErrorEstimator::ErrorMap
inherited

When calculating many error vectors at once, we need a data structure to hold them all.

Definition at line 121 of file error_estimator.h.

Constructor & Destructor Documentation

◆ KellyErrorEstimator() [1/3]

libMesh::KellyErrorEstimator::KellyErrorEstimator ( )

Constructor.

Responsible for initializing the _bc_function function pointer to nullptr. Defaults to H1 seminorm; changes to system norm are ignored.

Definition at line 41 of file kelly_error_estimator.C.

References libMesh::ErrorEstimator::error_norm, and libMesh::H1_SEMINORM.

41  :
43  _bc_function(nullptr)
44 {
46 }
std::pair< bool, Real >(* _bc_function)(const System &system, const Point &p, const std::string &var_name)
Pointer to function that provides BC information.
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...

◆ KellyErrorEstimator() [2/3]

libMesh::KellyErrorEstimator::KellyErrorEstimator ( const KellyErrorEstimator )
delete

This class cannot be (default) copy constructed/assigned because its base class has unique_ptr members.

◆ KellyErrorEstimator() [3/3]

libMesh::KellyErrorEstimator::KellyErrorEstimator ( KellyErrorEstimator &&  )
default

Defaulted move ctor, move assignment operator, and destructor.

◆ ~KellyErrorEstimator()

virtual libMesh::KellyErrorEstimator::~KellyErrorEstimator ( )
virtualdefault

Member Function Documentation

◆ attach_flux_bc_function()

void libMesh::KellyErrorEstimator::attach_flux_bc_function ( std::pair< bool, Real >   fptrconst System &system, const Point &p, const std::string &var_name)

Register a user function to use in computing the flux BCs.

Definition at line 207 of file kelly_error_estimator.C.

References _bc_function, fptr(), and libMesh::JumpErrorEstimator::integrate_boundary_sides.

210 {
211  _bc_function = fptr;
212 
213  // We may be turning boundary side integration on or off
214  if (fptr)
216  else
217  integrate_boundary_sides = false;
218 }
std::pair< bool, Real >(* _bc_function)(const System &system, const Point &p, const std::string &var_name)
Pointer to function that provides BC information.
Number fptr(const Point &p, const Parameters &, const std::string &libmesh_dbg_var(sys_name), const std::string &unknown_name)
Definition: projection.C:80
bool integrate_boundary_sides
A boolean flag, by default false, to be set to true if integrations with boundary_side_integration() ...

◆ boundary_side_integration()

bool libMesh::KellyErrorEstimator::boundary_side_integration ( )
overrideprotectedvirtual

The function which calculates a normal derivative jump based error term on a boundary side.

Returns
true if the flux bc function is in fact defined on the current side.

Reimplemented from libMesh::JumpErrorEstimator.

Definition at line 135 of file kelly_error_estimator.C.

References _bc_function, libMesh::Elem::dim(), libMesh::ErrorEstimator::error_norm, libMesh::JumpErrorEstimator::fine_context, libMesh::JumpErrorEstimator::fine_error, libMesh::FEGenericBase< OutputType >::get_dphi(), libMesh::FEAbstract::get_JxW(), libMesh::FEAbstract::get_normals(), libMesh::FEAbstract::get_xyz(), libMesh::Elem::hmax(), libMesh::FEAbstract::n_quadrature_points(), libMesh::TensorTools::norm_sq(), libMesh::Real, libMesh::JumpErrorEstimator::var, and libMesh::SystemNorm::weight().

136 {
137  const Elem & fine_elem = fine_context->get_elem();
138 
139  FEBase * fe_fine = nullptr;
140  fine_context->get_side_fe( var, fe_fine, fine_elem.dim() );
141 
142  const std::string & var_name =
143  fine_context->get_system().variable_name(var);
144 
145  std::vector<std::vector<RealGradient>> dphi_fine = fe_fine->get_dphi();
146  std::vector<Point> face_normals = fe_fine->get_normals();
147  std::vector<Real> JxW_face = fe_fine->get_JxW();
148  std::vector<Point> qface_point = fe_fine->get_xyz();
149 
150  // The reinitialization also recomputes the locations of
151  // the quadrature points on the side. By checking if the
152  // first quadrature point on the side is on a flux boundary
153  // for a particular variable, we will determine if the whole
154  // element is on a flux boundary (assuming quadrature points
155  // are strictly contained in the side).
156  if (this->_bc_function(fine_context->get_system(),
157  qface_point[0], var_name).first)
158  {
159  const Real h = fine_elem.hmax();
160 
161  // The number of quadrature points
162  const unsigned int n_qp = fe_fine->n_quadrature_points();
163 
164  // The error contribution from this face
165  Real error = 1.e-30;
166 
167  // loop over the integration points on the face.
168  for (unsigned int qp=0; qp<n_qp; qp++)
169  {
170  // Value of the imposed flux BC at this quadrature point.
171  const std::pair<bool,Real> flux_bc =
172  this->_bc_function(fine_context->get_system(),
173  qface_point[qp], var_name);
174 
175  // Be sure the BC function still thinks we're on the
176  // flux boundary.
177  libmesh_assert_equal_to (flux_bc.first, true);
178 
179  // The solution gradient from each element
180  Gradient grad_fine = fine_context->side_gradient(var, qp);
181 
182  // The difference between the desired BC and the approximate solution.
183  const Number jump = flux_bc.second - grad_fine*face_normals[qp];
184 
185  // The flux jump squared. If using complex numbers,
186  // TensorTools::norm_sq(z) returns |z|^2, where |z| is the modulus of z.
187  const Real jump2 = TensorTools::norm_sq(jump);
188 
189  // Integrate the error on the face. The error is
190  // scaled by an additional power of h, where h is
191  // the maximum side length for the element. This
192  // arises in the definition of the indicator.
193  error += JxW_face[qp]*jump2;
194 
195  } // End quadrature point loop
196 
197  fine_error = error*h*error_norm.weight(var);
198 
199  return true;
200  } // end if side on flux boundary
201  return false;
202 }
std::pair< bool, Real >(* _bc_function)(const System &system, const Point &p, const std::string &var_name)
Pointer to function that provides BC information.
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
std::unique_ptr< FEMContext > fine_context
Context objects for integrating on the fine and coarse elements sharing a face.
NumberVectorValue Gradient
FEGenericBase< Real > FEBase
Real fine_error
The fine and coarse error values to be set by each side_integration();.
Real weight(unsigned int var) const
Definition: system_norm.C:134
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
auto norm_sq(const T &a) -> decltype(std::norm(a))
Definition: tensor_tools.h:104
unsigned int var
The variable number currently being evaluated.

◆ coarse_n_flux_faces_increment()

float libMesh::JumpErrorEstimator::coarse_n_flux_faces_increment ( )
protectedinherited

A utility function to correctly increase n_flux_faces for the coarse element.

Definition at line 488 of file jump_error_estimator.C.

References libMesh::JumpErrorEstimator::coarse_context, dim, and libMesh::JumpErrorEstimator::fine_context.

Referenced by libMesh::JumpErrorEstimator::estimate_error().

489 {
490  // Keep track of the number of internal flux sides found on each
491  // element
492  unsigned short dim = coarse_context->get_elem().dim();
493 
494  const unsigned int divisor =
495  1 << (dim-1)*(fine_context->get_elem().level() -
496  coarse_context->get_elem().level());
497 
498  // With a difference of n levels between fine and coarse elements,
499  // we compute a fractional flux face for the coarse element by adding:
500  // 1/2^n in 2D
501  // 1/4^n in 3D
502  // each time. This code will get hit 2^n times in 2D and 4^n
503  // times in 3D so that the final flux face count for the coarse
504  // element will be an integer value.
505 
506  return 1.0f / static_cast<float>(divisor);
507 }
std::unique_ptr< FEMContext > fine_context
Context objects for integrating on the fine and coarse elements sharing a face.
unsigned int dim
std::unique_ptr< FEMContext > coarse_context

◆ estimate_error()

void libMesh::JumpErrorEstimator::estimate_error ( const System system,
ErrorVector error_per_cell,
const NumericVector< Number > *  solution_vector = nullptr,
bool  estimate_parent_error = false 
)
overridevirtualinherited

This function uses the derived class's jump error estimate formula to estimate the error on each cell.

The estimated error is output in the vector error_per_cell

Conventions for assigning the direction of the normal:

  • e & f are global element ids

Case (1.) Elements are at the same level, e<f Compute the flux jump on the face and add it as a contribution to error_per_cell[e] and error_per_cell[f]


| | |

f
e —> n

Case (2.) The neighbor is at a higher level. Compute the flux jump on e's face and add it as a contribution to error_per_cell[e] and error_per_cell[f]


| | | | | | e |—> n | | | | | |--------—| f |


Implements libMesh::ErrorEstimator.

Definition at line 56 of file jump_error_estimator.C.

References libMesh::Elem::active(), libMesh::ElemInternal::active_family_tree_by_neighbor(), libMesh::JumpErrorEstimator::boundary_side_integration(), libMesh::Elem::child_ref_range(), libMesh::JumpErrorEstimator::coarse_context, libMesh::JumpErrorEstimator::coarse_error, libMesh::JumpErrorEstimator::coarse_n_flux_faces_increment(), libMesh::FEGenericBase< OutputType >::coarsened_dof_values(), dim, libMesh::ErrorEstimator::error_norm, libMesh::ErrorVectorReal, libMesh::JumpErrorEstimator::fine_context, libMesh::JumpErrorEstimator::fine_error, libMesh::NumericVector< T >::get(), libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMesh::FEAbstract::get_xyz(), libMesh::DofObject::id(), libMesh::index_range(), libMesh::JumpErrorEstimator::init_context(), libMesh::JumpErrorEstimator::integrate_boundary_sides, libMesh::JumpErrorEstimator::internal_side_integration(), libMesh::Elem::level(), libMesh::libmesh_ignore(), mesh, n_vars, libMesh::System::n_vars(), libMesh::Elem::neighbor_ptr(), libMesh::Elem::parent(), libMesh::ErrorEstimator::reduce_error(), libMesh::JumpErrorEstimator::reinit_sides(), libMesh::SCALAR, libMesh::JumpErrorEstimator::scale_by_n_flux_faces, libMesh::Elem::side_index_range(), libMesh::DenseVector< T >::size(), libMesh::System::solution, std::sqrt(), libMesh::NumericVector< T >::swap(), libMesh::JumpErrorEstimator::use_unweighted_quadrature_rules, libMesh::JumpErrorEstimator::var, and libMesh::SystemNorm::weight().

Referenced by assemble_and_solve(), and main().

60 {
61  LOG_SCOPE("estimate_error()", "JumpErrorEstimator");
62 
99  // This parameter is not used when !LIBMESH_ENABLE_AMR.
100  libmesh_ignore(estimate_parent_error);
101 
102  // The current mesh
103  const MeshBase & mesh = system.get_mesh();
104 
105  // The number of variables in the system
106  const unsigned int n_vars = system.n_vars();
107 
108  // The DofMap for this system
109 #ifdef LIBMESH_ENABLE_AMR
110  const DofMap & dof_map = system.get_dof_map();
111 #endif
112 
113  // Resize the error_per_cell vector to be
114  // the number of elements, initialize it to 0.
115  error_per_cell.resize (mesh.max_elem_id());
116  std::fill (error_per_cell.begin(), error_per_cell.end(), 0.);
117 
118  // Declare a vector of floats which is as long as
119  // error_per_cell above, and fill with zeros. This vector will be
120  // used to keep track of the number of edges (faces) on each active
121  // element which are either:
122  // 1) an internal edge
123  // 2) an edge on a Neumann boundary for which a boundary condition
124  // function has been specified.
125  // The error estimator can be scaled by the number of flux edges (faces)
126  // which the element actually has to obtain a more uniform measure
127  // of the error. Use floats instead of ints since in case 2 (above)
128  // f gets 1/2 of a flux face contribution from each of his
129  // neighbors
130  std::vector<float> n_flux_faces;
132  n_flux_faces.resize(error_per_cell.size(), 0);
133 
134  // Prepare current_local_solution to localize a non-standard
135  // solution vector if necessary
136  if (solution_vector && solution_vector != system.solution.get())
137  {
138  NumericVector<Number> * newsol =
139  const_cast<NumericVector<Number> *>(solution_vector);
140  System & sys = const_cast<System &>(system);
141  newsol->swap(*sys.solution);
142  sys.update();
143  }
144 
145  // We don't use full elem_jacobian or subjacobians here.
146  fine_context = std::make_unique<FEMContext>
147  (system, nullptr, /* allocate_local_matrices = */ false);
148  coarse_context = std::make_unique<FEMContext>
149  (system, nullptr, /* allocate_local_matrices = */ false);
150 
151  // Don't overintegrate - we're evaluating differences of FE values,
152  // not products of them.
154  fine_context->use_unweighted_quadrature_rules(system.extra_quadrature_order);
155 
156  // Loop over all the variables we've been requested to find jumps in, to
157  // pre-request
158  for (var=0; var<n_vars; var++)
159  {
160  // Skip variables which aren't part of our norm,
161  // as well as SCALAR variables, which have no jumps
162  if (error_norm.weight(var) == 0.0 ||
163  system.variable_type(var).family == SCALAR)
164  continue;
165 
166  // FIXME: Need to generalize this to vector-valued elements. [PB]
167  FEBase * side_fe = nullptr;
168 
169  const std::set<unsigned char> & elem_dims =
170  fine_context->elem_dimensions();
171 
172  for (const auto & dim : elem_dims)
173  {
174  fine_context->get_side_fe( var, side_fe, dim );
175 
176  side_fe->get_xyz();
177  }
178  }
179 
180  this->init_context(*fine_context);
181  this->init_context(*coarse_context);
182 
183  // Iterate over all the active elements in the mesh
184  // that live on this processor.
185  for (const auto & e : mesh.active_local_element_ptr_range())
186  {
187  const dof_id_type e_id = e->id();
188 
189 #ifdef LIBMESH_ENABLE_AMR
190 
191  if (e->infinite())
192  {
193  libmesh_warning("Warning: Jumps on the border of infinite elements are ignored."
194  << std::endl);
195  continue;
196  }
197 
198  // See if the parent of element e has been examined yet;
199  // if not, we may want to compute the estimator on it
200  const Elem * parent = e->parent();
201 
202  // We only can compute and only need to compute on
203  // parents with all active children
204  bool compute_on_parent = true;
205  if (!parent || !estimate_parent_error)
206  compute_on_parent = false;
207  else
208  for (auto & child : parent->child_ref_range())
209  if (!child.active())
210  compute_on_parent = false;
211 
212  if (compute_on_parent &&
213  !error_per_cell[parent->id()])
214  {
215  // Compute a projection onto the parent
216  DenseVector<Number> Uparent;
218  (*(system.solution), dof_map, parent, Uparent, false);
219 
220  // Loop over the neighbors of the parent
221  for (auto n_p : parent->side_index_range())
222  {
223  if (parent->neighbor_ptr(n_p) != nullptr) // parent has a neighbor here
224  {
225  // Find the active neighbors in this direction
226  std::vector<const Elem *> active_neighbors;
227  parent->neighbor_ptr(n_p)->
228  active_family_tree_by_neighbor(active_neighbors,
229  parent);
230  // Compute the flux to each active neighbor
231  for (std::size_t a=0,
232  n_active_neighbors = active_neighbors.size();
233  a != n_active_neighbors; ++a)
234  {
235  const Elem * f = active_neighbors[a];
236 
237  if (f ->infinite()) // don't take infinite elements into account
238  continue;
239 
240  // FIXME - what about when f->level <
241  // parent->level()??
242  if (f->level() >= parent->level())
243  {
244  fine_context->pre_fe_reinit(system, f);
245  coarse_context->pre_fe_reinit(system, parent);
246  libmesh_assert_equal_to
247  (coarse_context->get_elem_solution().size(),
248  Uparent.size());
249  coarse_context->get_elem_solution() = Uparent;
250 
251  this->reinit_sides();
252 
253  // Loop over all significant variables in the system
254  for (var=0; var<n_vars; var++)
255  if (error_norm.weight(var) != 0.0 &&
256  system.variable_type(var).family != SCALAR)
257  {
259 
260  error_per_cell[fine_context->get_elem().id()] +=
261  static_cast<ErrorVectorReal>(fine_error);
262  error_per_cell[coarse_context->get_elem().id()] +=
263  static_cast<ErrorVectorReal>(coarse_error);
264  }
265 
266  // Keep track of the number of internal flux
267  // sides found on each element
269  {
270  n_flux_faces[fine_context->get_elem().id()]++;
271  n_flux_faces[coarse_context->get_elem().id()] +=
273  }
274  }
275  }
276  }
277  else if (integrate_boundary_sides)
278  {
279  fine_context->pre_fe_reinit(system, parent);
280  libmesh_assert_equal_to
281  (fine_context->get_elem_solution().size(),
282  Uparent.size());
283  fine_context->get_elem_solution() = Uparent;
284  fine_context->side = cast_int<unsigned char>(n_p);
285  fine_context->side_fe_reinit();
286 
287  // If we find a boundary flux for any variable,
288  // let's just count it as a flux face for all
289  // variables. Otherwise we'd need to keep track of
290  // a separate n_flux_faces and error_per_cell for
291  // every single var.
292  bool found_boundary_flux = false;
293 
294  for (var=0; var<n_vars; var++)
295  if (error_norm.weight(var) != 0.0 &&
296  system.variable_type(var).family != SCALAR)
297  {
298  if (this->boundary_side_integration())
299  {
300  error_per_cell[fine_context->get_elem().id()] +=
301  static_cast<ErrorVectorReal>(fine_error);
302  found_boundary_flux = true;
303  }
304  }
305 
306  if (scale_by_n_flux_faces && found_boundary_flux)
307  n_flux_faces[fine_context->get_elem().id()]++;
308  }
309  }
310  }
311 #endif // #ifdef LIBMESH_ENABLE_AMR
312 
313  // If we do any more flux integration, e will be the fine element
314  fine_context->pre_fe_reinit(system, e);
315 
316  // Loop over the neighbors of element e
317  for (auto n_e : e->side_index_range())
318  {
319  if ((e->neighbor_ptr(n_e) != nullptr) ||
321  {
322  fine_context->side = cast_int<unsigned char>(n_e);
323  fine_context->side_fe_reinit();
324  }
325 
326  // e is not on the boundary (infinite elements are treated as boundary)
327  if (e->neighbor_ptr(n_e) != nullptr
328  && !e->neighbor_ptr(n_e) ->infinite())
329  {
330 
331  const Elem * f = e->neighbor_ptr(n_e);
332  const dof_id_type f_id = f->id();
333 
334  // Compute flux jumps if we are in case 1 or case 2.
335  if ((f->active() && (f->level() == e->level()) && (e_id < f_id))
336  || (f->level() < e->level()))
337  {
338  // f is now the coarse element
339  coarse_context->pre_fe_reinit(system, f);
340 
341  this->reinit_sides();
342 
343  // Loop over all significant variables in the system
344  for (var=0; var<n_vars; var++)
345  if (error_norm.weight(var) != 0.0 &&
346  system.variable_type(var).family != SCALAR)
347  {
349 
350  error_per_cell[fine_context->get_elem().id()] +=
351  static_cast<ErrorVectorReal>(fine_error);
352  error_per_cell[coarse_context->get_elem().id()] +=
353  static_cast<ErrorVectorReal>(coarse_error);
354  }
355 
356  // Keep track of the number of internal flux
357  // sides found on each element
359  {
360  n_flux_faces[fine_context->get_elem().id()]++;
361  n_flux_faces[coarse_context->get_elem().id()] +=
363  }
364  } // end if (case1 || case2)
365  } // if (e->neighbor(n_e) != nullptr)
366 
367  // Otherwise, e is on the boundary. If it happens to
368  // be on a Dirichlet boundary, we need not do anything.
369  // On the other hand, if e is on a Neumann (flux) boundary
370  // with grad(u).n = g, we need to compute the additional residual
371  // (h * \int |g - grad(u_h).n|^2 dS)^(1/2).
372  // We can only do this with some knowledge of the boundary
373  // conditions, i.e. the user must have attached an appropriate
374  // BC function.
375  else if (integrate_boundary_sides)
376  {
377  bool found_boundary_flux = false;
378 
379  for (var=0; var<n_vars; var++)
380  if (error_norm.weight(var) != 0.0 &&
381  system.variable_type(var).family != SCALAR)
382  if (this->boundary_side_integration())
383  {
384  error_per_cell[fine_context->get_elem().id()] +=
385  static_cast<ErrorVectorReal>(fine_error);
386  found_boundary_flux = true;
387  }
388 
389  if (scale_by_n_flux_faces && found_boundary_flux)
390  n_flux_faces[fine_context->get_elem().id()]++;
391  } // end if (e->neighbor_ptr(n_e) == nullptr)
392  } // end loop over neighbors
393  } // End loop over active local elements
394 
395 
396  // Each processor has now computed the error contributions
397  // for its local elements. We need to sum the vector
398  // and then take the square-root of each component. Note
399  // that we only need to sum if we are running on multiple
400  // processors, and we only need to take the square-root
401  // if the value is nonzero. There will in general be many
402  // zeros for the inactive elements.
403 
404  // First sum the vector of estimated error values
405  this->reduce_error(error_per_cell, system.comm());
406 
407  // Compute the square-root of each component.
408  for (auto i : index_range(error_per_cell))
409  if (error_per_cell[i] != 0.)
410  error_per_cell[i] = std::sqrt(error_per_cell[i]);
411 
412 
413  if (this->scale_by_n_flux_faces)
414  {
415  // Sum the vector of flux face counts
416  this->reduce_error(n_flux_faces, system.comm());
417 
418  // Sanity check: Make sure the number of flux faces is
419  // always an integer value
420 #ifdef DEBUG
421  for (const auto & val : n_flux_faces)
422  libmesh_assert_equal_to (val, static_cast<float>(static_cast<unsigned int>(val)));
423 #endif
424 
425  // Scale the error by the number of flux faces for each element
426  for (auto i : index_range(n_flux_faces))
427  {
428  if (n_flux_faces[i] == 0.0) // inactive or non-local element
429  continue;
430 
431  error_per_cell[i] /= static_cast<ErrorVectorReal>(n_flux_faces[i]);
432  }
433  }
434 
435  // If we used a non-standard solution before, now is the time to fix
436  // the current_local_solution
437  if (solution_vector && solution_vector != system.solution.get())
438  {
439  NumericVector<Number> * newsol =
440  const_cast<NumericVector<Number> *>(solution_vector);
441  System & sys = const_cast<System &>(system);
442  newsol->swap(*sys.solution);
443  sys.update();
444  }
445 }
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
void active_family_tree_by_neighbor(T elem, std::vector< T > &family, T neighbor_in, bool reset=true)
std::unique_ptr< FEMContext > fine_context
Context objects for integrating on the fine and coarse elements sharing a face.
unsigned int dim
MeshBase & mesh
virtual void internal_side_integration()=0
The function, to be implemented by derived classes, which calculates an error term on an internal sid...
ADRealEigenVector< T, D, asd > sqrt(const ADRealEigenVector< T, D, asd > &)
Definition: type_vector.h:53
DIE A HORRIBLE DEATH HERE typedef float ErrorVectorReal
static void coarsened_dof_values(const NumericVector< Number > &global_vector, const DofMap &dof_map, const Elem *coarse_elem, DenseVector< Number > &coarse_dofs, const unsigned int var, const bool use_old_dof_indices=false)
Creates a local projection on coarse_elem, based on the DoF values in global_vector for it&#39;s children...
Definition: fe_base.C:1012
void reinit_sides()
A utility function to reinit the finite element data on elements sharing a side.
void libmesh_ignore(const Args &...)
virtual bool boundary_side_integration()
The function, to be implemented by derived classes, which calculates an error term on a boundary side...
unsigned int n_vars
void reduce_error(std::vector< ErrorVectorReal > &error_per_cell, const Parallel::Communicator &comm) const
This method takes the local error contributions in error_per_cell from each processor and combines th...
std::unique_ptr< FEMContext > coarse_context
FEGenericBase< Real > FEBase
Real fine_error
The fine and coarse error values to be set by each side_integration();.
Real weight(unsigned int var) const
Definition: system_norm.C:134
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...
bool integrate_boundary_sides
A boolean flag, by default false, to be set to true if integrations with boundary_side_integration() ...
template class LIBMESH_EXPORT NumericVector< Number >
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111
virtual void init_context(FEMContext &c)
An initialization function, to give derived classes a chance to request specific data from the FE obj...
uint8_t dof_id_type
Definition: id_types.h:67

◆ estimate_errors() [1/2]

void libMesh::ErrorEstimator::estimate_errors ( const EquationSystems equation_systems,
ErrorVector error_per_cell,
const std::map< const System *, SystemNorm > &  error_norms,
const std::map< const System *, const NumericVector< Number > *> *  solution_vectors = nullptr,
bool  estimate_parent_error = false 
)
virtualinherited

This virtual function can be redefined in derived classes, but by default computes the sum of the error_per_cell for each system in the equation_systems.

Currently this function ignores the error_norm member variable, and uses the function argument error_norms instead.

This function is named estimate_errors instead of estimate_error because otherwise C++ can get confused.

Reimplemented in libMesh::UniformRefinementEstimator.

Definition at line 47 of file error_estimator.C.

References libMesh::ErrorEstimator::error_norm, libMesh::ErrorEstimator::estimate_error(), libMesh::EquationSystems::get_system(), libMesh::index_range(), libMesh::make_range(), and libMesh::EquationSystems::n_systems().

52 {
53  SystemNorm old_error_norm = this->error_norm;
54 
55  // Sum the error values from each system
56  for (auto s : make_range(equation_systems.n_systems()))
57  {
58  ErrorVector system_error_per_cell;
59  const System & sys = equation_systems.get_system(s);
60  if (error_norms.find(&sys) == error_norms.end())
61  this->error_norm = old_error_norm;
62  else
63  this->error_norm = error_norms.find(&sys)->second;
64 
65  const NumericVector<Number> * solution_vector = nullptr;
66  if (solution_vectors &&
67  solution_vectors->find(&sys) != solution_vectors->end())
68  solution_vector = solution_vectors->find(&sys)->second;
69 
70  this->estimate_error(sys, system_error_per_cell,
71  solution_vector, estimate_parent_error);
72 
73  if (s)
74  {
75  libmesh_assert_equal_to (error_per_cell.size(), system_error_per_cell.size());
76  for (auto i : index_range(error_per_cell))
77  error_per_cell[i] += system_error_per_cell[i];
78  }
79  else
80  error_per_cell = system_error_per_cell;
81  }
82 
83  // Restore our old state before returning
84  this->error_norm = old_error_norm;
85 }
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
virtual void estimate_error(const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=nullptr, bool estimate_parent_error=false)=0
This pure virtual function must be redefined in derived classes to compute the error for each active ...
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
template class LIBMESH_EXPORT NumericVector< Number >
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111

◆ estimate_errors() [2/2]

void libMesh::ErrorEstimator::estimate_errors ( const EquationSystems equation_systems,
ErrorMap errors_per_cell,
const std::map< const System *, const NumericVector< Number > *> *  solution_vectors = nullptr,
bool  estimate_parent_error = false 
)
virtualinherited

This virtual function can be redefined in derived classes, but by default it calls estimate_error repeatedly to calculate the requested error vectors.

FIXME: This is a default implementation - derived classes should reimplement it for efficiency.

Currently this function ignores the error_norm.weight() values because it calculates each variable's error individually, unscaled.

The user selects which errors get computed by filling a map with error vectors: If errors_per_cell[&system][v] exists, it will be filled with the error values in variable v of system

Reimplemented in libMesh::UniformRefinementEstimator.

Definition at line 93 of file error_estimator.C.

References libMesh::ErrorEstimator::error_norm, libMesh::ErrorEstimator::estimate_error(), libMesh::EquationSystems::get_system(), libMesh::make_range(), libMesh::EquationSystems::n_systems(), n_vars, libMesh::System::n_vars(), and libMesh::SystemNorm::type().

97 {
98  SystemNorm old_error_norm = this->error_norm;
99 
100  // Find the requested error values from each system
101  for (auto s : make_range(equation_systems.n_systems()))
102  {
103  const System & sys = equation_systems.get_system(s);
104 
105  unsigned int n_vars = sys.n_vars();
106 
107  for (unsigned int v = 0; v != n_vars; ++v)
108  {
109  // Only fill in ErrorVectors the user asks for
110  if (!errors_per_cell.count(std::make_pair(&sys, v)))
111  continue;
112 
113  // Calculate error in only one variable
114  std::vector<Real> weights(n_vars, 0.0);
115  weights[v] = 1.0;
116  this->error_norm =
117  SystemNorm(std::vector<FEMNormType>(n_vars, old_error_norm.type(v)),
118  weights);
119 
120  const NumericVector<Number> * solution_vector = nullptr;
121  if (solution_vectors &&
122  solution_vectors->find(&sys) != solution_vectors->end())
123  solution_vector = solution_vectors->find(&sys)->second;
124 
125  this->estimate_error
126  (sys, *errors_per_cell[std::make_pair(&sys, v)],
127  solution_vector, estimate_parent_error);
128  }
129  }
130 
131  // Restore our old state before returning
132  this->error_norm = old_error_norm;
133 }
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
virtual void estimate_error(const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=nullptr, bool estimate_parent_error=false)=0
This pure virtual function must be redefined in derived classes to compute the error for each active ...
unsigned int n_vars
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
template class LIBMESH_EXPORT NumericVector< Number >

◆ init_context()

void libMesh::KellyErrorEstimator::init_context ( FEMContext c)
overrideprotectedvirtual

An initialization function, for requesting specific data from the FE objects.

Reimplemented from libMesh::JumpErrorEstimator.

Definition at line 59 of file kelly_error_estimator.C.

References libMesh::JumpErrorEstimator::coarse_context, dim, libMesh::FEMContext::elem_dimensions(), libMesh::ErrorEstimator::error_norm, libMesh::FEGenericBase< OutputType >::get_dphi(), libMesh::FEAbstract::get_normals(), libMesh::FEMContext::get_side_fe(), libMesh::DiffContext::n_vars(), n_vars, and libMesh::SystemNorm::weight().

60 {
61  const unsigned int n_vars = c.n_vars();
62  for (unsigned int v=0; v<n_vars; v++)
63  {
64  // Possibly skip this variable
65  if (error_norm.weight(v) == 0.0) continue;
66 
67  // FIXME: Need to generalize this to vector-valued elements. [PB]
68  FEBase * side_fe = nullptr;
69 
70  const std::set<unsigned char> & elem_dims =
71  c.elem_dimensions();
72 
73  for (const auto & dim : elem_dims)
74  {
75  c.get_side_fe( v, side_fe, dim );
76 
77  // We'll need gradients on both sides for flux jump computation
78  side_fe->get_dphi();
79 
80  // But we only need normal vectors from one side
81  if (&c != coarse_context.get())
82  side_fe->get_normals();
83  }
84  }
85 }
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
unsigned int dim
unsigned int n_vars
std::unique_ptr< FEMContext > coarse_context
FEGenericBase< Real > FEBase
Real weight(unsigned int var) const
Definition: system_norm.C:134

◆ internal_side_integration()

void libMesh::KellyErrorEstimator::internal_side_integration ( )
overrideprotectedvirtual

The function which calculates a normal derivative jump based error term on an internal side.

Implements libMesh::JumpErrorEstimator.

Definition at line 90 of file kelly_error_estimator.C.

References libMesh::JumpErrorEstimator::coarse_context, libMesh::JumpErrorEstimator::coarse_error, libMesh::Elem::dim(), libMesh::ErrorEstimator::error_norm, libMesh::JumpErrorEstimator::fine_context, libMesh::JumpErrorEstimator::fine_error, libMesh::FEGenericBase< OutputType >::get_dphi(), libMesh::FEAbstract::get_JxW(), libMesh::FEAbstract::get_normals(), libMesh::Elem::hmax(), libMesh::FEAbstract::n_quadrature_points(), libMesh::TensorTools::norm_sq(), libMesh::Real, libMesh::JumpErrorEstimator::var, and libMesh::SystemNorm::weight().

91 {
92  const Elem & coarse_elem = coarse_context->get_elem();
93  const Elem & fine_elem = fine_context->get_elem();
94 
95  FEBase * fe_fine = nullptr;
96  fine_context->get_side_fe( var, fe_fine, fine_elem.dim() );
97 
98  FEBase * fe_coarse = nullptr;
99  coarse_context->get_side_fe( var, fe_coarse, fine_elem.dim() );
100 
101  Real error = 1.e-30;
102  unsigned int n_qp = fe_fine->n_quadrature_points();
103 
104  std::vector<std::vector<RealGradient>> dphi_coarse = fe_coarse->get_dphi();
105  std::vector<std::vector<RealGradient>> dphi_fine = fe_fine->get_dphi();
106  std::vector<Point> face_normals = fe_fine->get_normals();
107  std::vector<Real> JxW_face = fe_fine->get_JxW();
108 
109  for (unsigned int qp=0; qp != n_qp; ++qp)
110  {
111  // Calculate solution gradients on fine and coarse elements
112  // at this quadrature point
113  Gradient
114  grad_fine = fine_context->side_gradient(var, qp),
115  grad_coarse = coarse_context->side_gradient(var, qp);
116 
117  // Find the jump in the normal derivative
118  // at this quadrature point
119  const Number jump = (grad_fine - grad_coarse)*face_normals[qp];
120  const Real jump2 = TensorTools::norm_sq(jump);
121 
122  // Accumulate the jump integral
123  error += JxW_face[qp] * jump2;
124  }
125 
126  // Add the h-weighted jump integral to each error term
127  fine_error =
128  error * fine_elem.hmax() * error_norm.weight(var);
129  coarse_error =
130  error * coarse_elem.hmax() * error_norm.weight(var);
131 }
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
std::unique_ptr< FEMContext > fine_context
Context objects for integrating on the fine and coarse elements sharing a face.
NumberVectorValue Gradient
std::unique_ptr< FEMContext > coarse_context
FEGenericBase< Real > FEBase
Real fine_error
The fine and coarse error values to be set by each side_integration();.
Real weight(unsigned int var) const
Definition: system_norm.C:134
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
auto norm_sq(const T &a) -> decltype(std::norm(a))
Definition: tensor_tools.h:104
unsigned int var
The variable number currently being evaluated.

◆ operator=() [1/2]

KellyErrorEstimator& libMesh::KellyErrorEstimator::operator= ( const KellyErrorEstimator )
delete

◆ operator=() [2/2]

KellyErrorEstimator& libMesh::KellyErrorEstimator::operator= ( KellyErrorEstimator &&  )
default

◆ reduce_error()

void libMesh::ErrorEstimator::reduce_error ( std::vector< ErrorVectorReal > &  error_per_cell,
const Parallel::Communicator comm 
) const
protectedinherited

This method takes the local error contributions in error_per_cell from each processor and combines them to get the global error vector.

Definition at line 32 of file error_estimator.C.

References TIMPI::Communicator::sum().

Referenced by libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::PatchRecoveryErrorEstimator::estimate_error(), libMesh::WeightedPatchRecoveryErrorEstimator::estimate_error(), libMesh::JumpErrorEstimator::estimate_error(), libMesh::AdjointRefinementEstimator::estimate_error(), and libMesh::ExactErrorEstimator::estimate_error().

34 {
35  // This function must be run on all processors at once
36  // parallel_object_only();
37 
38  // Each processor has now computed the error contributions
39  // for its local elements. We may need to sum the vector to
40  // recover the error for each element.
41 
42  comm.sum(error_per_cell);
43 }

◆ reinit_sides()

void libMesh::JumpErrorEstimator::reinit_sides ( )
protectedinherited

A utility function to reinit the finite element data on elements sharing a side.

Definition at line 450 of file jump_error_estimator.C.

References libMesh::JumpErrorEstimator::coarse_context, dim, libMesh::ErrorEstimator::error_norm, libMesh::JumpErrorEstimator::fine_context, libMesh::FEMap::inverse_map(), n_vars, libMesh::FEAbstract::reinit(), libMesh::SCALAR, and libMesh::SystemNorm::weight().

Referenced by libMesh::JumpErrorEstimator::estimate_error().

451 {
452  fine_context->side_fe_reinit();
453 
454  unsigned short dim = fine_context->get_elem().dim();
455  libmesh_assert_equal_to(dim, coarse_context->get_elem().dim());
456 
457  FEBase * fe_fine = nullptr;
458  fine_context->get_side_fe( 0, fe_fine, dim );
459 
460  // Get the physical locations of the fine element quadrature points
461  std::vector<Point> qface_point = fe_fine->get_xyz();
462 
463  // Find the master quadrature point locations on the coarse element
464  FEBase * fe_coarse = nullptr;
465  coarse_context->get_side_fe( 0, fe_coarse, dim );
466 
467  std::vector<Point> qp_coarse;
468 
469  FEMap::inverse_map (coarse_context->get_elem().dim(),
470  &coarse_context->get_elem(), qface_point,
471  qp_coarse);
472 
473  // The number of variables in the system
474  const unsigned int n_vars = fine_context->n_vars();
475 
476  // Calculate all coarse element shape functions at those locations
477  for (unsigned int v=0; v<n_vars; v++)
478  if (error_norm.weight(v) != 0.0 &&
479  fine_context->get_system().variable_type(v).family != SCALAR)
480  {
481  coarse_context->get_side_fe( v, fe_coarse, dim );
482  fe_coarse->reinit (&coarse_context->get_elem(), &qp_coarse);
483  }
484 }
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
std::unique_ptr< FEMContext > fine_context
Context objects for integrating on the fine and coarse elements sharing a face.
unsigned int dim
static Point inverse_map(const unsigned int dim, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true, const bool extra_checks=true)
Definition: fe_map.C:1626
unsigned int n_vars
std::unique_ptr< FEMContext > coarse_context
FEGenericBase< Real > FEBase
Real weight(unsigned int var) const
Definition: system_norm.C:134

◆ type()

ErrorEstimatorType libMesh::KellyErrorEstimator::type ( ) const
overridevirtual
Returns
The type for the ErrorEstimator subclass.

Implements libMesh::ErrorEstimator.

Definition at line 51 of file kelly_error_estimator.C.

References libMesh::KELLY.

52 {
53  return KELLY;
54 }

Member Data Documentation

◆ _bc_function

std::pair<bool,Real>(* libMesh::KellyErrorEstimator::_bc_function) (const System &system, const Point &p, const std::string &var_name)
protected

Pointer to function that provides BC information.

Definition at line 119 of file kelly_error_estimator.h.

Referenced by attach_flux_bc_function(), and boundary_side_integration().

◆ coarse_context

std::unique_ptr<FEMContext> libMesh::JumpErrorEstimator::coarse_context
protectedinherited

◆ coarse_error

Real libMesh::JumpErrorEstimator::coarse_error
protectedinherited

◆ error_norm

SystemNorm libMesh::ErrorEstimator::error_norm
inherited

When estimating the error in a single system, the error_norm is used to control the scaling and norm choice for each variable.

Not all estimators will support all norm choices. The default scaling is for all variables to be weighted equally. The default norm choice depends on the error estimator.

Part of this functionality was supported via component_scale and sobolev_order in older libMesh versions, and a small part was supported via component_mask in even older versions. Hopefully the encapsulation here will allow us to avoid changing this API again.

Definition at line 158 of file error_estimator.h.

Referenced by libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::AdjointRefinementEstimator::AdjointRefinementEstimator(), libMesh::DiscontinuityMeasure::boundary_side_integration(), boundary_side_integration(), libMesh::DiscontinuityMeasure::DiscontinuityMeasure(), libMesh::JumpErrorEstimator::estimate_error(), libMesh::AdjointResidualErrorEstimator::estimate_error(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::ErrorEstimator::estimate_errors(), libMesh::ExactErrorEstimator::ExactErrorEstimator(), libMesh::ExactErrorEstimator::find_squared_element_error(), libMesh::LaplacianErrorEstimator::init_context(), libMesh::DiscontinuityMeasure::init_context(), init_context(), libMesh::LaplacianErrorEstimator::internal_side_integration(), libMesh::DiscontinuityMeasure::internal_side_integration(), internal_side_integration(), KellyErrorEstimator(), libMesh::LaplacianErrorEstimator::LaplacianErrorEstimator(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::PatchRecoveryErrorEstimator(), libMesh::JumpErrorEstimator::reinit_sides(), and libMesh::UniformRefinementEstimator::UniformRefinementEstimator().

◆ fine_context

std::unique_ptr<FEMContext> libMesh::JumpErrorEstimator::fine_context
protectedinherited

◆ fine_error

Real libMesh::JumpErrorEstimator::fine_error
protectedinherited

◆ integrate_boundary_sides

bool libMesh::JumpErrorEstimator::integrate_boundary_sides
protectedinherited

A boolean flag, by default false, to be set to true if integrations with boundary_side_integration() should be performed.

Definition at line 152 of file jump_error_estimator.h.

Referenced by libMesh::DiscontinuityMeasure::attach_essential_bc_function(), attach_flux_bc_function(), and libMesh::JumpErrorEstimator::estimate_error().

◆ scale_by_n_flux_faces

bool libMesh::JumpErrorEstimator::scale_by_n_flux_faces
inherited

This boolean flag allows you to scale the error indicator result for each element by the number of "flux faces" the element actually has.

This tends to weight more evenly cells which are on the boundaries and thus have fewer contributions to their flux. The value is initialized to false, simply set it to true if you want to use the feature.

Definition at line 99 of file jump_error_estimator.h.

Referenced by libMesh::JumpErrorEstimator::estimate_error().

◆ use_unweighted_quadrature_rules

bool libMesh::JumpErrorEstimator::use_unweighted_quadrature_rules
inherited

This boolean flag allows you to use "unweighted" quadrature rules (sized to exactly integrate unweighted shape functions in master element space) rather than "default" quadrature rules (sized to exactly integrate polynomials of one higher degree than mass matrix terms).

The results with the former, lower-order rules will be somewhat less accurate in many cases but will be much cheaper to compute.

The value is initialized to false, simply set it to true if you want to use the feature.

Definition at line 113 of file jump_error_estimator.h.

Referenced by libMesh::JumpErrorEstimator::estimate_error(), and main().

◆ var

unsigned int libMesh::JumpErrorEstimator::var
protectedinherited

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