libMesh
exact_error_estimator.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 // C++ includes
20 #include <algorithm> // for std::fill
21 #include <cstdlib> // *must* precede <cmath> for proper std:abs() on PGI, Sun Studio CC
22 #include <cmath> // for sqrt
23 
24 
25 // Local Includes
26 #include "libmesh/libmesh_common.h"
27 #include "libmesh/exact_error_estimator.h"
28 #include "libmesh/dof_map.h"
29 #include "libmesh/equation_systems.h"
30 #include "libmesh/error_vector.h"
31 #include "libmesh/fe_base.h"
32 #include "libmesh/libmesh_logging.h"
33 #include "libmesh/elem.h"
34 #include "libmesh/mesh_base.h"
35 #include "libmesh/mesh_function.h"
36 #include "libmesh/numeric_vector.h"
37 #include "libmesh/quadrature.h"
38 #include "libmesh/system.h"
39 #include "libmesh/tensor_tools.h"
40 
41 namespace libMesh
42 {
43 
44 //-----------------------------------------------------------------
45 // ErrorEstimator implementations
47  const Parameters & parameters,
48  const std::string & sys_name,
49  const std::string & unknown_name))
50 {
53 
54  // We're not using a fine grid solution
56 
57  // We're not using user-provided functors
58  this->clear_functors();
59 }
60 
61 
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 }
79 
80 
81 void ExactErrorEstimator::attach_exact_value (unsigned int sys_num,
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 }
90 
91 
93  const Parameters & parameters,
94  const std::string & sys_name,
95  const std::string & unknown_name))
96 {
99 
100  // We're not using a fine grid solution
102 
103  // We're not using user-provided functors
104  this->clear_functors();
105 }
106 
107 
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 }
125 
126 
127 void ExactErrorEstimator::attach_exact_deriv (unsigned int sys_num,
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 }
136 
137 
138 
139 
141  const Parameters & parameters,
142  const std::string & sys_name,
143  const std::string & unknown_name))
144 {
145  libmesh_assert(hptr);
146  _exact_hessian = hptr;
147 
148  // We're not using a fine grid solution
150 
151  // We're not using user-provided functors
152  this->clear_functors();
153 }
154 
155 
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 }
173 
174 
175 void ExactErrorEstimator::attach_exact_hessian (unsigned int sys_num,
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 }
184 
185 
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 }
199 
201  ErrorVector & error_per_cell,
202  const NumericVector<Number> * solution_vector,
203  bool estimate_parent_error)
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;
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>
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;
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 }
411 
412 
413 
415  const std::string & var_name,
416  const Elem * elem,
417  const DenseVector<Number> & Uelem,
418  FEBase * fe,
419  MeshFunction * fine_values) const
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 }
544 
545 
546 
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 }
561 
562 
563 
564 } // namespace libMesh
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:178
Real time
For time-dependent problems, this is the time t at the beginning of the current timestep.
Definition: system.h:1545
void attach_exact_derivs(std::vector< FunctionBase< Gradient > * > g)
Clone and attach arbitrary functors which compute the exact gradients of the EquationSystems&#39; solutio...
This is the EquationSystems class.
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_n...
const std::string & variable_name(const unsigned int i) const
Definition: system.h:2134
This class provides the ability to map between arbitrary, user-defined strings and several data types...
Definition: parameters.h:63
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.
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:1697
Real norm_sq() const
Definition: type_tensor.h:1270
int _extra_order
Extra order to use for quadrature rule.
unsigned int dim
The ErrorVector is a specialization of the StatisticsVector for error data computed on a finite eleme...
Definition: error_vector.h:50
ImplicitSystem & sys
void attach_exact_hessians(std::vector< FunctionBase< Tensor > * > h)
Clone and attach arbitrary functors which compute the exact second derivatives of the EquationSystems...
Real weight(unsigned int var) const
Definition: system_norm.h:292
const Elem * parent() const
Definition: elem.h:2346
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...
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
MeshBase & mesh
UniquePtr< NumericVector< Number > > current_local_solution
All the values I need to compute my contribution to the simulation at hand.
Definition: system.h:1535
const class libmesh_nullptr_t libmesh_nullptr
Gradient gradient(const Point &p, const Real time=0.)
This class defines a vector in LIBMESH_DIM dimensional Real or Complex space.
const unsigned int n_vars
Definition: tecplot_io.C:68
The libMesh namespace provides an interface to certain functionality in the library.
virtual unsigned int size() const libmesh_override
Definition: dense_vector.h:87
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 a...
DIE A HORRIBLE DEATH HERE typedef float ErrorVectorReal
const std::string & name() const
Definition: system.h:1998
This is the MeshBase class.
Definition: mesh_base.h:68
SimpleRange< ChildRefIter > child_ref_range()
Returns a range with all children of a parent element, usable in range-based for loops.
Definition: elem.h:1699
Number fptr(const Point &p, const Parameters &, const std::string &libmesh_dbg_var(sys_name), const std::string &unknown_name)
Definition: projection.C:78
libmesh_assert(j)
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
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
virtual SimpleRange< element_iterator > active_local_element_ptr_range()=0
This class handles the numbering of degrees of freedom on a mesh.
Definition: dof_map.h:167
unsigned short int variable_number(const std::string &var) const
Definition: system.C:1263
virtual dof_id_type max_elem_id() const =0
void update_global_solution(std::vector< Number > &global_soln) const
Fill the input vector global_soln so that it contains the global solution on all processors.
Definition: system.C:662
const MeshBase & get_mesh() const
Definition: system.h:2014
const DofMap & get_dof_map() const
Definition: system.h:2030
const std::vector< Real > & get_JxW() const
Definition: fe_abstract.h:237
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
const std::vector< std::vector< OutputShape > > & get_phi() const
Definition: fe_base.h:208
UniquePtr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition: system.h:1523
Real norm_sq() const
Definition: type_vector.h:940
void clear_functors()
Helper method for cleanup.
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 &)
const std::vector< std::vector< OutputGradient > > & get_dphi() const
Definition: fe_base.h:216
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Number current_solution(const dof_id_type global_dof_number) const
Definition: system.C:192
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...
virtual void swap(NumericVector< T > &v)
Swaps the contents of this with v.
unsigned int number() const
Definition: system.h:2006
const std::vector< std::vector< OutputTensor > > & get_d2phi() const
Definition: fe_base.h:290
UniquePtr< QBase > default_quadrature_rule(const unsigned int dim, const int extraorder=0) const
Definition: fe_type.C:30
const Parallel::Communicator & comm() const
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
Gradient gptr(const Point &p, const Parameters &, const std::string &libmesh_dbg_var(sys_name), const std::string &unknown_name)
Definition: projection.C:93
void attach_exact_values(std::vector< FunctionBase< Number > * > f)
Clone and attach arbitrary functors which compute the exact values of the EquationSystems&#39; solutions ...
virtual UniquePtr< FunctionBase< Output > > clone() const =0
const EquationSystems & get_equation_systems() const
Definition: system.h:712
Parameters parameters
Data structure holding arbitrary parameters.
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.
unsigned int variable_scalar_number(const std::string &var, unsigned int component) const
Definition: system.h:2145
const T_sys & get_system(const std::string &name) const
unsigned int n_vars() const
Definition: system.h:2086
FEMNormType type(unsigned int var) const
Definition: system_norm.h:268
This class provides function-like objects for data distributed over a mesh.
Definition: mesh_function.h:53
void attach_reference_solution(EquationSystems *es_fine)
Attach function similar to system.h which allows the user to attach a second EquationSystems object w...
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 solutio...
dof_id_type id() const
Definition: dof_object.h:632
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...
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
Tensor hessian(const Point &p, const Real time=0.)
const std::vector< Point > & get_xyz() const
Definition: fe_abstract.h:230
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.
This class forms the foundation from which generic finite elements may be derived.
This class defines a tensor in LIBMESH_DIM dimensional Real or Complex space.
virtual void reinit(const Elem *elem, const std::vector< Point > *const pts=libmesh_nullptr, const std::vector< Real > *const weights=libmesh_nullptr)=0
This is at the core of this class.
uint8_t dof_id_type
Definition: id_types.h:64
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Fills the vector di with the global degree of freedom indices for the element.
Definition: dof_map.C:1917
static UniquePtr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.