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 : #ifndef LIBMESH_MESH_TRIANGLE_HOLES_H
19 : #define LIBMESH_MESH_TRIANGLE_HOLES_H
20 :
21 : #include "libmesh/libmesh_config.h"
22 :
23 : // Local includes
24 : #include "libmesh/boundary_info.h" // invalid_id
25 : #include "libmesh/point.h"
26 : #include "libmesh/triangulator_interface.h"
27 : #include "libmesh/vector_value.h"
28 :
29 : // C++ includes
30 :
31 : namespace libMesh
32 : {
33 :
34 : /**
35 : * An abstract class for defining a 2-dimensional hole. We assume that
36 : * the connectivity of the hole is implicit in the numbering of the points,
37 : * controlled by segment_indices vector.
38 : * The size of segment_indices is equal to the number of segments plus one.
39 : * Each segment has segment_indices[i+1]-segment_indices[i] connected points,
40 : * with node segment_indices[i] is connected to node segment_indices[i+1],
41 : * node segment_indices[i+1] is connected to node segment_indices[i+2], etc,
42 : * and the last node "wraps around" to connect back to node segment_indices[i].
43 : *
44 : * \author John W. Peterson
45 : * \date 2011
46 : * \brief Class for parameterizing 2D holes to be meshed with Triangle.
47 : */
48 : class TriangulatorInterface::Hole
49 : {
50 : public:
51 : /**
52 : * Constructor
53 : */
54 1678 : Hole() = default;
55 :
56 : /**
57 : * Destructor
58 : */
59 0 : virtual ~Hole() = default;
60 :
61 : /**
62 : * The number of geometric points which define the hole.
63 : */
64 : virtual unsigned int n_points() const = 0;
65 :
66 : /**
67 : * The number of geometric midpoints along *each* of the sides
68 : * defining the hole.
69 : */
70 0 : virtual unsigned int n_midpoints() const { return 0; }
71 :
72 : /**
73 : * Return the nth point defining the hole.
74 : */
75 : virtual Point point(const unsigned int n) const = 0;
76 :
77 : /**
78 : * Return the midpoint \p m along the side \p n defining the hole.
79 : */
80 0 : virtual Point midpoint(const unsigned int /* m */,
81 : const unsigned int /* n */) const
82 0 : { libmesh_error(); /* by default holes are polygonal */ }
83 :
84 : /**
85 : * Return an (arbitrary) point which lies inside the hole.
86 : */
87 : virtual Point inside() const = 0;
88 :
89 : /**
90 : * Return true iff \p p lies inside the hole.
91 : *
92 : * This method currently does not take any higher-order hole
93 : * geometry into account, but treats the hole as a polygon.
94 : */
95 : bool contains(Point p) const;
96 :
97 : /**
98 : * Return the area of the hole
99 : *
100 : * This method currently does not take any higher-order hole
101 : * geometry into account, but treats the hole as a polygon.
102 : */
103 : Real area() const;
104 :
105 : /**
106 : * Return a vector with right-hand-rule orientation and length of
107 : * twice area() squared. This is useful for determining
108 : * orientation of non-planar or non-counter-clockwise holes.
109 : *
110 : * This method currently does not take any higher-order hole
111 : * geometry into account, but treats the hole as a polygon.
112 : */
113 : RealGradient areavec() const;
114 :
115 : /**
116 : * Starting indices of points for a hole with multiple disconnected boundaries.
117 : */
118 30 : virtual std::vector<unsigned int> segment_indices() const
119 : {
120 : // default to only one enclosing boundary
121 0 : std::vector<unsigned int> seg;
122 30 : seg.push_back(0);
123 30 : seg.push_back(n_points());
124 30 : return seg;
125 : }
126 :
127 : /**
128 : * Set whether or not a triangulator is allowed to refine the
129 : * hole boundary when refining the mesh interior. This is true by
130 : * default, but may be set to false to make the hole boundary more
131 : * predictable (and so easier to stitch to other meshes) later.
132 : */
133 0 : virtual void set_refine_boundary_allowed (bool refine_bdy_allowed)
134 0 : { _refine_bdy_allowed = refine_bdy_allowed; }
135 :
136 : /**
137 : * Get whether or not the triangulation is allowed to refine the
138 : * mesh boundary when refining the interior. True by default.
139 : */
140 2760 : virtual bool refine_boundary_allowed () const
141 2760 : { return _refine_bdy_allowed; }
142 :
143 : protected:
144 :
145 : /**
146 : * Helper function for contains(), also useful for MeshedHole::inside()
147 : */
148 : std::vector<Real> find_ray_intersections(Point ray_start,
149 : Point ray_target) const;
150 :
151 : /**
152 : * Calculate an inside point based on our boundary
153 : */
154 : Point calculate_inside_point() const;
155 :
156 : /**
157 : * Whether to allow boundary refinement. True by default; specified
158 : * here so we can use the default constructor.
159 : */
160 : bool _refine_bdy_allowed = true;
161 : };
162 :
163 :
164 :
165 :
166 :
167 : /**
168 : * A concrete instantiation of the Hole class that describes polygonal
169 : * (triangular, square, pentagonal, ...) holes.
170 : */
171 : class TriangulatorInterface::PolygonHole : public TriangulatorInterface::Hole
172 : {
173 : public:
174 : /**
175 : * Constructor specifying the center, radius, and number of
176 : * points which comprise the hole. The points will all lie
177 : * on a circle of radius r.
178 : */
179 : PolygonHole(const Point & center, Real radius, unsigned int n_points);
180 :
181 : virtual unsigned int n_points() const override;
182 :
183 : virtual Point point(const unsigned int n) const override;
184 :
185 : virtual Point inside() const override;
186 :
187 : private:
188 : /**
189 : * (x,y) location of the center of the hole
190 : */
191 : Point _center;
192 :
193 : /**
194 : * circular hole radius
195 : */
196 : Real _radius;
197 :
198 : /**
199 : * number of points used to describe the hole. The actual
200 : * points can be generated knowing the center and radius.
201 : * For example, n_points=3 would generate a triangular hole.
202 : */
203 : unsigned int _n_points;
204 : };
205 :
206 :
207 :
208 : /**
209 : * A way to translate and/or rotate an existing hole; perhaps to tile
210 : * it in many places or to put it at an angle that the underlying hole
211 : * doesn't support.
212 : */
213 : class TriangulatorInterface::AffineHole : public TriangulatorInterface::Hole
214 : {
215 : public:
216 : /**
217 : * Constructor specifying the underlying hole, and the rotation
218 : * angle (in radians, done first) and translation (done second) with
219 : * which to transform it.
220 : */
221 : AffineHole(const Hole & underlying, Real angle, const Point & shift)
222 : : _underlying(underlying), _angle(angle), _shift(shift) {}
223 :
224 30033 : virtual unsigned int n_points() const override
225 30033 : { return _underlying.n_points(); }
226 :
227 : virtual Point point(const unsigned int n) const override;
228 :
229 : virtual Point inside() const override;
230 :
231 : private:
232 : /**
233 : * Rotate-and-shift equations
234 : */
235 : Point transform(const Point & p) const;
236 :
237 : /**
238 : * Hole to transform
239 : */
240 : const Hole & _underlying;
241 :
242 : /**
243 : * Angle to rotate (counter-clockwise) by
244 : */
245 : Real _angle;
246 :
247 : /**
248 : * (x,y) location to shift (0,0) to
249 : */
250 : Point _shift;
251 : };
252 :
253 :
254 :
255 :
256 :
257 : /**
258 : * Another concrete instantiation of the hole, this one should
259 : * be sufficiently general for most non-polygonal purposes. The user
260 : * supplies, at the time of construction, a reference to a vector
261 : * of Points which defines the hole (in order of connectivity) and
262 : * an arbitrary Point which lies inside the hole.
263 : */
264 : class TriangulatorInterface::ArbitraryHole : public TriangulatorInterface::Hole
265 : {
266 : public:
267 : /**
268 : * The fastest constructor requires a point which lies in the
269 : * interior of the hole and a reference to a vector of Points
270 : * defining the hole.
271 : */
272 : ArbitraryHole(const Point & center,
273 : std::vector<Point> points);
274 :
275 : ArbitraryHole(const Point & center,
276 : std::vector<Point> points,
277 : std::vector<unsigned int> segment_indices);
278 :
279 : /**
280 : * If we don't know a center a priori then we can calculate one
281 : */
282 : ArbitraryHole(std::vector<Point> points);
283 :
284 : /**
285 : * We can also construct an ArbitraryHole which just copies
286 : * a hole of any other type.
287 : */
288 : ArbitraryHole(const Hole & orig);
289 :
290 : ArbitraryHole(const ArbitraryHole & orig) = default;
291 :
292 : virtual unsigned int n_points() const override;
293 :
294 : virtual Point point(const unsigned int n) const override;
295 :
296 : virtual Point inside() const override;
297 :
298 : virtual std::vector<unsigned int> segment_indices() const override;
299 :
300 70 : const std::vector<Point> & get_points() const
301 : {
302 70 : return _points;
303 : }
304 :
305 621 : void set_points(std::vector<Point> points)
306 : {
307 639 : _points = std::move(points);
308 639 : _center = this->calculate_inside_point();
309 621 : }
310 :
311 : private:
312 : /**
313 : * arbitrary (x,y) location inside the hole
314 : */
315 : Point _center;
316 :
317 : /**
318 : * Reference to the vector of points which makes up
319 : * the hole.
320 : */
321 : std::vector<Point> _points;
322 :
323 : std::vector<unsigned int> _segment_indices;
324 : };
325 :
326 :
327 :
328 : /**
329 : * Another concrete instantiation of the hole, as general as
330 : * ArbitraryHole, but based on an existing 1D or 2D mesh.
331 : *
332 : * If ids are given, 2D edges on a boundary with a listed id or 1D
333 : * edges in a subdomain with a listed id will define the hole.
334 : *
335 : * If no ids are given, the hole will be defined by all 1D Edge
336 : * elements and all outward-facing 2D boundary edges.
337 : *
338 : * In either case, the hole definition should give a single connected
339 : * boundary, topologically a circle. The hole is defined when the
340 : * MeshedHole is constructed, and ignores any subsequent changes to
341 : * the input mesh.
342 : */
343 : class TriangulatorInterface::MeshedHole : public TriangulatorInterface::Hole
344 : {
345 : public:
346 : /**
347 : * The constructor requires a mesh defining the hole, and optionally
348 : * boundary+subdomain ids restricting the definition.
349 : */
350 : MeshedHole(const MeshBase & mesh,
351 : std::set<std::size_t> ids = {});
352 :
353 : virtual unsigned int n_points() const override;
354 :
355 : virtual unsigned int n_midpoints() const override;
356 :
357 : virtual Point point(const unsigned int n) const override;
358 :
359 : virtual Point midpoint(const unsigned int m,
360 : const unsigned int n) const override;
361 :
362 : virtual Point inside() const override;
363 :
364 : private:
365 : /**
366 : * An (x,y) location inside the hole. Cached because this is too
367 : * expensive to compute for an arbitrary input mesh unless we need
368 : * it for Triangle.
369 : */
370 : mutable Point _center;
371 :
372 : /**
373 : * The sorted vector of points which makes up the hole.
374 : */
375 : std::vector<Point> _points;
376 :
377 : /**
378 : * The sorted vector of midpoints in between points along the edges
379 : * of the hole. For a hole with m midpoints per edge, between
380 : * _points[n] and _points[n+1] lies _midpoints[n*m] through
381 : * _midpoints[n*m+m-1]
382 : */
383 : std::vector<Point> _midpoints;
384 : };
385 :
386 :
387 :
388 :
389 :
390 :
391 : /**
392 : * A class for defining a 2-dimensional region for Triangle.
393 : */
394 : class TriangulatorInterface::Region
395 : {
396 : public:
397 : /**
398 : * Constructor
399 : */
400 : Region(const Point & center,
401 : Real attribute = 0,
402 : Real max_area = -std::numeric_limits<Real>::max())
403 : : _center(center),
404 : _attribute(attribute),
405 : _max_area(max_area) {}
406 :
407 0 : Point inside() const { return _center; }
408 :
409 0 : Real attribute() const { return _attribute; }
410 :
411 0 : Real max_area() const { return _max_area; }
412 :
413 : private:
414 : /**
415 : * Arbitrary (x,y) location inside the region
416 : */
417 : const Point _center;
418 :
419 : /**
420 : * Attribute of the region
421 : * Default value for attribute is zero.
422 : */
423 : Real _attribute;
424 :
425 : /**
426 : * Maximum area that is allowed for all triangles in this region
427 : * Default negative maximum area means that no area constraint will be imposed.
428 : */
429 : Real _max_area;
430 : };
431 :
432 : } // namespace libMesh
433 :
434 : #endif // LIBMESH_MESH_TRIANGLE_HOLES_H
|
