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

This class implements an "error estimator" based on the difference between the approximate and exact solution. More...

#include <exact_error_estimator.h>

Inheritance diagram for libMesh::ExactErrorEstimator:
[legend]

Public Types

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

Public Member Functions

 ExactErrorEstimator ()
 Constructor. More...
 
 ~ExactErrorEstimator ()
 Destructor. More...
 
void attach_exact_values (std::vector< FunctionBase< Number > * > f)
 Clone and attach arbitrary functors which compute the exact values of the EquationSystems' solutions at any point. More...
 
void attach_exact_value (unsigned int sys_num, FunctionBase< Number > *f)
 Clone and attach an arbitrary functor which computes the exact value of the system sys_num solution at any point. More...
 
void attach_exact_value (Number fptr(const Point &p, const Parameters &Parameters, const std::string &sys_name, const std::string &unknown_name))
 Attach an arbitrary function which computes the exact value of the solution at any point. More...
 
void attach_exact_derivs (std::vector< FunctionBase< Gradient > * > g)
 Clone and attach arbitrary functors which compute the exact gradients of the EquationSystems' solutions at any point. More...
 
void attach_exact_deriv (unsigned int sys_num, FunctionBase< Gradient > *g)
 Clone and attach an arbitrary functor which computes the exact gradient of the system sys_num solution at any point. More...
 
void attach_exact_deriv (Gradient gptr(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name))
 Attach an arbitrary function which computes the exact gradient of the solution at any point. More...
 
void attach_exact_hessians (std::vector< FunctionBase< Tensor > * > h)
 Clone and attach arbitrary functors which compute the exact second derivatives of the EquationSystems' solutions at any point. More...
 
void attach_exact_hessian (unsigned int sys_num, FunctionBase< Tensor > *h)
 Clone and attach an arbitrary functor which computes the exact second derivatives of the system sys_num solution at any point. More...
 
void attach_exact_hessian (Tensor hptr(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name))
 Attach an arbitrary function which computes the exact second derivatives of the solution at any point. More...
 
void attach_reference_solution (EquationSystems *es_fine)
 Attach function similar to system.h which allows the user to attach a second EquationSystems object with a reference fine grid solution. More...
 
void extra_quadrature_order (const int extraorder)
 Increases or decreases the order of the quadrature rule used for numerical integration. More...
 
virtual void estimate_error (const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=libmesh_nullptr, bool estimate_parent_error=false) libmesh_override
 This function uses the exact solution function to estimate the error on each cell. More...
 
virtual ErrorEstimatorType type () const libmesh_override
 
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=libmesh_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=libmesh_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

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

void reduce_error (std::vector< ErrorVectorReal > &error_per_cell, const Parallel::Communicator &comm LIBMESH_CAN_DEFAULT_TO_COMMWORLD) 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...
 

Private Member Functions

Real find_squared_element_error (const System &system, const std::string &var_name, const Elem *elem, const DenseVector< Number > &Uelem, FEBase *fe, MeshFunction *fine_values) const
 Helper method for calculating on each element. More...
 
void clear_functors ()
 Helper method for cleanup. More...
 

Private Attributes

Number(* _exact_value )(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
 Function pointer to user-provided function which computes the exact value of the solution. More...
 
Gradient(* _exact_deriv )(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
 Function pointer to user-provided function which computes the exact derivative of the solution. More...
 
Tensor(* _exact_hessian )(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
 Function pointer to user-provided function which computes the exact hessian of the solution. More...
 
std::vector< FunctionBase< Number > * > _exact_values
 User-provided functors which compute the exact value of the solution for each system. More...
 
std::vector< FunctionBase< Gradient > * > _exact_derivs
 User-provided functors which compute the exact derivative of the solution for each system. More...
 
std::vector< FunctionBase< Tensor > * > _exact_hessians
 User-provided functors which compute the exact hessians of the solution for each system. More...
 
EquationSystems_equation_systems_fine
 Constant pointer to the EquationSystems object containing a fine grid solution. More...
 
int _extra_order
 Extra order to use for quadrature rule. More...
 

Detailed Description

This class implements an "error estimator" based on the difference between the approximate and exact solution.

In theory the quadrature error in this estimate should be much lower than the approximation error in other estimates, so this estimator can be used to calculate effectivity.

Author
Roy Stogner
Date
2006

Definition at line 67 of file exact_error_estimator.h.

Member Typedef Documentation

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

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

Definition at line 113 of file error_estimator.h.

Constructor & Destructor Documentation

libMesh::ExactErrorEstimator::ExactErrorEstimator ( )

Constructor.

Responsible for initializing the _bc_function function pointer to libmesh_nullptr, and defaulting the norm type to H1.

Definition at line 75 of file exact_error_estimator.h.

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

75  :
81  _extra_order(0)
82  { error_norm = H1; }
Number(* _exact_value)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
Function pointer to user-provided function which computes the exact value of the solution.
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
int _extra_order
Extra order to use for quadrature rule.
ErrorEstimator()
Constructor.
const class libmesh_nullptr_t libmesh_nullptr
EquationSystems * _equation_systems_fine
Constant pointer to the EquationSystems object containing a fine grid solution.
Tensor(* _exact_hessian)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
Function pointer to user-provided function which computes the exact hessian of the solution...
Gradient(* _exact_deriv)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
Function pointer to user-provided function which computes the exact derivative of the solution...
libMesh::ExactErrorEstimator::~ExactErrorEstimator ( )

Member Function Documentation

void libMesh::ExactErrorEstimator::attach_exact_deriv ( unsigned int  sys_num,
FunctionBase< Gradient > *  g 
)

Clone and attach an arbitrary functor which computes the exact gradient of the system sys_num solution at any point.

Definition at line 127 of file exact_error_estimator.C.

References _equation_systems_fine, _exact_derivs, _exact_hessian, attach_exact_hessian(), clear_functors(), libMesh::FunctionBase< Output >::clone(), libMesh::libmesh_assert(), and libmesh_nullptr.

Referenced by attach_exact_value(), main(), and ~ExactErrorEstimator().

129 {
130  if (_exact_derivs.size() <= sys_num)
131  _exact_derivs.resize(sys_num+1, libmesh_nullptr);
132 
133  if (g)
134  _exact_derivs[sys_num] = g->clone().release();
135 }
std::vector< FunctionBase< Gradient > * > _exact_derivs
User-provided functors which compute the exact derivative of the solution for each system...
const class libmesh_nullptr_t libmesh_nullptr
virtual UniquePtr< FunctionBase< Output > > clone() const =0
void libMesh::ExactErrorEstimator::attach_exact_deriv ( Gradient   gptrconst Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)

Attach an arbitrary function which computes the exact gradient of the solution at any point.

void libMesh::ExactErrorEstimator::attach_exact_derivs ( std::vector< FunctionBase< Gradient > * >  g)

Clone and attach arbitrary functors which compute the exact gradients of the EquationSystems' solutions at any point.

Definition at line 108 of file exact_error_estimator.C.

References _exact_derivs, and libmesh_nullptr.

Referenced by ~ExactErrorEstimator().

109 {
110  // Clear out any previous _exact_derivs entries, then add a new
111  // entry for each system.
112  for (std::size_t i=0; i != _exact_derivs.size(); ++i)
113  delete (_exact_derivs[i]);
114 
115  _exact_derivs.clear();
116  _exact_derivs.resize(g.size(), libmesh_nullptr);
117 
118  // We use clone() to get non-sliced copies of FunctionBase
119  // subclasses, but we don't currently put the resulting UniquePtrs
120  // into an STL container.
121  for (std::size_t i=0; i != g.size(); ++i)
122  if (g[i])
123  _exact_derivs[i] = g[i]->clone().release();
124 }
std::vector< FunctionBase< Gradient > * > _exact_derivs
User-provided functors which compute the exact derivative of the solution for each system...
const class libmesh_nullptr_t libmesh_nullptr
void libMesh::ExactErrorEstimator::attach_exact_hessian ( unsigned int  sys_num,
FunctionBase< Tensor > *  h 
)

Clone and attach an arbitrary functor which computes the exact second derivatives of the system sys_num solution at any point.

Definition at line 175 of file exact_error_estimator.C.

References _exact_hessians, libMesh::FunctionBase< Output >::clone(), and libmesh_nullptr.

Referenced by attach_exact_deriv(), and ~ExactErrorEstimator().

177 {
178  if (_exact_hessians.size() <= sys_num)
179  _exact_hessians.resize(sys_num+1, libmesh_nullptr);
180 
181  if (h)
182  _exact_hessians[sys_num] = h->clone().release();
183 }
std::vector< FunctionBase< Tensor > * > _exact_hessians
User-provided functors which compute the exact hessians of the solution for each system.
const class libmesh_nullptr_t libmesh_nullptr
void libMesh::ExactErrorEstimator::attach_exact_hessian ( Tensor   hptrconst Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)

Attach an arbitrary function which computes the exact second derivatives of the solution at any point.

void libMesh::ExactErrorEstimator::attach_exact_hessians ( std::vector< FunctionBase< Tensor > * >  h)

Clone and attach arbitrary functors which compute the exact second derivatives of the EquationSystems' solutions at any point.

Definition at line 156 of file exact_error_estimator.C.

References _exact_hessians, and libmesh_nullptr.

Referenced by ~ExactErrorEstimator().

157 {
158  // Clear out any previous _exact_hessians entries, then add a new
159  // entry for each system.
160  for (std::size_t i=0; i != _exact_hessians.size(); ++i)
161  delete (_exact_hessians[i]);
162 
163  _exact_hessians.clear();
164  _exact_hessians.resize(h.size(), libmesh_nullptr);
165 
166  // We use clone() to get non-sliced copies of FunctionBase
167  // subclasses, but we don't currently put the resulting UniquePtrs
168  // into an STL container.
169  for (std::size_t i=0; i != h.size(); ++i)
170  if (h[i])
171  _exact_hessians[i] = h[i]->clone().release();
172 }
std::vector< FunctionBase< Tensor > * > _exact_hessians
User-provided functors which compute the exact hessians of the solution for each system.
const class libmesh_nullptr_t libmesh_nullptr
void libMesh::ExactErrorEstimator::attach_exact_value ( unsigned int  sys_num,
FunctionBase< Number > *  f 
)

Clone and attach an arbitrary functor which computes the exact value of the system sys_num solution at any point.

Definition at line 81 of file exact_error_estimator.C.

References _equation_systems_fine, _exact_deriv, _exact_values, attach_exact_deriv(), clear_functors(), libMesh::FunctionBase< Output >::clone(), gptr(), libMesh::libmesh_assert(), and libmesh_nullptr.

Referenced by main(), and ~ExactErrorEstimator().

83 {
84  if (_exact_values.size() <= sys_num)
85  _exact_values.resize(sys_num+1, libmesh_nullptr);
86 
87  if (f)
88  _exact_values[sys_num] = f->clone().release();
89 }
std::vector< FunctionBase< Number > * > _exact_values
User-provided functors which compute the exact value of the solution for each system.
const class libmesh_nullptr_t libmesh_nullptr
virtual UniquePtr< FunctionBase< Output > > clone() const =0
void libMesh::ExactErrorEstimator::attach_exact_value ( Number   fptrconst Point &p, const Parameters &Parameters, const std::string &sys_name, const std::string &unknown_name)

Attach an arbitrary function which computes the exact value of the solution at any point.

void libMesh::ExactErrorEstimator::attach_exact_values ( std::vector< FunctionBase< Number > * >  f)

Clone and attach arbitrary functors which compute the exact values of the EquationSystems' solutions at any point.

Definition at line 62 of file exact_error_estimator.C.

References _exact_values, and libmesh_nullptr.

Referenced by ~ExactErrorEstimator().

63 {
64  // Clear out any previous _exact_values entries, then add a new
65  // entry for each system.
66  for (std::size_t i=0; i != _exact_values.size(); ++i)
67  delete (_exact_values[i]);
68 
69  _exact_values.clear();
70  _exact_values.resize(f.size(), libmesh_nullptr);
71 
72  // We use clone() to get non-sliced copies of FunctionBase
73  // subclasses, but we don't currently put the resulting UniquePtrs
74  // into an STL container.
75  for (std::size_t i=0; i != f.size(); ++i)
76  if (f[i])
77  _exact_values[i] = f[i]->clone().release();
78 }
std::vector< FunctionBase< Number > * > _exact_values
User-provided functors which compute the exact value of the solution for each system.
const class libmesh_nullptr_t libmesh_nullptr
void libMesh::ExactErrorEstimator::attach_reference_solution ( EquationSystems es_fine)

Attach function similar to system.h which allows the user to attach a second EquationSystems object with a reference fine grid solution.

Definition at line 186 of file exact_error_estimator.C.

References _equation_systems_fine, _exact_deriv, _exact_hessian, _exact_value, clear_functors(), libMesh::libmesh_assert(), and libmesh_nullptr.

Referenced by ~ExactErrorEstimator().

187 {
188  libmesh_assert(es_fine);
189  _equation_systems_fine = es_fine;
190 
191  // If we're using a fine grid solution, we're not using exact value
192  // function pointers or functors.
196 
197  this->clear_functors();
198 }
Number(* _exact_value)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
Function pointer to user-provided function which computes the exact value of the solution.
const class libmesh_nullptr_t libmesh_nullptr
libmesh_assert(j)
void clear_functors()
Helper method for cleanup.
EquationSystems * _equation_systems_fine
Constant pointer to the EquationSystems object containing a fine grid solution.
Tensor(* _exact_hessian)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
Function pointer to user-provided function which computes the exact hessian of the solution...
Gradient(* _exact_deriv)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
Function pointer to user-provided function which computes the exact derivative of the solution...
void libMesh::ExactErrorEstimator::clear_functors ( )
private

Helper method for cleanup.

Definition at line 547 of file exact_error_estimator.C.

References _exact_derivs, _exact_hessians, and _exact_values.

Referenced by attach_exact_deriv(), attach_exact_value(), and attach_reference_solution().

548 {
549  // delete will clean up any cloned functors and no-op on any NULL
550  // pointers
551 
552  for (std::size_t i=0; i != _exact_values.size(); ++i)
553  delete (_exact_values[i]);
554 
555  for (std::size_t i=0; i != _exact_derivs.size(); ++i)
556  delete (_exact_derivs[i]);
557 
558  for (std::size_t i=0; i != _exact_hessians.size(); ++i)
559  delete (_exact_hessians[i]);
560 }
std::vector< FunctionBase< Number > * > _exact_values
User-provided functors which compute the exact value of the solution for each system.
std::vector< FunctionBase< Tensor > * > _exact_hessians
User-provided functors which compute the exact hessians of the solution for each system.
std::vector< FunctionBase< Gradient > * > _exact_derivs
User-provided functors which compute the exact derivative of the solution for each system...
void libMesh::ExactErrorEstimator::estimate_error ( const System system,
ErrorVector error_per_cell,
const NumericVector< Number > *  solution_vector = libmesh_nullptr,
bool  estimate_parent_error = false 
)
virtual

This function uses the exact solution function to estimate the error on each cell.

The estimated error is output in the vector error_per_cell

Implements libMesh::ErrorEstimator.

Definition at line 200 of file exact_error_estimator.C.

References _equation_systems_fine, _exact_derivs, _exact_hessians, _exact_values, _extra_order, libMesh::MeshBase::active_local_element_ptr_range(), libMesh::FEGenericBase< OutputType >::build(), libMesh::NumericVector< T >::build(), libMesh::Elem::child_ref_range(), libMesh::FEGenericBase< OutputType >::coarsened_dof_values(), libMesh::ParallelObject::comm(), libMesh::System::current_local_solution, libMesh::System::current_solution(), libMesh::FEType::default_quadrature_rule(), dim, libMesh::DofMap::dof_indices(), libMesh::ErrorEstimator::error_norm, libMesh::ErrorVectorReal, find_squared_element_error(), libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMesh::EquationSystems::get_system(), libMesh::DofObject::id(), libMesh::libmesh_ignore(), libMesh::MeshBase::max_elem_id(), mesh, libMesh::MeshBase::mesh_dimension(), n_vars, libMesh::System::n_vars(), libMesh::System::name(), libMesh::Elem::parent(), libMesh::ErrorEstimator::reduce_error(), libMesh::SERIAL, libMesh::System::solution, libMesh::NumericVector< T >::swap(), libMesh::sys, libMesh::System::update_global_solution(), libMesh::System::variable_name(), libMesh::System::variable_number(), libMesh::DofMap::variable_type(), and libMesh::SystemNorm::weight().

Referenced by extra_quadrature_order(), and main().

204 {
205  // Ignore the fact that this variable is unused when !LIBMESH_ENABLE_AMR
206  libmesh_ignore(estimate_parent_error);
207 
208  // The current mesh
209  const MeshBase & mesh = system.get_mesh();
210 
211  // The dimensionality of the mesh
212  const unsigned int dim = mesh.mesh_dimension();
213 
214  // The number of variables in the system
215  const unsigned int n_vars = system.n_vars();
216 
217  // The DofMap for this system
218  const DofMap & dof_map = system.get_dof_map();
219 
220  // Resize the error_per_cell vector to be
221  // the number of elements, initialize it to 0.
222  error_per_cell.resize (mesh.max_elem_id());
223  std::fill (error_per_cell.begin(), error_per_cell.end(), 0.);
224 
225  // Prepare current_local_solution to localize a non-standard
226  // solution vector if necessary
227  if (solution_vector && solution_vector != system.solution.get())
228  {
229  NumericVector<Number> * newsol =
230  const_cast<NumericVector<Number> *>(solution_vector);
231  System & sys = const_cast<System &>(system);
232  newsol->swap(*sys.solution);
233  sys.update();
234  }
235 
236  // Loop over all the variables in the system
237  for (unsigned int var=0; var<n_vars; var++)
238  {
239  // Possibly skip this variable
240  if (error_norm.weight(var) == 0.0) continue;
241 
242  // The (string) name of this variable
243  const std::string & var_name = system.variable_name(var);
244 
245  // The type of finite element to use for this variable
246  const FEType & fe_type = dof_map.variable_type (var);
247 
248  UniquePtr<FEBase> fe (FEBase::build (dim, fe_type));
249 
250  // Build an appropriate Gaussian quadrature rule
251  UniquePtr<QBase> qrule =
252  fe_type.default_quadrature_rule (dim,
253  _extra_order);
254 
255  fe->attach_quadrature_rule (qrule.get());
256 
257  // Prepare a global solution and a MeshFunction of the fine system if we need one
258  UniquePtr<MeshFunction> fine_values;
259  UniquePtr<NumericVector<Number>> fine_soln = NumericVector<Number>::build(system.comm());
261  {
262  const System & fine_system = _equation_systems_fine->get_system(system.name());
263 
264  std::vector<Number> global_soln;
265  // FIXME - we're assuming that the fine system solution gets
266  // used even when a different vector is used for the coarse
267  // system
268  fine_system.update_global_solution(global_soln);
269  fine_soln->init
270  (cast_int<numeric_index_type>(global_soln.size()), true,
271  SERIAL);
272  (*fine_soln) = global_soln;
273 
274  fine_values = UniquePtr<MeshFunction>
275  (new MeshFunction(*_equation_systems_fine,
276  *fine_soln,
277  fine_system.get_dof_map(),
278  fine_system.variable_number(var_name)));
279  fine_values->init();
280  } else {
281  // Initialize functors if we're using them
282  for (std::size_t i=0; i != _exact_values.size(); ++i)
283  if (_exact_values[i])
284  _exact_values[i]->init();
285 
286  for (std::size_t i=0; i != _exact_derivs.size(); ++i)
287  if (_exact_derivs[i])
288  _exact_derivs[i]->init();
289 
290  for (std::size_t i=0; i != _exact_hessians.size(); ++i)
291  if (_exact_hessians[i])
292  _exact_hessians[i]->init();
293  }
294 
295  // Request the data we'll need to compute with
296  fe->get_JxW();
297  fe->get_phi();
298  fe->get_dphi();
299 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
300  fe->get_d2phi();
301 #endif
302  fe->get_xyz();
303 
304 #ifdef LIBMESH_ENABLE_AMR
305  // If we compute on parent elements, we'll want to do so only
306  // once on each, so we need to keep track of which we've done.
307  std::vector<bool> computed_var_on_parent;
308 
309  if (estimate_parent_error)
310  computed_var_on_parent.resize(error_per_cell.size(), false);
311 #endif
312 
313  // TODO: this ought to be threaded (and using subordinate
314  // MeshFunction objects in each thread rather than a single
315  // master)
316 
317  // Iterate over all the active elements in the mesh
318  // that live on this processor.
319  for (const auto & elem : mesh.active_local_element_ptr_range())
320  {
321  const dof_id_type e_id = elem->id();
322 
323 #ifdef LIBMESH_ENABLE_AMR
324  // See if the parent of element e has been examined yet;
325  // if not, we may want to compute the estimator on it
326  const Elem * parent = elem->parent();
327 
328  // We only can compute and only need to compute on
329  // parents with all active children
330  bool compute_on_parent = true;
331  if (!parent || !estimate_parent_error)
332  compute_on_parent = false;
333  else
334  for (auto & child : parent->child_ref_range())
335  if (!child.active())
336  compute_on_parent = false;
337 
338  if (compute_on_parent &&
339  !computed_var_on_parent[parent->id()])
340  {
341  computed_var_on_parent[parent->id()] = true;
342 
343  // Compute a projection onto the parent
344  DenseVector<Number> Uparent;
345  FEBase::coarsened_dof_values(*(system.current_local_solution),
346  dof_map, parent, Uparent,
347  var, false);
348 
349  error_per_cell[parent->id()] +=
350  static_cast<ErrorVectorReal>
351  (find_squared_element_error(system, var_name,
352  parent, Uparent,
353  fe.get(),
354  fine_values.get()));
355  }
356 #endif
357 
358  // Get the local to global degree of freedom maps
359  std::vector<dof_id_type> dof_indices;
360  dof_map.dof_indices (elem, dof_indices, var);
361  const unsigned int n_dofs =
362  cast_int<unsigned int>(dof_indices.size());
363  DenseVector<Number> Uelem(n_dofs);
364  for (unsigned int i=0; i != n_dofs; ++i)
365  Uelem(i) = system.current_solution(dof_indices[i]);
366 
367  error_per_cell[e_id] +=
368  static_cast<ErrorVectorReal>
369  (find_squared_element_error(system, var_name, elem,
370  Uelem, fe.get(),
371  fine_values.get()));
372 
373  } // End loop over active local elements
374  } // End loop over variables
375 
376 
377 
378  // Each processor has now computed the error contributions
379  // for its local elements. We need to sum the vector
380  // and then take the square-root of each component. Note
381  // that we only need to sum if we are running on multiple
382  // processors, and we only need to take the square-root
383  // if the value is nonzero. There will in general be many
384  // zeros for the inactive elements.
385 
386  // First sum the vector of estimated error values
387  this->reduce_error(error_per_cell, system.comm());
388 
389  // Compute the square-root of each component.
390  {
391  LOG_SCOPE("std::sqrt()", "ExactErrorEstimator");
392  for (std::size_t i=0; i<error_per_cell.size(); i++)
393  if (error_per_cell[i] != 0.)
394  {
395  libmesh_assert_greater (error_per_cell[i], 0.);
396  error_per_cell[i] = std::sqrt(error_per_cell[i]);
397  }
398  }
399 
400  // If we used a non-standard solution before, now is the time to fix
401  // the current_local_solution
402  if (solution_vector && solution_vector != system.solution.get())
403  {
404  NumericVector<Number> * newsol =
405  const_cast<NumericVector<Number> *>(solution_vector);
406  System & sys = const_cast<System &>(system);
407  newsol->swap(*sys.solution);
408  sys.update();
409  }
410 }
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
std::vector< FunctionBase< Number > * > _exact_values
User-provided functors which compute the exact value of the solution for each system.
int _extra_order
Extra order to use for quadrature rule.
unsigned int dim
ImplicitSystem & sys
Real weight(unsigned int var) const
Definition: system_norm.h:292
std::vector< FunctionBase< Tensor > * > _exact_hessians
User-provided functors which compute the exact hessians of the solution for each system.
std::vector< FunctionBase< Gradient > * > _exact_derivs
User-provided functors which compute the exact derivative of the solution for each system...
MeshBase & mesh
const unsigned int n_vars
Definition: tecplot_io.C:68
DIE A HORRIBLE DEATH HERE typedef float ErrorVectorReal
static UniquePtr< NumericVector< T > > build(const Parallel::Communicator &comm, const SolverPackage solver_package=libMesh::default_solver_package())
Builds a NumericVector on the processors in communicator comm using the linear solver package specifi...
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:802
void reduce_error(std::vector< ErrorVectorReal > &error_per_cell, const Parallel::Communicator &comm LIBMESH_CAN_DEFAULT_TO_COMMWORLD) const
This method takes the local error contributions in error_per_cell from each processor and combines th...
void libmesh_ignore(const T &)
virtual void update()
Update the local values to reflect the solution on neighboring processors.
Definition: system.C:425
EquationSystems * _equation_systems_fine
Constant pointer to the EquationSystems object containing a fine grid solution.
virtual void get(const std::vector< numeric_index_type > &index, T *values) const
Access multiple components at once.
const T_sys & get_system(const std::string &name) const
Real find_squared_element_error(const System &system, const std::string &var_name, const Elem *elem, const DenseVector< Number > &Uelem, FEBase *fe, MeshFunction *fine_values) const
Helper method for calculating on each element.
uint8_t dof_id_type
Definition: id_types.h:64
static UniquePtr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
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 = libmesh_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 48 of file error_estimator.C.

References libMesh::ErrorEstimator::error_norm, libMesh::ErrorEstimator::estimate_error(), libMesh::EquationSystems::get_system(), libmesh_nullptr, libMesh::EquationSystems::n_systems(), and libMesh::sys.

Referenced by libMesh::ErrorEstimator::~ErrorEstimator().

53 {
54  SystemNorm old_error_norm = this->error_norm;
55 
56  // Sum the error values from each system
57  for (unsigned int s = 0; s != equation_systems.n_systems(); ++s)
58  {
59  ErrorVector system_error_per_cell;
60  const System & sys = equation_systems.get_system(s);
61  if (error_norms.find(&sys) == error_norms.end())
62  this->error_norm = old_error_norm;
63  else
64  this->error_norm = error_norms.find(&sys)->second;
65 
66  const NumericVector<Number> * solution_vector = libmesh_nullptr;
67  if (solution_vectors &&
68  solution_vectors->find(&sys) != solution_vectors->end())
69  solution_vector = solution_vectors->find(&sys)->second;
70 
71  this->estimate_error(sys, system_error_per_cell,
72  solution_vector, estimate_parent_error);
73 
74  if (s)
75  {
76  libmesh_assert_equal_to (error_per_cell.size(), system_error_per_cell.size());
77  for (std::size_t i=0; i != error_per_cell.size(); ++i)
78  error_per_cell[i] += system_error_per_cell[i];
79  }
80  else
81  error_per_cell = system_error_per_cell;
82  }
83 
84  // Restore our old state before returning
85  this->error_norm = old_error_norm;
86 }
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
ImplicitSystem & sys
const class libmesh_nullptr_t libmesh_nullptr
virtual void estimate_error(const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=libmesh_nullptr, bool estimate_parent_error=false)=0
This pure virtual function must be redefined in derived classes to compute the error for each active ...
void libMesh::ErrorEstimator::estimate_errors ( const EquationSystems equation_systems,
ErrorMap errors_per_cell,
const std::map< const System *, const NumericVector< Number > * > *  solution_vectors = libmesh_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 94 of file error_estimator.C.

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

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

Increases or decreases the order of the quadrature rule used for numerical integration.

Definition at line 167 of file exact_error_estimator.h.

References _extra_order, estimate_error(), and libmesh_nullptr.

168  { _extra_order = extraorder; }
int _extra_order
Extra order to use for quadrature rule.
Real libMesh::ExactErrorEstimator::find_squared_element_error ( const System system,
const std::string &  var_name,
const Elem elem,
const DenseVector< Number > &  Uelem,
FEBase fe,
MeshFunction fine_values 
) const
private

Helper method for calculating on each element.

Definition at line 414 of file exact_error_estimator.C.

References _equation_systems_fine, _exact_deriv, _exact_derivs, _exact_hessian, _exact_hessians, _exact_value, _exact_values, libMesh::ErrorEstimator::error_norm, libMesh::FEGenericBase< OutputType >::get_d2phi(), libMesh::FEGenericBase< OutputType >::get_dphi(), libMesh::System::get_equation_systems(), libMesh::FEAbstract::get_JxW(), libMesh::FEGenericBase< OutputType >::get_phi(), libMesh::FEAbstract::get_xyz(), libMesh::MeshFunction::gradient(), libMesh::H1, libMesh::H1_SEMINORM, libMesh::H2, libMesh::H2_SEMINORM, libMesh::MeshFunction::hessian(), libMesh::L2, libMesh::System::name(), libMesh::TensorTools::norm_sq(), libMesh::TypeVector< T >::norm_sq(), libMesh::TypeTensor< T >::norm_sq(), libMesh::System::number(), libMesh::EquationSystems::parameters, libMesh::Real, libMesh::FEAbstract::reinit(), libMesh::DenseVector< T >::size(), libMesh::System::time, libMesh::SystemNorm::type(), libMesh::System::variable_number(), and libMesh::System::variable_scalar_number().

Referenced by estimate_error().

420 {
421  // The (string) name of this system
422  const std::string & sys_name = system.name();
423  const unsigned int sys_num = system.number();
424 
425  const unsigned int var = system.variable_number(var_name);
426  const unsigned int var_component =
427  system.variable_scalar_number(var, 0);
428 
429  const Parameters & parameters = system.get_equation_systems().parameters;
430 
431  // reinitialize the element-specific data
432  // for the current element
433  fe->reinit (elem);
434 
435  // Get the data we need to compute with
436  const std::vector<Real> & JxW = fe->get_JxW();
437  const std::vector<std::vector<Real>> & phi_values = fe->get_phi();
438  const std::vector<std::vector<RealGradient>> & dphi_values = fe->get_dphi();
439  const std::vector<Point> & q_point = fe->get_xyz();
440 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
441  const std::vector<std::vector<RealTensor>> & d2phi_values = fe->get_d2phi();
442 #endif
443 
444  // The number of shape functions
445  const unsigned int n_sf =
446  cast_int<unsigned int>(Uelem.size());
447 
448  // The number of quadrature points
449  const unsigned int n_qp =
450  cast_int<unsigned int>(JxW.size());
451 
452  Real error_val = 0;
453 
454  // Begin the loop over the Quadrature points.
455  //
456  for (unsigned int qp=0; qp<n_qp; qp++)
457  {
458  // Real u_h = 0.;
459  // RealGradient grad_u_h;
460 
461  Number u_h = 0.;
462 
463  Gradient grad_u_h;
464 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
465  Tensor grad2_u_h;
466 #endif
467 
468  // Compute solution values at the current
469  // quadrature point. This requires a sum
470  // over all the shape functions evaluated
471  // at the quadrature point.
472  for (unsigned int i=0; i<n_sf; i++)
473  {
474  // Values from current solution.
475  u_h += phi_values[i][qp]*Uelem(i);
476  grad_u_h += dphi_values[i][qp]*Uelem(i);
477 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
478  grad2_u_h += d2phi_values[i][qp]*Uelem(i);
479 #endif
480  }
481 
482  // Compute the value of the error at this quadrature point
483  if (error_norm.type(var) == L2 ||
484  error_norm.type(var) == H1 ||
485  error_norm.type(var) == H2)
486  {
487  Number val_error = u_h;
488  if (_exact_value)
489  val_error -= _exact_value(q_point[qp],parameters,sys_name,var_name);
490  else if (_exact_values.size() > sys_num && _exact_values[sys_num])
491  val_error -= _exact_values[sys_num]->
492  component(var_component, q_point[qp], system.time);
493  else if (_equation_systems_fine)
494  val_error -= (*fine_values)(q_point[qp]);
495 
496  // Add the squares of the error to each contribution
497  error_val += JxW[qp]*TensorTools::norm_sq(val_error);
498  }
499 
500  // Compute the value of the error in the gradient at this
501  // quadrature point
502  if (error_norm.type(var) == H1 ||
503  error_norm.type(var) == H1_SEMINORM ||
504  error_norm.type(var) == H2)
505  {
506  Gradient grad_error = grad_u_h;
507  if (_exact_deriv)
508  grad_error -= _exact_deriv(q_point[qp],parameters,sys_name,var_name);
509  else if (_exact_derivs.size() > sys_num && _exact_derivs[sys_num])
510  grad_error -= _exact_derivs[sys_num]->
511  component(var_component, q_point[qp], system.time);
512  else if (_equation_systems_fine)
513  grad_error -= fine_values->gradient(q_point[qp]);
514 
515  error_val += JxW[qp]*grad_error.norm_sq();
516  }
517 
518 
519 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
520  // Compute the value of the error in the hessian at this
521  // quadrature point
522  if ((error_norm.type(var) == H2_SEMINORM ||
523  error_norm.type(var) == H2))
524  {
525  Tensor grad2_error = grad2_u_h;
526  if (_exact_hessian)
527  grad2_error -= _exact_hessian(q_point[qp],parameters,sys_name,var_name);
528  else if (_exact_hessians.size() > sys_num && _exact_hessians[sys_num])
529  grad2_error -= _exact_hessians[sys_num]->
530  component(var_component, q_point[qp], system.time);
531  else if (_equation_systems_fine)
532  grad2_error -= fine_values->hessian(q_point[qp]);
533 
534  error_val += JxW[qp]*grad2_error.norm_sq();
535  }
536 #endif
537 
538  } // end qp loop
539 
540  libmesh_assert_greater_equal (error_val, 0.);
541 
542  return error_val;
543 }
Number(* _exact_value)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
Function pointer to user-provided function which computes the exact value of the solution.
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
std::vector< FunctionBase< Number > * > _exact_values
User-provided functors which compute the exact value of the solution for each system.
std::vector< FunctionBase< Tensor > * > _exact_hessians
User-provided functors which compute the exact hessians of the solution for each system.
std::vector< FunctionBase< Gradient > * > _exact_derivs
User-provided functors which compute the exact derivative of the solution for each system...
virtual unsigned int size() const libmesh_override
Definition: dense_vector.h:87
NumberVectorValue Gradient
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
EquationSystems * _equation_systems_fine
Constant pointer to the EquationSystems object containing a fine grid solution.
Tensor(* _exact_hessian)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
Function pointer to user-provided function which computes the exact hessian of the solution...
NumberTensorValue Tensor
FEMNormType type(unsigned int var) const
Definition: system_norm.h:268
Gradient(* _exact_deriv)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
Function pointer to user-provided function which computes the exact derivative of the solution...
void libMesh::ErrorEstimator::reduce_error ( std::vector< ErrorVectorReal > &  error_per_cell,
const Parallel::Communicator &comm  LIBMESH_CAN_DEFAULT_TO_COMMWORLD 
) 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 33 of file error_estimator.C.

References libMesh::Parallel::Communicator::sum().

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

35 {
36  // This function must be run on all processors at once
37  // parallel_object_only();
38 
39  // Each processor has now computed the error contributions
40  // for its local elements. We may need to sum the vector to
41  // recover the error for each element.
42 
43  comm.sum(error_per_cell);
44 }
virtual ErrorEstimatorType libMesh::ExactErrorEstimator::type ( ) const
virtual
Returns
The type for the ErrorEstimator subclass.

Implements libMesh::ErrorEstimator.

Definition at line 189 of file exact_error_estimator.h.

References libMesh::EXACT.

Member Data Documentation

EquationSystems* libMesh::ExactErrorEstimator::_equation_systems_fine
private

Constant pointer to the EquationSystems object containing a fine grid solution.

Definition at line 243 of file exact_error_estimator.h.

Referenced by attach_exact_deriv(), attach_exact_value(), attach_reference_solution(), estimate_error(), and find_squared_element_error().

Gradient(* libMesh::ExactErrorEstimator::_exact_deriv) (const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
private

Function pointer to user-provided function which computes the exact derivative of the solution.

Definition at line 207 of file exact_error_estimator.h.

Referenced by attach_exact_value(), attach_reference_solution(), and find_squared_element_error().

std::vector<FunctionBase<Gradient> *> libMesh::ExactErrorEstimator::_exact_derivs
private

User-provided functors which compute the exact derivative of the solution for each system.

Definition at line 231 of file exact_error_estimator.h.

Referenced by attach_exact_deriv(), attach_exact_derivs(), clear_functors(), estimate_error(), and find_squared_element_error().

Tensor(* libMesh::ExactErrorEstimator::_exact_hessian) (const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
private

Function pointer to user-provided function which computes the exact hessian of the solution.

Definition at line 216 of file exact_error_estimator.h.

Referenced by attach_exact_deriv(), attach_reference_solution(), and find_squared_element_error().

std::vector<FunctionBase<Tensor> *> libMesh::ExactErrorEstimator::_exact_hessians
private

User-provided functors which compute the exact hessians of the solution for each system.

Definition at line 237 of file exact_error_estimator.h.

Referenced by attach_exact_hessian(), attach_exact_hessians(), clear_functors(), estimate_error(), and find_squared_element_error().

Number(* libMesh::ExactErrorEstimator::_exact_value) (const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
private

Function pointer to user-provided function which computes the exact value of the solution.

Definition at line 198 of file exact_error_estimator.h.

Referenced by attach_reference_solution(), and find_squared_element_error().

std::vector<FunctionBase<Number> *> libMesh::ExactErrorEstimator::_exact_values
private

User-provided functors which compute the exact value of the solution for each system.

Definition at line 225 of file exact_error_estimator.h.

Referenced by attach_exact_value(), attach_exact_values(), clear_functors(), estimate_error(), and find_squared_element_error().

int libMesh::ExactErrorEstimator::_extra_order
private

Extra order to use for quadrature rule.

Definition at line 263 of file exact_error_estimator.h.

Referenced by estimate_error(), and extra_quadrature_order().

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 150 of file error_estimator.h.

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


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