libMesh
exact_solution.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2024 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 #ifndef LIBMESH_EXACT_SOLUTION_H
19 #define LIBMESH_EXACT_SOLUTION_H
20 
21 
22 // Local Includes
23 #include "libmesh/libmesh_common.h" // for Number
24 
25 // C++ includes
26 #include <map>
27 #include <vector>
28 #include <memory>
29 #include <set>
30 
31 namespace libMesh
32 {
33 
34 
35 // Forward Declarations
36 class Point;
37 class EquationSystems;
38 class Parameters;
39 class Mesh;
40 template <typename Output> class FunctionBase;
41 enum FEMNormType : int;
42 
43 // Is there any way to simplify this?
44 // All we need are Tensor and Gradient. - RHS
45 template <typename T> class TensorValue;
46 template <typename T> class VectorValue;
51 
66 {
67 
68 public:
74  explicit
75  ExactSolution (const EquationSystems & es);
76 
83  ExactSolution(const ExactSolution &) = delete;
84  ExactSolution & operator= (const ExactSolution &) = delete;
85  ExactSolution & operator= (ExactSolution &&) = delete;
86 
93 
100  void set_excluded_subdomains(const std::set<subdomain_id_type> & excluded);
101 
107  void attach_reference_solution (const EquationSystems * es_fine);
108 
113  void attach_exact_values (const std::vector<FunctionBase<Number> *> & f);
114 
119  void attach_exact_value (unsigned int sys_num,
121 
126  typedef Number (*ValueFunctionPointer)(const Point & p,
127  const Parameters & Parameters,
128  const std::string & sys_name,
129  const std::string & unknown_name);
131 
136  void attach_exact_derivs (const std::vector<FunctionBase<Gradient> *> & g);
137 
142  void attach_exact_deriv (unsigned int sys_num,
144 
149  typedef Gradient (*GradientFunctionPointer)(const Point & p,
150  const Parameters & parameters,
151  const std::string & sys_name,
152  const std::string & unknown_name);
154 
159  void attach_exact_hessians (std::vector<FunctionBase<Tensor> *> h);
160 
165  void attach_exact_hessian (unsigned int sys_num,
167 
172  typedef Tensor (*HessianFunctionPointer)(const Point & p,
173  const Parameters & parameters,
174  const std::string & sys_name,
175  const std::string & unknown_name);
177 
185  void extra_quadrature_order (const int extraorder)
186  { _extra_order = extraorder; }
187 
195  void compute_error(std::string_view sys_name,
196  std::string_view unknown_name);
197 
205  Real l2_error(std::string_view sys_name,
206  std::string_view unknown_name);
207 
215  Real l1_error(std::string_view sys_name,
216  std::string_view unknown_name);
217 
230  Real l_inf_error(std::string_view sys_name,
231  std::string_view unknown_name);
232 
240  Real h1_error(std::string_view sys_name,
241  std::string_view unknown_name);
242 
253  Real hcurl_error(std::string_view sys_name,
254  std::string_view unknown_name);
255 
266  Real hdiv_error(std::string_view sys_name,
267  std::string_view unknown_name);
268 
276  Real h2_error(std::string_view sys_name,
277  std::string_view unknown_name);
278 
289  Real error_norm(std::string_view sys_name,
290  std::string_view unknown_name,
291  const FEMNormType & norm);
292 private:
293 
300  template<typename OutputShape>
301  void _compute_error(std::string_view sys_name,
302  std::string_view unknown_name,
303  std::vector<Real> & error_vals);
304 
311  std::vector<Real> & _check_inputs(std::string_view sys_name,
312  std::string_view unknown_name);
313 
318  std::vector<std::unique_ptr<FunctionBase<Number>>> _exact_values;
319 
324  std::vector<std::unique_ptr<FunctionBase<Gradient>>> _exact_derivs;
325 
330  std::vector<std::unique_ptr<FunctionBase<Tensor>>> _exact_hessians;
331 
340  typedef std::map<std::string, std::vector<Real>, std::less<>>
342 
349  std::map<std::string, SystemErrorMap, std::less<>> _errors;
350 
356 
362 
367 
372  std::set<subdomain_id_type> _excluded_subdomains;
373 };
374 
375 
376 
377 } // namespace libMesh
378 
379 
380 #endif // LIBMESH_EXACT_SOLUTION_H
std::map< std::string, std::vector< Real >, std::less<> > SystemErrorMap
Data structure which stores the errors: ||e|| = ||u - u_h|| ||grad(e)|| = ||grad(u) - grad(u_h)|| for...
This class handles the computation of the L2 and/or H1 error for the Systems in the EquationSystems o...
Real l_inf_error(std::string_view sys_name, std::string_view unknown_name)
This is the EquationSystems class.
Real h2_error(std::string_view sys_name, std::string_view unknown_name)
Real l1_error(std::string_view sys_name, std::string_view unknown_name)
Real hdiv_error(std::string_view sys_name, std::string_view unknown_name)
This class provides the ability to map between arbitrary, user-defined strings and several data types...
Definition: parameters.h:67
const EquationSystems & _equation_systems
Constant reference to the EquationSystems object used for the simulation.
std::map< std::string, SystemErrorMap, std::less<> > _errors
A map of SystemErrorMaps, which contains entries for each system in the EquationSystems object...
ExactSolution(const EquationSystems &es)
Constructor.
void attach_exact_values(const std::vector< FunctionBase< Number > *> &f)
Clone and attach arbitrary functors which compute the exact values of the EquationSystems&#39; solutions ...
int _extra_order
Extra order to use for quadrature rule.
FEMNormType
defines an enum for norms defined on vectors of finite element coefficients
Real hcurl_error(std::string_view sys_name, std::string_view unknown_name)
void extra_quadrature_order(const int extraorder)
Increases or decreases the order of the quadrature rule used for numerical integration.
std::vector< Real > & _check_inputs(std::string_view sys_name, std::string_view unknown_name)
This function is responsible for checking the validity of the sys_name and unknown_name inputs...
This class defines a vector in LIBMESH_DIM dimensional Real or Complex space.
The libMesh namespace provides an interface to certain functionality in the library.
void attach_exact_derivs(const std::vector< FunctionBase< Gradient > *> &g)
Clone and attach arbitrary functors which compute the exact gradients of the EquationSystems&#39; solutio...
const EquationSystems * _equation_systems_fine
Constant pointer to the EquationSystems object containing the fine grid solution. ...
Number fptr(const Point &p, const Parameters &, const std::string &libmesh_dbg_var(sys_name), const std::string &unknown_name)
Definition: projection.C:80
Number(* 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...
std::set< subdomain_id_type > _excluded_subdomains
Elements in a subdomain from this set are skipped during the error computation.
TensorValue< Number > NumberTensorValue
void attach_exact_hessians(std::vector< FunctionBase< Tensor > *> h)
Clone and attach arbitrary functors which compute the exact second derivatives of the EquationSystems...
ExactSolution & operator=(const ExactSolution &)=delete
void compute_error(std::string_view sys_name, std::string_view unknown_name)
Computes and stores the error in the solution value e = u-u_h, the gradient grad(e) = grad(u) - grad(...
NumberVectorValue Gradient
std::vector< std::unique_ptr< FunctionBase< Gradient > > > _exact_derivs
User-provided functors which compute the exact derivative of the solution for each system...
Real l2_error(std::string_view sys_name, std::string_view unknown_name)
auto norm(const T &a) -> decltype(std::abs(a))
Definition: tensor_tools.h:74
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...
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...
std::vector< std::unique_ptr< FunctionBase< Tensor > > > _exact_hessians
User-provided functors which compute the exact hessians of the solution for each system.
void _compute_error(std::string_view sys_name, std::string_view unknown_name, std::vector< Real > &error_vals)
This function computes the error (in the solution and its first derivative) for a single unknown in a...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
NumberTensorValue Tensor
void set_excluded_subdomains(const std::set< subdomain_id_type > &excluded)
The user can indicate that elements in certain subdomains should be excluded from the error calculati...
Tensor(* 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...
Gradient gptr(const Point &p, const Parameters &, const std::string &libmesh_dbg_var(sys_name), const std::string &unknown_name)
Definition: projection.C:95
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...
void attach_reference_solution(const EquationSystems *es_fine)
Attach function similar to system.h which allows the user to attach a second EquationSystems object w...
VectorValue< Number > NumberVectorValue
Gradient(* 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...
Real h1_error(std::string_view sys_name, std::string_view unknown_name)
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:39
std::vector< std::unique_ptr< FunctionBase< Number > > > _exact_values
User-provided functors which compute the exact value of the solution for each system.
void ErrorVector unsigned int
Definition: adjoints_ex3.C:360
This class defines a tensor in LIBMESH_DIM dimensional Real or Complex space.
Real error_norm(std::string_view sys_name, std::string_view unknown_name, const FEMNormType &norm)