libMesh
exact_error_estimator.h
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 
20 #ifndef LIBMESH_EXACT_ERROR_ESTIMATOR_H
21 #define LIBMESH_EXACT_ERROR_ESTIMATOR_H
22 
23 // Local Includes
24 #include "libmesh/error_estimator.h"
25 #include "libmesh/function_base.h"
26 
27 // C++ includes
28 #include <cstddef>
29 #include <string>
30 #include <vector>
31 
32 namespace libMesh
33 {
34 
35 // Forward Declarations
36 class Elem;
37 template <typename T> class FEGenericBase;
39 class MeshFunction;
40 class Point;
41 class Parameters;
42 
43 template <typename T> class DenseVector;
44 
45 // Is there any way to simplify this?
46 // All we need are Tensor and Gradient. - RHS
47 template <typename T> class TensorValue;
48 template <typename T> class VectorValue;
50 typedef NumberTensorValue Tensor;
52 typedef NumberVectorValue Gradient;
53 
54 
55 
68 {
69 public:
70 
81  _extra_order(0)
82  { error_norm = H1; }
83 
88 
93  void attach_exact_values (std::vector<FunctionBase<Number> *> f);
94 
99  void attach_exact_value (unsigned int sys_num,
101 
106  void attach_exact_value ( Number fptr(const Point & p,
107  const Parameters & Parameters,
108  const std::string & sys_name,
109  const std::string & unknown_name));
110 
115  void attach_exact_derivs (std::vector<FunctionBase<Gradient> *> g);
116 
121  void attach_exact_deriv (unsigned int sys_num,
123 
128  void attach_exact_deriv ( Gradient gptr(const Point & p,
129  const Parameters & parameters,
130  const std::string & sys_name,
131  const std::string & unknown_name));
132 
137  void attach_exact_hessians (std::vector<FunctionBase<Tensor> *> h);
138 
143  void attach_exact_hessian (unsigned int sys_num,
145 
150  void attach_exact_hessian ( Tensor hptr(const Point & p,
151  const Parameters & parameters,
152  const std::string & sys_name,
153  const std::string & unknown_name));
154 
161 
162 
167  void extra_quadrature_order (const int extraorder)
168  { _extra_order = extraorder; }
169 
170 
171  // Bring the base class functionality into the name lookup
172  // procedure. This allows for alternative calling formats
173  // defined in the base class. Thanks Wolfgang.
174  // GCC 2.95.3 cannot compile such code. Since it was not really
175  // essential to the functioning of this class, it's been removed.
176  // using ErrorEstimator::estimate_error;
177 
184  virtual void estimate_error (const System & system,
185  ErrorVector & error_per_cell,
186  const NumericVector<Number> * solution_vector = libmesh_nullptr,
187  bool estimate_parent_error = false) libmesh_override;
188 
189  virtual ErrorEstimatorType type() const libmesh_override
190  { return EXACT;}
191 
192 private:
193 
198  Number (* _exact_value) (const Point & p,
199  const Parameters & parameters,
200  const std::string & sys_name,
201  const std::string & unknown_name);
202 
207  Gradient (* _exact_deriv) (const Point & p,
208  const Parameters & parameters,
209  const std::string & sys_name,
210  const std::string & unknown_name);
211 
216  Tensor (* _exact_hessian) (const Point & p,
217  const Parameters & parameters,
218  const std::string & sys_name,
219  const std::string & unknown_name);
220 
225  std::vector<FunctionBase<Number> *> _exact_values;
226 
231  std::vector<FunctionBase<Gradient> *> _exact_derivs;
232 
237  std::vector<FunctionBase<Tensor> *> _exact_hessians;
238 
244 
248  Real find_squared_element_error (const System & system,
249  const std::string & var_name,
250  const Elem * elem,
251  const DenseVector<Number> & Uelem,
252  FEBase * fe,
253  MeshFunction * fine_values) const;
254 
258  void clear_functors ();
259 
264 };
265 
266 
267 } // namespace libMesh
268 
269 #endif // LIBMESH_EXACT_ERROR_ESTIMATOR_H
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...
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.
int _extra_order
Extra order to use for quadrature rule.
The ErrorVector is a specialization of the StatisticsVector for error data computed on a finite eleme...
Definition: error_vector.h:50
void attach_exact_hessians(std::vector< FunctionBase< Tensor > * > h)
Clone and attach arbitrary functors which compute the exact second derivatives of the EquationSystems...
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
const class libmesh_nullptr_t libmesh_nullptr
This class defines a vector in LIBMESH_DIM dimensional Real or Complex space.
void extra_quadrature_order(const int extraorder)
Increases or decreases the order of the quadrature rule used for numerical integration.
The libMesh namespace provides an interface to certain functionality in the library.
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...
Number fptr(const Point &p, const Parameters &, const std::string &libmesh_dbg_var(sys_name), const std::string &unknown_name)
Definition: projection.C:78
TensorValue< Number > NumberTensorValue
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
NumberVectorValue Gradient
This class holds functions that will estimate the error in a finite element solution on a given mesh...
virtual ErrorEstimatorType type() const libmesh_override
void clear_functors()
Helper method for cleanup.
FEGenericBase< Real > FEBase
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
EquationSystems * _equation_systems_fine
Constant pointer to the EquationSystems object containing a fine grid solution.
Tensor(* _exact_hessian)(const Point &p, const Parameters &parameters, const std::string &sys_name, const std::string &unknown_name)
Function pointer to user-provided function which computes the exact hessian of the solution...
NumberTensorValue Tensor
This class implements an "error estimator" based on the difference between the approximate and exact ...
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 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.
Defines a dense vector for use in Finite Element-type computations.
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...
VectorValue< Number > NumberVectorValue
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
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.