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 (const ExactErrorEstimator &)=delete
 
ExactErrorEstimatoroperator= (const ExactErrorEstimator &)=delete
 
 ExactErrorEstimator (ExactErrorEstimator &&)=default
 
ExactErrorEstimatoroperator= (ExactErrorEstimator &&)=default
 
virtual ~ExactErrorEstimator ()=default
 
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=nullptr, bool estimate_parent_error=false) override
 
virtual ErrorEstimatorType type () const 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=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=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) 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

ValueFunctionPointer _exact_value
 
GradientFunctionPointer _exact_deriv
 
HessianFunctionPointer _exact_hessian
 
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 57 of file exact_error_estimator.h.

Member Typedef Documentation

◆ ErrorMap

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

◆ GradientFunctionPointer

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 122 of file exact_error_estimator.h.

◆ HessianFunctionPointer

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 145 of file exact_error_estimator.h.

◆ ValueFunctionPointer

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 99 of file exact_error_estimator.h.

Constructor & Destructor Documentation

◆ ExactErrorEstimator() [1/3]

libMesh::ExactErrorEstimator::ExactErrorEstimator ( )

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

Definition at line 48 of file exact_error_estimator.C.

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

48  :
50  _exact_value(nullptr),
51  _exact_deriv(nullptr),
52  _exact_hessian(nullptr),
53  _equation_systems_fine(nullptr),
54  _extra_order(0)
55 {
56  error_norm = H1;
57 }
HessianFunctionPointer _exact_hessian
GradientFunctionPointer _exact_deriv

◆ ExactErrorEstimator() [2/3]

libMesh::ExactErrorEstimator::ExactErrorEstimator ( const ExactErrorEstimator )
delete

This class cannot be (default) copy constructed/assigned because it has containers of unique_ptrs. Explicitly deleting these functions is the best way to document this fact.

◆ ExactErrorEstimator() [3/3]

libMesh::ExactErrorEstimator::ExactErrorEstimator ( ExactErrorEstimator &&  )
default

Defaulted move ctor, move assignment operator, and destructor.

◆ ~ExactErrorEstimator()

virtual libMesh::ExactErrorEstimator::~ExactErrorEstimator ( )
virtualdefault

Member Function Documentation

◆ attach_exact_deriv() [1/2]

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

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

127 {
128  if (_exact_derivs.size() <= sys_num)
129  _exact_derivs.resize(sys_num+1);
130 
131  if (g)
132  _exact_derivs[sys_num] = g->clone();
133 }
std::vector< std::unique_ptr< FunctionBase< Gradient > > > _exact_derivs

◆ attach_exact_deriv() [2/2]

void libMesh::ExactErrorEstimator::attach_exact_deriv ( GradientFunctionPointer  gptr)

Definition at line 101 of file exact_error_estimator.C.

References _equation_systems_fine, _exact_deriv, and clear_functors().

102 {
103  libmesh_assert(gptr);
104  _exact_deriv = gptr;
105 
106  // We're not using a fine grid solution
107  _equation_systems_fine = nullptr;
108 
109  // We're not using user-provided functors
110  this->clear_functors();
111 }
GradientFunctionPointer _exact_deriv

◆ attach_exact_derivs()

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

References _exact_derivs.

115 {
116  // Automatically delete any previous _exact_derivs entries, then add a new
117  // entry for each system.
118  _exact_derivs.clear();
119 
120  for (auto ptr : g)
121  _exact_derivs.emplace_back(ptr ? ptr->clone() : nullptr);
122 }
std::vector< std::unique_ptr< FunctionBase< Gradient > > > _exact_derivs

◆ attach_exact_hessian() [1/2]

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

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

164 {
165  if (_exact_hessians.size() <= sys_num)
166  _exact_hessians.resize(sys_num+1);
167 
168  if (h)
169  _exact_hessians[sys_num] = h->clone();
170 }
std::vector< std::unique_ptr< FunctionBase< Tensor > > > _exact_hessians

◆ attach_exact_hessian() [2/2]

void libMesh::ExactErrorEstimator::attach_exact_hessian ( HessianFunctionPointer  hptr)

Definition at line 138 of file exact_error_estimator.C.

References _equation_systems_fine, _exact_hessian, and clear_functors().

139 {
140  libmesh_assert(hptr);
141  _exact_hessian = hptr;
142 
143  // We're not using a fine grid solution
144  _equation_systems_fine = nullptr;
145 
146  // We're not using user-provided functors
147  this->clear_functors();
148 }
HessianFunctionPointer _exact_hessian

◆ attach_exact_hessians()

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

References _exact_hessians.

152 {
153  // Automatically delete any previous _exact_hessians entries, then add a new
154  // entry for each system.
155  _exact_hessians.clear();
156 
157  for (auto ptr : h)
158  _exact_hessians.emplace_back(ptr ? ptr->clone() : nullptr);
159 }
std::vector< std::unique_ptr< FunctionBase< Tensor > > > _exact_hessians

◆ attach_exact_value() [1/2]

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

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

92 {
93  if (_exact_values.size() <= sys_num)
94  _exact_values.resize(sys_num+1);
95 
96  if (f)
97  _exact_values[sys_num] = f->clone();
98 }
virtual std::unique_ptr< FunctionBase< Output > > clone() const =0
std::vector< std::unique_ptr< FunctionBase< Number > > > _exact_values

◆ attach_exact_value() [2/2]

void libMesh::ExactErrorEstimator::attach_exact_value ( ValueFunctionPointer  fptr)

Definition at line 66 of file exact_error_estimator.C.

References _equation_systems_fine, _exact_value, and clear_functors().

67 {
68  libmesh_assert(fptr);
69  _exact_value = fptr;
70 
71  // We're not using a fine grid solution
72  _equation_systems_fine = nullptr;
73 
74  // We're not using user-provided functors
75  this->clear_functors();
76 }

◆ attach_exact_values()

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

References _exact_values.

80 {
81  // Automatically delete any previous _exact_values entries, then add a new
82  // entry for each system.
83  _exact_values.clear();
84 
85  for (auto ptr : f)
86  _exact_values.emplace_back(ptr ? ptr->clone() : nullptr);
87 }
std::vector< std::unique_ptr< FunctionBase< Number > > > _exact_values

◆ attach_reference_solution()

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

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

174 {
175  libmesh_assert(es_fine);
176  _equation_systems_fine = es_fine;
177 
178  // If we're using a fine grid solution, we're not using exact value
179  // function pointers or functors.
180  _exact_value = nullptr;
181  _exact_deriv = nullptr;
182  _exact_hessian = nullptr;
183 
184  this->clear_functors();
185 }
HessianFunctionPointer _exact_hessian
GradientFunctionPointer _exact_deriv

◆ clear_functors()

void libMesh::ExactErrorEstimator::clear_functors ( )
private

Helper method for cleanup

Definition at line 534 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().

535 {
536  _exact_values.clear();
537  _exact_derivs.clear();
538  _exact_hessians.clear();
539 }
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

◆ estimate_error()

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

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

References _equation_systems_fine, _exact_derivs, _exact_hessians, _exact_values, _extra_order, 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(), mesh, 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().

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

◆ 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(), 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 (unsigned int s = 0; s != equation_systems.n_systems(); ++s)
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 }
IntRange< std::size_t > index_range(const std::vector< T > &vec)
Definition: int_range.h:104
virtual void estimate_error(const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=nullptr, bool estimate_parent_error=false)=0

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

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 93 of file error_estimator.C.

References libMesh::ErrorEstimator::error_norm, libMesh::ErrorEstimator::estimate_error(), libMesh::EquationSystems::get_system(), 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 (unsigned int s = 0; s != equation_systems.n_systems(); ++s)
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.find(std::make_pair(&sys, v)) ==
111  errors_per_cell.end())
112  continue;
113 
114  // Calculate error in only one variable
115  std::vector<Real> weights(n_vars, 0.0);
116  weights[v] = 1.0;
117  this->error_norm =
118  SystemNorm(std::vector<FEMNormType>(n_vars, old_error_norm.type(v)),
119  weights);
120 
121  const NumericVector<Number> * solution_vector = nullptr;
122  if (solution_vectors &&
123  solution_vectors->find(&sys) != solution_vectors->end())
124  solution_vector = solution_vectors->find(&sys)->second;
125 
126  this->estimate_error
127  (sys, *errors_per_cell[std::make_pair(&sys, v)],
128  solution_vector, estimate_parent_error);
129  }
130  }
131 
132  // Restore our old state before returning
133  this->error_norm = old_error_norm;
134 }
virtual void estimate_error(const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=nullptr, bool estimate_parent_error=false)=0
const unsigned int n_vars
Definition: tecplot_io.C:69

◆ extra_quadrature_order()

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 163 of file exact_error_estimator.h.

References _extra_order.

164  { _extra_order = extraorder; }

◆ find_squared_element_error()

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 401 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().

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

◆ operator=() [1/2]

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

◆ operator=() [2/2]

ExactErrorEstimator& libMesh::ExactErrorEstimator::operator= ( ExactErrorEstimator &&  )
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 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().

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 }

◆ type()

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

Implements libMesh::ErrorEstimator.

Definition at line 60 of file exact_error_estimator.C.

References libMesh::EXACT.

61 {
62  return EXACT;
63 }

Member Data Documentation

◆ _equation_systems_fine

EquationSystems* libMesh::ExactErrorEstimator::_equation_systems_fine
private

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

Definition at line 229 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().

◆ _exact_deriv

GradientFunctionPointer libMesh::ExactErrorEstimator::_exact_deriv
private

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

Definition at line 199 of file exact_error_estimator.h.

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

◆ _exact_derivs

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 217 of file exact_error_estimator.h.

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

◆ _exact_hessian

HessianFunctionPointer libMesh::ExactErrorEstimator::_exact_hessian
private

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

Definition at line 205 of file exact_error_estimator.h.

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

◆ _exact_hessians

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 223 of file exact_error_estimator.h.

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

◆ _exact_value

ValueFunctionPointer libMesh::ExactErrorEstimator::_exact_value
private

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

Definition at line 193 of file exact_error_estimator.h.

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

◆ _exact_values

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 211 of file exact_error_estimator.h.

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

◆ _extra_order

int libMesh::ExactErrorEstimator::_extra_order
private

Extra order to use for quadrature rule

Definition at line 249 of file exact_error_estimator.h.

Referenced by estimate_error(), and extra_quadrature_order().

◆ 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 161 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::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::PatchRecoveryErrorEstimator(), libMesh::JumpErrorEstimator::reinit_sides(), and libMesh::UniformRefinementEstimator::UniformRefinementEstimator().


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