libMesh::ExactErrorEstimator Class Reference

#include <exact_error_estimator.h>

Inheritance diagram for libMesh::ExactErrorEstimator:

Public Types

typedef Number(* ValueFunctionPointer) (const Point &p, const Parameters &Parameters, const std::string &sys_name, const std::string &unknown_name)
 
typedef Gradient(* GradientFunctionPointer) (const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
 
typedef Tensor(* HessianFunctionPointer) (const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
 
typedef std::map< std::pair< const System *, unsigned int >, ErrorVector * > ErrorMap
 

Public Member Functions

 ExactErrorEstimator ()
 
 ~ExactErrorEstimator ()
 
void attach_exact_values (std::vector< FunctionBase< Number > * > f)
 
void attach_exact_value (unsigned int sys_num, FunctionBase< Number > *f)
 
void attach_exact_value (ValueFunctionPointer fptr)
 
void attach_exact_derivs (std::vector< FunctionBase< Gradient > * > g)
 
void attach_exact_deriv (unsigned int sys_num, FunctionBase< Gradient > *g)
 
void attach_exact_deriv (GradientFunctionPointer gptr)
 
void attach_exact_hessians (std::vector< FunctionBase< Tensor > * > h)
 
void attach_exact_hessian (unsigned int sys_num, FunctionBase< Tensor > *h)
 
void attach_exact_hessian (HessianFunctionPointer hptr)
 
void attach_reference_solution (EquationSystems *es_fine)
 
void extra_quadrature_order (const int extraorder)
 
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
 
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)
 
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)
 

Public Attributes

SystemNorm error_norm
 

Protected Member Functions

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

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
 
void clear_functors ()
 

Private Attributes

Number(* _exact_value )(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
 
Gradient(* _exact_deriv )(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
 
Tensor(* _exact_hessian )(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
 
std::vector< std::unique_ptr< FunctionBase< Number > > > _exact_values
 
std::vector< std::unique_ptr< FunctionBase< Gradient > > > _exact_derivs
 
std::vector< std::unique_ptr< FunctionBase< Tensor > > > _exact_hessians
 
EquationSystems_equation_systems_fine
 
int _extra_order
 

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.

typedef Gradient(* libMesh::ExactErrorEstimator::GradientFunctionPointer) (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.

Definition at line 129 of file exact_error_estimator.h.

typedef Tensor(* libMesh::ExactErrorEstimator::HessianFunctionPointer) (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.

Definition at line 152 of file exact_error_estimator.h.

typedef Number(* libMesh::ExactErrorEstimator::ValueFunctionPointer) (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.

Definition at line 106 of file exact_error_estimator.h.

Constructor & Destructor Documentation

libMesh::ExactErrorEstimator::ExactErrorEstimator ( )
inline

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)
const class libmesh_nullptr_t libmesh_nullptr
Tensor(* _exact_hessian)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
Gradient(* _exact_deriv)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
libMesh::ExactErrorEstimator::~ExactErrorEstimator ( )
inline

Destructor.

Definition at line 87 of file exact_error_estimator.h.

References attach_exact_value(), and attach_exact_values().

87 {}

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 105 of file exact_error_estimator.C.

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

107 {
108  if (_exact_derivs.size() <= sys_num)
109  _exact_derivs.resize(sys_num+1);
110 
111  if (g)
112  _exact_derivs[sys_num] = g->clone();
113 }
std::vector< std::unique_ptr< FunctionBase< Gradient > > > _exact_derivs
void libMesh::ExactErrorEstimator::attach_exact_deriv ( GradientFunctionPointer  gptr)

Definition at line 81 of file exact_error_estimator.C.

References _equation_systems_fine, _exact_deriv, clear_functors(), and libmesh_nullptr.

82 {
83  libmesh_assert(gptr);
84  _exact_deriv = gptr;
85 
86  // We're not using a fine grid solution
88 
89  // We're not using user-provided functors
90  this->clear_functors();
91 }
const class libmesh_nullptr_t libmesh_nullptr
Gradient(* _exact_deriv)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
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 94 of file exact_error_estimator.C.

References _exact_derivs, and libmesh_nullptr.

95 {
96  // Automatically delete any previous _exact_derivs entries, then add a new
97  // entry for each system.
98  _exact_derivs.clear();
99 
100  for (auto ptr : g)
101  _exact_derivs.emplace_back(ptr ? ptr->clone() : libmesh_nullptr);
102 }
const class libmesh_nullptr_t libmesh_nullptr
std::vector< std::unique_ptr< FunctionBase< Gradient > > > _exact_derivs
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 142 of file exact_error_estimator.C.

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

144 {
145  if (_exact_hessians.size() <= sys_num)
146  _exact_hessians.resize(sys_num+1);
147 
148  if (h)
149  _exact_hessians[sys_num] = h->clone();
150 }
std::vector< std::unique_ptr< FunctionBase< Tensor > > > _exact_hessians
void libMesh::ExactErrorEstimator::attach_exact_hessian ( HessianFunctionPointer  hptr)

Definition at line 118 of file exact_error_estimator.C.

References _equation_systems_fine, _exact_hessian, clear_functors(), and libmesh_nullptr.

119 {
120  libmesh_assert(hptr);
121  _exact_hessian = hptr;
122 
123  // We're not using a fine grid solution
125 
126  // We're not using user-provided functors
127  this->clear_functors();
128 }
const class libmesh_nullptr_t libmesh_nullptr
Tensor(* _exact_hessian)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
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 131 of file exact_error_estimator.C.

References _exact_hessians, and libmesh_nullptr.

132 {
133  // Automatically delete any previous _exact_hessians entries, then add a new
134  // entry for each system.
135  _exact_hessians.clear();
136 
137  for (auto ptr : h)
138  _exact_hessians.emplace_back(ptr ? ptr->clone() : libmesh_nullptr);
139 }
std::vector< std::unique_ptr< FunctionBase< Tensor > > > _exact_hessians
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 70 of file exact_error_estimator.C.

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

Referenced by ~ExactErrorEstimator().

72 {
73  if (_exact_values.size() <= sys_num)
74  _exact_values.resize(sys_num+1);
75 
76  if (f)
77  _exact_values[sys_num] = f->clone();
78 }
virtual std::unique_ptr< FunctionBase< Output > > clone() const =0
std::vector< std::unique_ptr< FunctionBase< Number > > > _exact_values
void libMesh::ExactErrorEstimator::attach_exact_value ( ValueFunctionPointer  fptr)

Definition at line 46 of file exact_error_estimator.C.

References _equation_systems_fine, _exact_value, clear_functors(), and libmesh_nullptr.

47 {
48  libmesh_assert(fptr);
49  _exact_value = fptr;
50 
51  // We're not using a fine grid solution
53 
54  // We're not using user-provided functors
55  this->clear_functors();
56 }
Number(* _exact_value)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
const class libmesh_nullptr_t libmesh_nullptr
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 59 of file exact_error_estimator.C.

References _exact_values, and libmesh_nullptr.

Referenced by ~ExactErrorEstimator().

60 {
61  // Automatically delete any previous _exact_values entries, then add a new
62  // entry for each system.
63  _exact_values.clear();
64 
65  for (auto ptr : f)
66  _exact_values.emplace_back(ptr ? ptr->clone() : libmesh_nullptr);
67 }
const class libmesh_nullptr_t libmesh_nullptr
std::vector< std::unique_ptr< FunctionBase< Number > > > _exact_values
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 153 of file exact_error_estimator.C.

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

154 {
155  libmesh_assert(es_fine);
156  _equation_systems_fine = es_fine;
157 
158  // If we're using a fine grid solution, we're not using exact value
159  // function pointers or functors.
163 
164  this->clear_functors();
165 }
Number(* _exact_value)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
const class libmesh_nullptr_t libmesh_nullptr
Tensor(* _exact_hessian)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
Gradient(* _exact_deriv)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
void libMesh::ExactErrorEstimator::clear_functors ( )
private

Helper method for cleanup

Definition at line 514 of file exact_error_estimator.C.

References _exact_derivs, _exact_hessians, and _exact_values.

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

515 {
516  _exact_values.clear();
517  _exact_derivs.clear();
518  _exact_hessians.clear();
519 }
std::vector< std::unique_ptr< FunctionBase< Tensor > > > _exact_hessians
std::vector< std::unique_ptr< FunctionBase< Gradient > > > _exact_derivs
std::vector< std::unique_ptr< FunctionBase< Number > > > _exact_values
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 167 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(), 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::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().

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

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 }
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
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.

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

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

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(), 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 }
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
void libMesh::ExactErrorEstimator::extra_quadrature_order ( const int  extraorder)
inline

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

Definition at line 170 of file exact_error_estimator.h.

References _extra_order, estimate_error(), and libmesh_nullptr.

171  { _extra_order = extraorder; }
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 381 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().

387 {
388  // The (string) name of this system
389  const std::string & sys_name = system.name();
390  const unsigned int sys_num = system.number();
391 
392  const unsigned int var = system.variable_number(var_name);
393  const unsigned int var_component =
394  system.variable_scalar_number(var, 0);
395 
396  const Parameters & parameters = system.get_equation_systems().parameters;
397 
398  // reinitialize the element-specific data
399  // for the current element
400  fe->reinit (elem);
401 
402  // Get the data we need to compute with
403  const std::vector<Real> & JxW = fe->get_JxW();
404  const std::vector<std::vector<Real>> & phi_values = fe->get_phi();
405  const std::vector<std::vector<RealGradient>> & dphi_values = fe->get_dphi();
406  const std::vector<Point> & q_point = fe->get_xyz();
407 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
408  const std::vector<std::vector<RealTensor>> & d2phi_values = fe->get_d2phi();
409 #endif
410 
411  // The number of shape functions
412  const unsigned int n_sf =
413  cast_int<unsigned int>(Uelem.size());
414 
415  // The number of quadrature points
416  const unsigned int n_qp =
417  cast_int<unsigned int>(JxW.size());
418 
419  Real error_val = 0;
420 
421  // Begin the loop over the Quadrature points.
422  //
423  for (unsigned int qp=0; qp<n_qp; qp++)
424  {
425  // Real u_h = 0.;
426  // RealGradient grad_u_h;
427 
428  Number u_h = 0.;
429 
430  Gradient grad_u_h;
431 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
432  Tensor grad2_u_h;
433 #endif
434 
435  // Compute solution values at the current
436  // quadrature point. This requires a sum
437  // over all the shape functions evaluated
438  // at the quadrature point.
439  for (unsigned int i=0; i<n_sf; i++)
440  {
441  // Values from current solution.
442  u_h += phi_values[i][qp]*Uelem(i);
443  grad_u_h += dphi_values[i][qp]*Uelem(i);
444 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
445  grad2_u_h += d2phi_values[i][qp]*Uelem(i);
446 #endif
447  }
448 
449  // Compute the value of the error at this quadrature point
450  if (error_norm.type(var) == L2 ||
451  error_norm.type(var) == H1 ||
452  error_norm.type(var) == H2)
453  {
454  Number val_error = u_h;
455  if (_exact_value)
456  val_error -= _exact_value(q_point[qp],parameters,sys_name,var_name);
457  else if (_exact_values.size() > sys_num && _exact_values[sys_num])
458  val_error -= _exact_values[sys_num]->
459  component(var_component, q_point[qp], system.time);
460  else if (_equation_systems_fine)
461  val_error -= (*fine_values)(q_point[qp]);
462 
463  // Add the squares of the error to each contribution
464  error_val += JxW[qp]*TensorTools::norm_sq(val_error);
465  }
466 
467  // Compute the value of the error in the gradient at this
468  // quadrature point
469  if (error_norm.type(var) == H1 ||
470  error_norm.type(var) == H1_SEMINORM ||
471  error_norm.type(var) == H2)
472  {
473  Gradient grad_error = grad_u_h;
474  if (_exact_deriv)
475  grad_error -= _exact_deriv(q_point[qp],parameters,sys_name,var_name);
476  else if (_exact_derivs.size() > sys_num && _exact_derivs[sys_num])
477  grad_error -= _exact_derivs[sys_num]->
478  component(var_component, q_point[qp], system.time);
479  else if (_equation_systems_fine)
480  grad_error -= fine_values->gradient(q_point[qp]);
481 
482  error_val += JxW[qp]*grad_error.norm_sq();
483  }
484 
485 
486 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
487  // Compute the value of the error in the hessian at this
488  // quadrature point
489  if ((error_norm.type(var) == H2_SEMINORM ||
490  error_norm.type(var) == H2))
491  {
492  Tensor grad2_error = grad2_u_h;
493  if (_exact_hessian)
494  grad2_error -= _exact_hessian(q_point[qp],parameters,sys_name,var_name);
495  else if (_exact_hessians.size() > sys_num && _exact_hessians[sys_num])
496  grad2_error -= _exact_hessians[sys_num]->
497  component(var_component, q_point[qp], system.time);
498  else if (_equation_systems_fine)
499  grad2_error -= fine_values->hessian(q_point[qp]);
500 
501  error_val += JxW[qp]*grad2_error.norm_sq();
502  }
503 #endif
504 
505  } // end qp loop
506 
507  libmesh_assert_greater_equal (error_val, 0.);
508 
509  return error_val;
510 }
Number(* _exact_value)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
std::vector< std::unique_ptr< FunctionBase< Tensor > > > _exact_hessians
virtual unsigned int size() const libmesh_override
Definition: dense_vector.h:87
std::vector< std::unique_ptr< FunctionBase< Gradient > > > _exact_derivs
NumberVectorValue Gradient
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Tensor(* _exact_hessian)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
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)
std::vector< std::unique_ptr< FunctionBase< Number > > > _exact_values
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
inlinevirtual
Returns
The type for the ErrorEstimator subclass.

Implements libMesh::ErrorEstimator.

Definition at line 192 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 246 of file exact_error_estimator.h.

Referenced by attach_exact_deriv(), attach_exact_hessian(), 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 210 of file exact_error_estimator.h.

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

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

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

Definition at line 234 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 219 of file exact_error_estimator.h.

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

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

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

Definition at line 240 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 201 of file exact_error_estimator.h.

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

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

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

Definition at line 228 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 266 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(), libMesh::PatchRecoveryErrorEstimator::PatchRecoveryErrorEstimator(), libMesh::JumpErrorEstimator::reinit_sides(), and libMesh::UniformRefinementEstimator::UniformRefinementEstimator().


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