Line data Source code
1 : // The libMesh Finite Element Library.
2 : // Copyright (C) 2002-2026 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"
29 : #include "libmesh/parallel_object.h"
30 :
31 : // C++ includes
32 : #include <cstddef>
33 : #include <vector>
34 : #include <memory>
35 :
36 : namespace libMesh
37 : {
38 :
39 :
40 : // Forward Declarations
41 : template <typename T> class DenseVector;
42 : class EquationSystems;
43 : template <typename T> class NumericVector;
44 : class DofMap;
45 : class PointLocatorBase;
46 :
47 : /**
48 : * This class provides function-like objects for data
49 : * distributed over a mesh.
50 : *
51 : * \author Daniel Dreyer
52 : * \date 2003
53 : */
54 124 : class MeshFunction : public FunctionBase<Number>,
55 : public ParallelObject
56 : {
57 : public:
58 :
59 : /**
60 : * Constructor for mesh based functions with vectors
61 : * as return value. Optionally takes a master function.
62 : * If the MeshFunction is to be evaluated outside of the
63 : * local partition of the mesh, then both the mesh in
64 : * \p eqn_systems and the coefficient vector \p vec
65 : * should be serialized.
66 : */
67 : MeshFunction (const EquationSystems & eqn_systems,
68 : const NumericVector<Number> & vec,
69 : const DofMap & dof_map,
70 : std::vector<unsigned int> vars,
71 : const FunctionBase<Number> * master=nullptr);
72 :
73 : /**
74 : * Constructor for mesh based functions with a number
75 : * as return value. Optionally takes a master function.
76 : * If the MeshFunction is to be evaluated outside of the
77 : * local partition of the mesh, then both the mesh in
78 : * \p eqn_systems and the coefficient vector \p vec
79 : * should be serialized.
80 : */
81 : MeshFunction (const EquationSystems & eqn_systems,
82 : const NumericVector<Number> & vec,
83 : const DofMap & dof_map,
84 : const unsigned int var,
85 : const FunctionBase<Number> * master=nullptr);
86 :
87 : /**
88 : * A regular copy constructor.
89 : */
90 : MeshFunction (const MeshFunction & mf);
91 :
92 : /**
93 : * Special functions.
94 : * - This class contains a unique_ptr so it can't be default copy constructed.
95 : * - This class contains const references so it can't be default copy/move assigned.
96 : * - The destructor is defaulted out-of-line.
97 : */
98 : MeshFunction (MeshFunction &&) = default;
99 : MeshFunction & operator= (const MeshFunction &) = delete;
100 : MeshFunction & operator= (MeshFunction &&) = delete;
101 :
102 : /**
103 : * Destructor.
104 : */
105 : ~MeshFunction ();
106 :
107 : /**
108 : * Override the FunctionBase::init() member function.
109 : */
110 : virtual void init () override;
111 :
112 : /**
113 : * Clears the function.
114 : */
115 : virtual void clear () override;
116 :
117 : /**
118 : * \returns A new copy of the function.
119 : *
120 : * The new copy uses the original as a master function to enable
121 : * simultaneous evaluations of the copies in different threads.
122 : *
123 : * \note This implies the copy should not be used after the
124 : * original is destroyed.
125 : */
126 : virtual std::unique_ptr<FunctionBase<Number>> clone () const override;
127 :
128 : /**
129 : * \returns The value of variable 0 at point \p p and for \p time,
130 : * which defaults to zero.
131 : */
132 : virtual Number operator() (const Point & p,
133 : const Real time=0.) override;
134 :
135 : /**
136 : * \returns A map of values of variable 0 at point
137 : * \p p and for \p time.
138 : *
139 : * The std::map is from element to Number and accounts for
140 : * doubly-defined values on faces if discontinuous variables are
141 : * used.
142 : */
143 : std::map<const Elem *, Number> discontinuous_value (const Point & p,
144 : const Real time=0.);
145 :
146 : /**
147 : * \returns The first derivatives of variable 0 at point
148 : * \p p and for \p time, which defaults to zero.
149 : */
150 : Gradient gradient (const Point & p,
151 : const Real time=0.);
152 :
153 : /**
154 : * \returns A map of first derivatives (gradients) of variable 0 at point
155 : * \p p and for \p time.
156 : * map is from element to Gradient and accounts for double defined
157 : * values on faces if the gradient is discontinuous
158 : */
159 : std::map<const Elem *, Gradient> discontinuous_gradient (const Point & p,
160 : const Real time=0.);
161 :
162 : #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
163 : /**
164 : * \returns The second derivatives of variable 0 at point
165 : * \p p and for \p time, which defaults to zero.
166 : */
167 : Tensor hessian (const Point & p,
168 : const Real time=0.);
169 : #endif
170 :
171 : /**
172 : * Computes values at coordinate \p p and for time \p time, which
173 : * defaults to zero, optionally restricting the point to the
174 : * MeshFunction subdomain_ids. This is useful in cases where there
175 : * are multiple dimensioned elements, for example.
176 : */
177 : virtual void operator() (const Point & p,
178 : const Real time,
179 : DenseVector<Number> & output) override;
180 :
181 : /**
182 : * Computes values at coordinate \p p and for time \p time,
183 : * restricting the point to the passed subdomain_ids, which
184 : * parameter overrides the internal subdomain_ids.
185 : */
186 : void operator() (const Point & p,
187 : const Real time,
188 : DenseVector<Number> & output,
189 : const std::set<subdomain_id_type> * subdomain_ids);
190 :
191 : /**
192 : * Similar to operator() with the same parameter list, but with the difference
193 : * that multiple values on faces are explicitly permitted. This is useful for
194 : * discontinuous shape functions that are evaluated on faces.
195 : */
196 : void discontinuous_value (const Point & p,
197 : const Real time,
198 : std::map<const Elem *, DenseVector<Number>> & output);
199 :
200 : /**
201 : * Similar to operator() with the same parameter list, but with the difference
202 : * that multiple values on faces are explicitly permitted. This is useful for
203 : * discontinuous shape functions that are evaluated on faces.
204 : */
205 : void discontinuous_value (const Point & p,
206 : const Real time,
207 : std::map<const Elem *, DenseVector<Number>> & output,
208 : const std::set<subdomain_id_type> * subdomain_ids);
209 :
210 : /**
211 : * Computes gradients at coordinate \p p and for time \p time, which
212 : * defaults to zero, optionally restricting the point to the
213 : * MeshFunction subdomain_ids.
214 : */
215 : void gradient (const Point & p,
216 : const Real time,
217 : std::vector<Gradient> & output);
218 :
219 : /**
220 : * Computes gradients at coordinate \p p and for time \p time, which
221 : * defaults to zero, optionally restricting the point to the passed
222 : * subdomain_ids, which parameter overrides the internal
223 : * subdomain_ids.
224 : */
225 : void gradient (const Point & p,
226 : const Real time,
227 : std::vector<Gradient> & output,
228 : const std::set<subdomain_id_type> * subdomain_ids);
229 :
230 : /**
231 : * Similar to gradient, but with the difference
232 : * that multiple values on faces are explicitly permitted. This is useful for
233 : * evaluating gradients on faces where the values to the left and right are different.
234 : */
235 : void discontinuous_gradient (const Point & p,
236 : const Real time,
237 : std::map<const Elem *, std::vector<Gradient>> & output);
238 :
239 : /**
240 : * Similar to gradient, but with the difference
241 : * that multiple values on faces are explicitly permitted. This is useful for
242 : * evaluating gradients on faces where the values to the left and right are different.
243 : */
244 : void discontinuous_gradient (const Point & p,
245 : const Real time,
246 : std::map<const Elem *, std::vector<Gradient>> & output,
247 : const std::set<subdomain_id_type> * subdomain_ids);
248 :
249 : /**
250 : * Computes gradients at coordinate \p p and for time \p time, which
251 : * defaults to zero, optionally restricting the point to the
252 : * MeshFunction subdomain_ids.
253 : */
254 : void hessian (const Point & p,
255 : const Real time,
256 : std::vector<Tensor> & output);
257 :
258 : /**
259 : * Computes gradients at coordinate \p p and for time \p time, which
260 : * defaults to zero, optionally restricting the point to the passed
261 : * subdomain_ids. This is useful in cases where there are multiple
262 : * dimensioned elements, for example.
263 : */
264 : void hessian (const Point & p,
265 : const Real time,
266 : std::vector<Tensor> & output,
267 : const std::set<subdomain_id_type> * subdomain_ids);
268 :
269 : /**
270 : * \returns The current \p PointLocator object, for use elsewhere.
271 : *
272 : * \note The \p MeshFunction object must be initialized before this
273 : * is called.
274 : */
275 : const PointLocatorBase & get_point_locator () const;
276 : PointLocatorBase & get_point_locator ();
277 :
278 : /**
279 : * Enables out-of-mesh mode. In this mode, if asked for a point
280 : * that is not contained in any element, the \p MeshFunction will
281 : * return the given \p value instead of crashing. This mode is off
282 : * per default. If you use a master mesh function and you want to
283 : * enable this mode, you will have to enable it for the master mesh
284 : * function as well and for all mesh functions that have the same
285 : * master mesh function. You may, however, specify different
286 : * values.
287 : */
288 : void enable_out_of_mesh_mode(const DenseVector<Number> & value);
289 :
290 : /**
291 : * Enables out-of-mesh mode. In this mode, if asked for a point
292 : * that is not contained in any element, the \p MeshFunction will
293 : * return the given \p value instead of crashing. This mode is off
294 : * per default. If you use a master mesh function and you want to
295 : * enable this mode, you will have to enable it for the master mesh
296 : * function as well and for all mesh functions that have the same
297 : * master mesh function. You may, however, specify different
298 : * values.
299 : */
300 : void enable_out_of_mesh_mode(const Number & value);
301 :
302 : /**
303 : * Disables out-of-mesh mode. This is also the default.
304 : */
305 : void disable_out_of_mesh_mode();
306 :
307 : /**
308 : * We may want to specify a tolerance for the PointLocator to use,
309 : * since in some cases the point we want to evaluate at might be
310 : * slightly outside the mesh (due to numerical rounding issues,
311 : * for example).
312 : */
313 : void set_point_locator_tolerance(Real tol);
314 :
315 : /**
316 : * Turn off the user-specified PointLocator tolerance.
317 : */
318 : void unset_point_locator_tolerance();
319 :
320 : /**
321 : * Choose a default list of subdomain ids to be searched for points.
322 : * If the provided list pointer is null or if no list has been
323 : * provided, then all subdomain ids are searched. This list can be
324 : * overridden on a per-evaluation basis by using the method
325 : * overrides with a similar argument.
326 : */
327 : void set_subdomain_ids(const std::set<subdomain_id_type> * subdomain_ids);
328 :
329 : protected:
330 :
331 : /**
332 : * Helper function to reduce code duplication
333 : */
334 : const Elem * find_element(const Point & p,
335 : const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
336 :
337 : /**
338 : * \returns All elements that are close to a point \p p.
339 : *
340 : * This is similar to \p find_element() but covers cases where \p p
341 : * is on the boundary.
342 : */
343 : std::set<const Elem *> find_elements(const Point & p,
344 : const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
345 :
346 : /**
347 : * Helper function that is called by MeshFunction::find_element()
348 : * and MeshFunction::find_elements() to ensure that Elems found by
349 : * the PointLocator are actually "evaluable" on the processor where
350 : * they are found.
351 : */
352 : const Elem * check_found_elem(const Elem * element, const Point & p) const;
353 :
354 : /**
355 : * Helper function for finding a gradient as evaluated from a
356 : * specific element
357 : */
358 : void _gradient_on_elem (const Point & p,
359 : const Elem * element,
360 : std::vector<Gradient> & output);
361 :
362 : /**
363 : * The equation systems handler, from which
364 : * the data are gathered.
365 : */
366 : const EquationSystems & _eqn_systems;
367 :
368 : /**
369 : * A reference to the vector that holds the data
370 : * that is to be interpolated.
371 : */
372 : const NumericVector<Number> & _vector;
373 :
374 : /**
375 : * Need access to the \p DofMap of the other system.
376 : */
377 : const DofMap & _dof_map;
378 :
379 : /**
380 : * The indices of the variables within the other system
381 : * for which data are to be gathered.
382 : */
383 : const std::vector<unsigned int> _system_vars;
384 :
385 : /**
386 : * A point locator is needed to locate the
387 : * points in the mesh.
388 : */
389 : std::unique_ptr<PointLocatorBase> _point_locator;
390 :
391 : /**
392 : * A default set of subdomain ids in which to search for points.
393 : */
394 : std::unique_ptr<std::set<subdomain_id_type>> _subdomain_ids;
395 :
396 : /**
397 : * \p true if out-of-mesh mode is enabled. See \p
398 : * enable_out_of_mesh_mode() for more details. Default is \p false.
399 : */
400 : bool _out_of_mesh_mode;
401 :
402 : /**
403 : * Value to return outside the mesh if out-of-mesh mode is enabled.
404 : * See \p enable_out_of_mesh_mode() for more details.
405 : */
406 : DenseVector<Number> _out_of_mesh_value;
407 : };
408 :
409 :
410 : } // namespace libMesh
411 :
412 :
413 : #endif // LIBMESH_MESH_FUNCTION_H
|
