mesh_function.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2018 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_MESH_FUNCTION_H
21 #define LIBMESH_MESH_FUNCTION_H
22 
23 // Local Includes
24 #include "libmesh/function_base.h"
25 #include "libmesh/dense_vector.h"
26 #include "libmesh/vector_value.h"
27 #include "libmesh/tensor_value.h"
28 #include "libmesh/tree_base.h"
30 
31 // C++ includes
32 #include <cstddef>
33 #include <vector>
34 
35 namespace libMesh
36 {
37 
38 
39 // Forward Declarations
40 template <typename T> class DenseVector;
41 class EquationSystems;
42 template <typename T> class NumericVector;
43 class DofMap;
44 class PointLocatorBase;
45 
53 class MeshFunction : public FunctionBase<Number>,
54  public ParallelObject
55 {
56 public:
57 
66  MeshFunction (const EquationSystems & eqn_systems,
67  const NumericVector<Number> & vec,
68  const DofMap & dof_map,
69  const std::vector<unsigned int> & vars,
70  const FunctionBase<Number> * master=nullptr);
71 
80  MeshFunction (const EquationSystems & eqn_systems,
81  const NumericVector<Number> & vec,
82  const DofMap & dof_map,
83  const unsigned int var,
84  const FunctionBase<Number> * master=nullptr);
85 
91  MeshFunction (MeshFunction &&) = delete;
92  MeshFunction (const MeshFunction &) = delete;
93 
97  MeshFunction & operator= (const MeshFunction &) = delete;
98  MeshFunction & operator= (MeshFunction &&) = delete;
99 
103  ~MeshFunction ();
104 
110  virtual void init () override { this->init(Trees::NODES); }
111 
116  void init (const Trees::BuildType point_locator_build_type);
117 
121  virtual void clear () override;
122 
132  virtual std::unique_ptr<FunctionBase<Number>> clone () const override;
133 
138  Number operator() (const Point & p,
139  const Real time=0.) override;
140 
149  std::map<const Elem *, Number> discontinuous_value (const Point & p,
150  const Real time=0.);
151 
156  Gradient gradient (const Point & p,
157  const Real time=0.);
158 
165  std::map<const Elem *, Gradient> discontinuous_gradient (const Point & p,
166  const Real time=0.);
167 
168 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
169 
173  Tensor hessian (const Point & p,
174  const Real time=0.);
175 #endif
176 
183  void operator() (const Point & p,
184  const Real time,
185  DenseVector<Number> & output) override;
186 
192  void operator() (const Point & p,
193  const Real time,
194  DenseVector<Number> & output,
195  const std::set<subdomain_id_type> * subdomain_ids);
196 
202  void discontinuous_value (const Point & p,
203  const Real time,
204  std::map<const Elem *, DenseVector<Number>> & output);
205 
211  void discontinuous_value (const Point & p,
212  const Real time,
213  std::map<const Elem *, DenseVector<Number>> & output,
214  const std::set<subdomain_id_type> * subdomain_ids);
215 
222  void gradient (const Point & p,
223  const Real time,
224  std::vector<Gradient> & output,
225  const std::set<subdomain_id_type> * subdomain_ids = nullptr);
226 
232  void discontinuous_gradient (const Point & p,
233  const Real time,
234  std::map<const Elem *, std::vector<Gradient>> & output);
235 
241  void discontinuous_gradient (const Point & p,
242  const Real time,
243  std::map<const Elem *, std::vector<Gradient>> & output,
244  const std::set<subdomain_id_type> * subdomain_ids);
245 
252  void hessian (const Point & p,
253  const Real time,
254  std::vector<Tensor> & output,
255  const std::set<subdomain_id_type> * subdomain_ids = nullptr);
256 
263  const PointLocatorBase & get_point_locator (void) const;
264 
276 
287  void enable_out_of_mesh_mode(const Number & value);
288 
292  void disable_out_of_mesh_mode(void);
293 
301 
306 
307 protected:
308 
312  const Elem * find_element(const Point & p,
313  const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
314 
321  std::set<const Elem *> find_elements(const Point & p,
322  const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
323 
329 
335 
339  const DofMap & _dof_map;
340 
345  const std::vector<unsigned int> _system_vars;
346 
352 
358 
364 };
365 
366 
367 
368 
369 // ------------------------------------------------------------
370 // MeshFunction inline methods
371 
372 
373 } // namespace libMesh
374 
375 
376 #endif // LIBMESH_MESH_FUNCTION_H
Manages multiples systems of equations.
const EquationSystems & _eqn_systems
virtual std::unique_ptr< FunctionBase< Number > > clone() const override
const NumericVector< Number > & _vector
PointLocatorBase * _point_locator
MeshFunction & operator=(const MeshFunction &)=delete
MeshFunction(const EquationSystems &eqn_systems, const NumericVector< Number > &vec, const DofMap &dof_map, const std::vector< unsigned int > &vars, const FunctionBase< Number > *master=nullptr)
Definition: mesh_function.C:43
const DofMap & _dof_map
const Elem * find_element(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
The base class for all geometric element types.
Definition: elem.h:100
Gradient gradient(const Point &p, const Real time=0.)
void disable_out_of_mesh_mode(void)
Manages the degrees of freedom (DOFs) in a simulation.
Definition: dof_map.h:176
DenseVector< Number > _out_of_mesh_value
void enable_out_of_mesh_mode(const DenseVector< Number > &value)
virtual void clear() override
std::set< const Elem * > find_elements(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
std::map< const Elem *, Number > discontinuous_value(const Point &p, const Real time=0.)
void unset_point_locator_tolerance()
An object whose state is distributed along a set of processors.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const PointLocatorBase & get_point_locator(void) const
Number operator()(const Point &p, const Real time=0.) override
static const bool value
Definition: xdr_io.C:109
BuildType
Base class for different Tree types.
Definition: tree_base.h:55
void set_point_locator_tolerance(Real tol)
Base class for functors that can be evaluated at a point and (optionally) time.
std::map< const Elem *, Gradient > discontinuous_gradient(const Point &p, const Real time=0.)
A geometric point in (x,y,z) space.
Definition: point.h:38
Tensor hessian(const Point &p, const Real time=0.)
virtual void init() override
const std::vector< unsigned int > _system_vars