Line data Source code
1 : // The libMesh Finite Element Library.
2 : // Copyright (C) 2002-2025 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 156 : 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 : #ifdef LIBMESH_ENABLE_DEPRECATED
113 : /**
114 : * The actual initialization process. Takes an optional argument which
115 : * specifies the method to use when building a \p PointLocator
116 : *
117 : * \deprecated The input argument is not used (was it ever?) to
118 : * control the PointLocator type which is built, so one should
119 : * instead call the version of init() taking no args.
120 : */
121 : void init (const Trees::BuildType point_locator_build_type);
122 : #endif // LIBMESH_ENABLE_DEPRECATED
123 :
124 : /**
125 : * Clears the function.
126 : */
127 : virtual void clear () override;
128 :
129 : /**
130 : * \returns A new copy of the function.
131 : *
132 : * The new copy uses the original as a master function to enable
133 : * simultaneous evaluations of the copies in different threads.
134 : *
135 : * \note This implies the copy should not be used after the
136 : * original is destroyed.
137 : */
138 : virtual std::unique_ptr<FunctionBase<Number>> clone () const override;
139 :
140 : /**
141 : * \returns The value of variable 0 at point \p p and for \p time,
142 : * which defaults to zero.
143 : */
144 : virtual Number operator() (const Point & p,
145 : const Real time=0.) override;
146 :
147 : /**
148 : * \returns A map of values of variable 0 at point
149 : * \p p and for \p time.
150 : *
151 : * The std::map is from element to Number and accounts for
152 : * doubly-defined values on faces if discontinuous variables are
153 : * used.
154 : */
155 : std::map<const Elem *, Number> discontinuous_value (const Point & p,
156 : const Real time=0.);
157 :
158 : /**
159 : * \returns The first derivatives of variable 0 at point
160 : * \p p and for \p time, which defaults to zero.
161 : */
162 : Gradient gradient (const Point & p,
163 : const Real time=0.);
164 :
165 : /**
166 : * \returns A map of first derivatives (gradients) of variable 0 at point
167 : * \p p and for \p time.
168 : * map is from element to Gradient and accounts for double defined
169 : * values on faces if the gradient is discontinuous
170 : */
171 : std::map<const Elem *, Gradient> discontinuous_gradient (const Point & p,
172 : const Real time=0.);
173 :
174 : #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
175 : /**
176 : * \returns The second derivatives of variable 0 at point
177 : * \p p and for \p time, which defaults to zero.
178 : */
179 : Tensor hessian (const Point & p,
180 : const Real time=0.);
181 : #endif
182 :
183 : /**
184 : * Computes values at coordinate \p p and for time \p time, which
185 : * defaults to zero, optionally restricting the point to the
186 : * MeshFunction subdomain_ids. This is useful in cases where there
187 : * are multiple dimensioned elements, for example.
188 : */
189 : virtual void operator() (const Point & p,
190 : const Real time,
191 : DenseVector<Number> & output) override;
192 :
193 : /**
194 : * Computes values at coordinate \p p and for time \p time,
195 : * restricting the point to the passed subdomain_ids, which
196 : * parameter overrides the internal subdomain_ids.
197 : */
198 : void operator() (const Point & p,
199 : const Real time,
200 : DenseVector<Number> & output,
201 : const std::set<subdomain_id_type> * subdomain_ids);
202 :
203 : /**
204 : * Similar to operator() with the same parameter list, but with the difference
205 : * that multiple values on faces are explicitly permitted. This is useful for
206 : * discontinuous shape functions that are evaluated on faces.
207 : */
208 : void discontinuous_value (const Point & p,
209 : const Real time,
210 : std::map<const Elem *, DenseVector<Number>> & output);
211 :
212 : /**
213 : * Similar to operator() with the same parameter list, but with the difference
214 : * that multiple values on faces are explicitly permitted. This is useful for
215 : * discontinuous shape functions that are evaluated on faces.
216 : */
217 : void discontinuous_value (const Point & p,
218 : const Real time,
219 : std::map<const Elem *, DenseVector<Number>> & output,
220 : const std::set<subdomain_id_type> * subdomain_ids);
221 :
222 : /**
223 : * Computes gradients at coordinate \p p and for time \p time, which
224 : * defaults to zero, optionally restricting the point to the
225 : * MeshFunction subdomain_ids.
226 : */
227 : void gradient (const Point & p,
228 : const Real time,
229 : std::vector<Gradient> & output);
230 :
231 : /**
232 : * Computes gradients at coordinate \p p and for time \p time, which
233 : * defaults to zero, optionally restricting the point to the passed
234 : * subdomain_ids, which parameter overrides the internal
235 : * subdomain_ids.
236 : */
237 : void gradient (const Point & p,
238 : const Real time,
239 : std::vector<Gradient> & output,
240 : const std::set<subdomain_id_type> * subdomain_ids);
241 :
242 : /**
243 : * Similar to gradient, but with the difference
244 : * that multiple values on faces are explicitly permitted. This is useful for
245 : * evaluating gradients on faces where the values to the left and right are different.
246 : */
247 : void discontinuous_gradient (const Point & p,
248 : const Real time,
249 : std::map<const Elem *, std::vector<Gradient>> & output);
250 :
251 : /**
252 : * Similar to gradient, but with the difference
253 : * that multiple values on faces are explicitly permitted. This is useful for
254 : * evaluating gradients on faces where the values to the left and right are different.
255 : */
256 : void discontinuous_gradient (const Point & p,
257 : const Real time,
258 : std::map<const Elem *, std::vector<Gradient>> & output,
259 : const std::set<subdomain_id_type> * subdomain_ids);
260 :
261 : /**
262 : * Computes gradients at coordinate \p p and for time \p time, which
263 : * defaults to zero, optionally restricting the point to the
264 : * MeshFunction subdomain_ids.
265 : */
266 : void hessian (const Point & p,
267 : const Real time,
268 : std::vector<Tensor> & output);
269 :
270 : /**
271 : * Computes gradients at coordinate \p p and for time \p time, which
272 : * defaults to zero, optionally restricting the point to the passed
273 : * subdomain_ids. This is useful in cases where there are multiple
274 : * dimensioned elements, for example.
275 : */
276 : void hessian (const Point & p,
277 : const Real time,
278 : std::vector<Tensor> & output,
279 : const std::set<subdomain_id_type> * subdomain_ids);
280 :
281 : /**
282 : * \returns The current \p PointLocator object, for use elsewhere.
283 : *
284 : * \note The \p MeshFunction object must be initialized before this
285 : * is called.
286 : */
287 : const PointLocatorBase & get_point_locator () const;
288 : PointLocatorBase & get_point_locator ();
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 DenseVector<Number> & value);
301 :
302 : /**
303 : * Enables out-of-mesh mode. In this mode, if asked for a point
304 : * that is not contained in any element, the \p MeshFunction will
305 : * return the given \p value instead of crashing. This mode is off
306 : * per default. If you use a master mesh function and you want to
307 : * enable this mode, you will have to enable it for the master mesh
308 : * function as well and for all mesh functions that have the same
309 : * master mesh function. You may, however, specify different
310 : * values.
311 : */
312 : void enable_out_of_mesh_mode(const Number & value);
313 :
314 : /**
315 : * Disables out-of-mesh mode. This is also the default.
316 : */
317 : void disable_out_of_mesh_mode();
318 :
319 : /**
320 : * We may want to specify a tolerance for the PointLocator to use,
321 : * since in some cases the point we want to evaluate at might be
322 : * slightly outside the mesh (due to numerical rounding issues,
323 : * for example).
324 : */
325 : void set_point_locator_tolerance(Real tol);
326 :
327 : /**
328 : * Turn off the user-specified PointLocator tolerance.
329 : */
330 : void unset_point_locator_tolerance();
331 :
332 : /**
333 : * Choose a default list of subdomain ids to be searched for points.
334 : * If the provided list pointer is null or if no list has been
335 : * provided, then all subdomain ids are searched. This list can be
336 : * overridden on a per-evaluation basis by using the method
337 : * overrides with a similar argument.
338 : */
339 : void set_subdomain_ids(const std::set<subdomain_id_type> * subdomain_ids);
340 :
341 : protected:
342 :
343 : /**
344 : * Helper function to reduce code duplication
345 : */
346 : const Elem * find_element(const Point & p,
347 : const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
348 :
349 : /**
350 : * \returns All elements that are close to a point \p p.
351 : *
352 : * This is similar to \p find_element() but covers cases where \p p
353 : * is on the boundary.
354 : */
355 : std::set<const Elem *> find_elements(const Point & p,
356 : const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
357 :
358 : /**
359 : * Helper function that is called by MeshFunction::find_element()
360 : * and MeshFunction::find_elements() to ensure that Elems found by
361 : * the PointLocator are actually "evaluable" on the processor where
362 : * they are found.
363 : */
364 : const Elem * check_found_elem(const Elem * element, const Point & p) const;
365 :
366 : /**
367 : * Helper function for finding a gradient as evaluated from a
368 : * specific element
369 : */
370 : void _gradient_on_elem (const Point & p,
371 : const Elem * element,
372 : std::vector<Gradient> & output);
373 :
374 : /**
375 : * The equation systems handler, from which
376 : * the data are gathered.
377 : */
378 : const EquationSystems & _eqn_systems;
379 :
380 : /**
381 : * A reference to the vector that holds the data
382 : * that is to be interpolated.
383 : */
384 : const NumericVector<Number> & _vector;
385 :
386 : /**
387 : * Need access to the \p DofMap of the other system.
388 : */
389 : const DofMap & _dof_map;
390 :
391 : /**
392 : * The indices of the variables within the other system
393 : * for which data are to be gathered.
394 : */
395 : const std::vector<unsigned int> _system_vars;
396 :
397 : /**
398 : * A point locator is needed to locate the
399 : * points in the mesh.
400 : */
401 : std::unique_ptr<PointLocatorBase> _point_locator;
402 :
403 : /**
404 : * A default set of subdomain ids in which to search for points.
405 : */
406 : std::unique_ptr<std::set<subdomain_id_type>> _subdomain_ids;
407 :
408 : /**
409 : * \p true if out-of-mesh mode is enabled. See \p
410 : * enable_out_of_mesh_mode() for more details. Default is \p false.
411 : */
412 : bool _out_of_mesh_mode;
413 :
414 : /**
415 : * Value to return outside the mesh if out-of-mesh mode is enabled.
416 : * See \p enable_out_of_mesh_mode() for more details.
417 : */
418 : DenseVector<Number> _out_of_mesh_value;
419 : };
420 :
421 :
422 : } // namespace libMesh
423 :
424 :
425 : #endif // LIBMESH_MESH_FUNCTION_H
|
