Loading [MathJax]/extensions/tex2jax.js
ug4
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
anisotropy_util_impl.h
Go to the documentation of this file.
1/*
2 * Copyright (c) 2009-2015: G-CSC, Goethe University Frankfurt
3 * Author: Markus Breit
4 * Date: 2018-05-25
5 *
6 * This file is part of UG4.
7 *
8 * UG4 is free software: you can redistribute it and/or modify it under the
9 * terms of the GNU Lesser General Public License version 3 (as published by the
10 * Free Software Foundation) with the following additional attribution
11 * requirements (according to LGPL/GPL v3 §7):
12 *
13 * (1) The following notice must be displayed in the Appropriate Legal Notices
14 * of covered and combined works: "Based on UG4 (www.ug4.org/license)".
15 *
16 * (2) The following notice must be displayed at a prominent place in the
17 * terminal output of covered works: "Based on UG4 (www.ug4.org/license)".
18 *
19 * (3) The following bibliography is recommended for citation and must be
20 * preserved in all covered files:
21 * "Reiter, S., Vogel, A., Heppner, I., Rupp, M., and Wittum, G. A massively
22 * parallel geometric multigrid solver on hierarchically distributed grids.
23 * Computing and visualization in science 16, 4 (2013), 151-164"
24 * "Vogel, A., Reiter, S., Rupp, M., Nägel, A., and Wittum, G. UG4 -- a novel
25 * flexible software system for simulating pde based models on high performance
26 * computers. Computing and visualization in science 16, 4 (2013), 165-179"
27 *
28 * This program is distributed in the hope that it will be useful,
29 * but WITHOUT ANY WARRANTY; without even the implied warranty of
30 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
31 * GNU Lesser General Public License for more details.
32 */
33
34#include "lib_grid/algorithms/element_side_util.h" // for GetOpposingSide
35#include "lib_grid/grid/grid_util.h" // for CompareVertices
36
37
38namespace ug {
39
40
41template <typename TAAPos>
43(
44 Edge* elem,
45 const TAAPos& aaPos,
46 number thresholdRatio
47)
48{
49 return ISOTROPIC;
50}
51
52
53
54template <typename TAAPos>
56(
58 const TAAPos& aaPos,
59 number thresholdRatio
60)
61{
62 // check whether elem is anisotropic
63 number sideLength02 = VertexDistance(q->vertex(0), q->vertex(1), aaPos)
64 + VertexDistance(q->vertex(2), q->vertex(3), aaPos);
65
66 number sideLength13 = VertexDistance(q->vertex(1), q->vertex(2), aaPos)
67 + VertexDistance(q->vertex(3), q->vertex(0), aaPos);
68
69 if (sideLength02 < thresholdRatio * sideLength13)
70 return QUAD_SHORTX;
71 if (sideLength13 < thresholdRatio * sideLength02)
72 return QUAD_SHORTY;
73
74 return ISOTROPIC;
75}
76
77
78
79template <typename TAAPos>
81(
82 Face* elem,
83 const TAAPos& aaPos,
84 number thresholdRatio
85)
86{
87 Quadrilateral* q = dynamic_cast<Quadrilateral*>(elem);
88 if (!q)
89 return ISOTROPIC;
90
91 return is_anisotropic(q, aaPos, thresholdRatio);
92}
93
94
95
96
97template <typename TAAPos>
99(
100 Prism* p,
101 const TAAPos& aaPos,
102 number thresholdRatio
103)
104{
105 // check whether elem is anisotropic
106 number length = VertexDistance(p->vertex(0), p->vertex(3), aaPos)
107 + VertexDistance(p->vertex(1), p->vertex(4), aaPos)
108 + VertexDistance(p->vertex(2), p->vertex(5), aaPos);
109 number width = VertexDistance(p->vertex(0), p->vertex(1), aaPos)
110 + VertexDistance(p->vertex(1), p->vertex(2), aaPos)
111 + VertexDistance(p->vertex(2), p->vertex(0), aaPos)
112 + VertexDistance(p->vertex(3), p->vertex(4), aaPos)
113 + VertexDistance(p->vertex(4), p->vertex(5), aaPos)
114 + VertexDistance(p->vertex(5), p->vertex(3), aaPos);
115
116 number ratio = width ? 2.0*length/width : std::numeric_limits<number>::max();
117
118 // flat case
119 if (ratio < thresholdRatio)
120 return PRISM_FLAT;
121
122 // long case
123 if (ratio*thresholdRatio > 1)
124 return PRISM_LONG;
125
126 return ISOTROPIC;
127}
128
129
130
131template <typename TAAPos>
133(
134 Hexahedron* hex,
135 const TAAPos& aaPos,
136 number thresholdRatio
137)
138{
139 number length1 = VertexDistance(hex->vertex(0), hex->vertex(1), aaPos)
140 + VertexDistance(hex->vertex(2), hex->vertex(3), aaPos)
141 + VertexDistance(hex->vertex(4), hex->vertex(5), aaPos)
142 + VertexDistance(hex->vertex(6), hex->vertex(7), aaPos);
143 number length2 = VertexDistance(hex->vertex(0), hex->vertex(3), aaPos)
144 + VertexDistance(hex->vertex(1), hex->vertex(2), aaPos)
145 + VertexDistance(hex->vertex(4), hex->vertex(7), aaPos)
146 + VertexDistance(hex->vertex(5), hex->vertex(6), aaPos);
147 number length3 = VertexDistance(hex->vertex(0), hex->vertex(4), aaPos)
148 + VertexDistance(hex->vertex(1), hex->vertex(5), aaPos)
149 + VertexDistance(hex->vertex(2), hex->vertex(6), aaPos)
150 + VertexDistance(hex->vertex(3), hex->vertex(7), aaPos);
151
152 bool shortx = false;
153 bool shorty = false;
154 bool shortz = false;
155
156 if (length1 < thresholdRatio * length2 || length1 < thresholdRatio * length3)
157 shortx = true;
158
159 if (length2 < thresholdRatio * length1 || length2 < thresholdRatio * length3)
160 shorty = true;
161
162 if (length3 < thresholdRatio * length1 || length3 < thresholdRatio * length2)
163 shortz = true;
164
165 if (shortx)
166 {
167 if (shorty)
168 return HEX_SHORTXY;
169 if (shortz)
170 return HEX_SHORTXZ;
171 return HEX_SHORTX;
172 }
173
174 if (shorty)
175 {
176 if (shortz)
177 return HEX_SHORTYZ;
178 return HEX_SHORTY;
179 }
180
181 if (shortz)
182 return HEX_SHORTZ;
183
184 return ISOTROPIC;
185}
186
187
188
189template <typename TAAPos>
191(
192 Volume* elem,
193 const TAAPos& aaPos,
194 number thresholdRatio
195)
196{
197 // treat prism case
198 Prism* prism = dynamic_cast<Prism*>(elem);
199 if (prism)
200 return is_anisotropic(prism, aaPos, thresholdRatio);
201
202 // treat hexahedron case
203 Hexahedron* hex = dynamic_cast<Hexahedron*>(elem);
204 if (hex)
205 return is_anisotropic(hex, aaPos, thresholdRatio);
206
207 // other cases do not exist
208 return ISOTROPIC;
209}
210
211
212
213template <typename TAAPos>
215(
216 Edge* elem,
217 Grid& grid,
218 const TAAPos& aaPos,
219 number thresholdRatio,
220 std::vector<Vertex*>& sidesOut
221)
222{
223 return ISOTROPIC;
224}
225
226
227template <typename TAAPos>
229(
230 Face* elem,
231 Grid& grid,
232 const TAAPos& aaPos,
233 number thresholdRatio,
234 std::vector<Edge*>& sidesOut
235)
236{
237 // check whether this is a quadrilateral
238 Quadrilateral* q = dynamic_cast<Quadrilateral*>(elem);
239 if (!q)
240 return ISOTROPIC;
241
242 // check whether element is anisotropic (and which case)
243 AnisotropyState state = is_anisotropic(q, aaPos, thresholdRatio);
244 if (state == ISOTROPIC)
245 return state;
246
247 if (state == QUAD_SHORTX)
248 {
250 grid.associated_elements_sorted(assEd, q);
251 sidesOut.push_back(assEd[1]);
252 sidesOut.push_back(assEd[3]);
253
254 return state;
255 }
256
257 if (state == QUAD_SHORTY)
258 {
260 grid.associated_elements_sorted(assEd, q);
261 sidesOut.push_back(assEd[0]);
262 sidesOut.push_back(assEd[2]);
263
264 return state;
265 }
266
267 return state;
268}
269
270
271
272template <typename TAAPos>
274(
275 Volume* elem,
276 Grid& grid,
277 const TAAPos& aaPos,
278 number thresholdRatio,
279 std::vector<Face*>& sidesOut
280)
281{
282 // treat prism case
283 Prism* prism = dynamic_cast<Prism*>(elem);
284 if (prism)
285 {
286 AnisotropyState state = is_anisotropic(prism, aaPos, thresholdRatio);
287 if (state == ISOTROPIC)
288 return state;
289
290 // flat case
291 if (state == PRISM_FLAT)
292 {
293 Face* side = grid.get_side(prism, 0);
294 if (side)
295 sidesOut.push_back(side);
296
297 side = grid.get_side(prism, 4);
298 if (side)
299 sidesOut.push_back(side);
300 }
301
302 // long case
303 if (state == PRISM_LONG)
304 {
305 Face* side = grid.get_side(prism, 1);
306 if (side)
307 sidesOut.push_back(side);
308
309 side = grid.get_side(prism, 2);
310 if (side)
311 sidesOut.push_back(side);
312
313 side = grid.get_side(prism, 3);
314 if (side)
315 sidesOut.push_back(side);
316 }
317
318 return state;
319 }
320
321
322 // treat heaxahedron case
323 Hexahedron* hex = dynamic_cast<Hexahedron*>(elem);
324 if (hex)
325 {
326 AnisotropyState state = is_anisotropic(hex, aaPos, thresholdRatio);
327 if (state == ISOTROPIC)
328 return state;
329
330
331 // short in x direction
332 if (state == HEX_SHORTX || state == HEX_SHORTXY || state == HEX_SHORTXZ)
333 {
334 Face* side = grid.get_side(hex, 2);
335 if (side)
336 sidesOut.push_back(side);
337
338 side = grid.get_side(hex, 4);
339 if (side)
340 sidesOut.push_back(side);
341 }
342
343 // short in y direction
344 if (state == HEX_SHORTY || state == HEX_SHORTXY || state == HEX_SHORTYZ)
345 {
346 Face* side = grid.get_side(hex, 1);
347 if (side)
348 sidesOut.push_back(side);
349
350 side = grid.get_side(hex, 3);
351 if (side)
352 sidesOut.push_back(side);
353 }
354
355 // short in z direction
356 if (state == HEX_SHORTZ || state == HEX_SHORTXZ || state == HEX_SHORTYZ)
357 {
358 Face* side = grid.get_side(hex, 0);
359 if (side)
360 sidesOut.push_back(side);
361
362 side = grid.get_side(hex, 5);
363 if (side)
364 sidesOut.push_back(side);
365 }
366
367 return state;
368 }
369
370
371 // other elements cannot be anisotropic
372 return ISOTROPIC;
373}
374
375
376
377
378
379template <typename TAAPos>
381(
382 Edge* elem,
383 Grid& grid,
384 const TAAPos& aaPos,
385 number thresholdRatio,
386 std::vector<Edge*>& longEdges
387)
388{
389 return ISOTROPIC;
390}
391
392
393template <typename TAAPos>
395(
396 Face* elem,
397 Grid& grid,
398 const TAAPos& aaPos,
399 number thresholdRatio,
400 std::vector<Edge*>& longEdges
401)
402{
403 return close_sides_of_anisotropic_elem(elem, grid, aaPos, thresholdRatio, longEdges);
404}
405
406
407template <typename TAAPos>
409(
410 Volume* elem,
411 Grid& grid,
412 const TAAPos& aaPos,
413 number thresholdRatio,
414 std::vector<Edge*>& longEdges
415)
416{
417 // treat prism case
418 Prism* prism = dynamic_cast<Prism*>(elem);
419 if (prism)
420 {
421 AnisotropyState state = is_anisotropic(prism, aaPos, thresholdRatio);
422 if (state == ISOTROPIC)
423 return state;
424
425 // flat case
426 if (state == PRISM_FLAT)
427 {
429 grid.associated_elements_sorted(assEd, prism);
430
431 UG_COND_THROW(assEd.size() < 9, "Prism needs to have 9 edges, but only "
432 << assEd.size() << " found.");
433
434 longEdges.reserve(longEdges.size() + 6);
435 longEdges.push_back(assEd[0]);
436 longEdges.push_back(assEd[1]);
437 longEdges.push_back(assEd[2]);
438 longEdges.push_back(assEd[6]);
439 longEdges.push_back(assEd[7]);
440 longEdges.push_back(assEd[8]);
441
442 return state;
443 }
444
445 // long case
446 if (state == PRISM_LONG)
447 {
449 grid.associated_elements_sorted(assEd, prism);
450
451 UG_COND_THROW(assEd.size() < 9, "Prism needs to have 9 edges, but only "
452 << assEd.size() << " found.");
453
454 longEdges.reserve(longEdges.size() + 3);
455 longEdges.push_back(assEd[3]);
456 longEdges.push_back(assEd[4]);
457 longEdges.push_back(assEd[5]);
458
459 return state;
460 }
461
462 return state;
463 }
464
465
466 // treat heaxahedron case
467 Hexahedron* hex = dynamic_cast<Hexahedron*>(elem);
468 if (hex)
469 {
470 AnisotropyState state = is_anisotropic(hex, aaPos, thresholdRatio);
471 if (state == ISOTROPIC)
472 return state;
473
474
475 // short in x direction
476 if (state == HEX_SHORTX)
477 {
479 grid.associated_elements_sorted(assEd, hex);
480
481 UG_COND_THROW(assEd.size() < 12, "Hexahedron needs to have 12 edges, but only "
482 << assEd.size() << " found.");
483
484 longEdges.reserve(longEdges.size() + 8);
485 longEdges.push_back(assEd[1]);
486 longEdges.push_back(assEd[3]);
487 longEdges.push_back(assEd[4]);
488 longEdges.push_back(assEd[5]);
489 longEdges.push_back(assEd[6]);
490 longEdges.push_back(assEd[7]);
491 longEdges.push_back(assEd[8]);
492 longEdges.push_back(assEd[11]);
493
494 return state;
495 }
496
497 // short in y direction
498 if (state == HEX_SHORTY)
499 {
501 grid.associated_elements_sorted(assEd, hex);
502
503 UG_COND_THROW(assEd.size() < 12, "Hexahedron needs to have 12 edges, but only "
504 << assEd.size() << " found.");
505
506 longEdges.reserve(longEdges.size() + 8);
507 longEdges.push_back(assEd[0]);
508 longEdges.push_back(assEd[2]);
509 longEdges.push_back(assEd[4]);
510 longEdges.push_back(assEd[5]);
511 longEdges.push_back(assEd[6]);
512 longEdges.push_back(assEd[7]);
513 longEdges.push_back(assEd[8]);
514 longEdges.push_back(assEd[10]);
515
516 return state;
517 }
518
519 // short in z direction
520 if (state == HEX_SHORTZ)
521 {
523 grid.associated_elements_sorted(assEd, hex);
524
525 UG_COND_THROW(assEd.size() < 12, "Hexahedron needs to have 12 edges, but only "
526 << assEd.size() << " found.");
527
528 longEdges.reserve(longEdges.size() + 8);
529 longEdges.push_back(assEd[0]);
530 longEdges.push_back(assEd[1]);
531 longEdges.push_back(assEd[2]);
532 longEdges.push_back(assEd[3]);
533 longEdges.push_back(assEd[8]);
534 longEdges.push_back(assEd[9]);
535 longEdges.push_back(assEd[10]);
536 longEdges.push_back(assEd[11]);
537
538 return state;
539 }
540
541 // short in x and y direction
542 if (state == HEX_SHORTXY)
543 {
545 grid.associated_elements_sorted(assEd, hex);
546
547 UG_COND_THROW(assEd.size() < 12, "Hexahedron needs to have 12 edges, but only "
548 << assEd.size() << " found.");
549
550 longEdges.reserve(longEdges.size() + 4);
551 longEdges.push_back(assEd[4]);
552 longEdges.push_back(assEd[5]);
553 longEdges.push_back(assEd[6]);
554 longEdges.push_back(assEd[7]);
555
556 return state;
557 }
558
559 // short in x and z direction
560 if (state == HEX_SHORTXZ)
561 {
563 grid.associated_elements_sorted(assEd, hex);
564
565 UG_COND_THROW(assEd.size() < 12, "Hexahedron needs to have 12 edges, but only "
566 << assEd.size() << " found.");
567
568 longEdges.reserve(longEdges.size() + 4);
569 longEdges.push_back(assEd[1]);
570 longEdges.push_back(assEd[3]);
571 longEdges.push_back(assEd[9]);
572 longEdges.push_back(assEd[11]);
573
574 return state;
575 }
576
577 // short in y and z direction
578 if (state == HEX_SHORTYZ)
579 {
581 grid.associated_elements_sorted(assEd, hex);
582
583 UG_COND_THROW(assEd.size() < 12, "Hexahedron needs to have 12 edges, but only "
584 << assEd.size() << " found.");
585
586 longEdges.reserve(longEdges.size() + 4);
587 longEdges.push_back(assEd[0]);
588 longEdges.push_back(assEd[2]);
589 longEdges.push_back(assEd[8]);
590 longEdges.push_back(assEd[10]);
591
592 return state;
593 }
594
595 return state;
596 }
597
598
599 // other elements cannot be anisotropic
600 return ISOTROPIC;
601}
602
603} // namespace ug
604
parameterString p
virtual Vertex * vertex(size_t index) const
Definition grid_objects_2d.h:265
Base-class for edges.
Definition grid_base_objects.h:397
Faces are 2-dimensional objects.
Definition grid_base_objects.h:510
Manages the elements of a grid and their interconnection.
Definition grid.h:132
Vertex::side * get_side(Vertex *obj, size_t side)
This method returns the i-th side of an Edge, Face or Volume.
Definition grid.cpp:1170
void associated_elements_sorted(traits< Vertex >::secure_container &elemsOut, TElem *e)
Puts all elements of type TAss which are contained in 'e' into elemsOut in the reference elements ord...
Definition grid_impl.hpp:503
traits< Edge >::secure_container SecureEdgeContainer
Container to store associated edges.
Definition grid.h:164
A volume element with 6 quadrilateral sides.
Definition grid_objects_3d.h:227
virtual Vertex * vertex(size_t index) const
Definition grid_objects_3d.h:243
A volume element with 2 triangle and 3 quadrilateral sides.
Definition grid_objects_3d.h:360
a face with four points.
Definition grid_objects_2d.h:323
Volumes are 3-dimensional objects.
Definition grid_base_objects.h:754
UG_API number VertexDistance(Vertex *v0, Vertex *v1, TAAPos &aaPos)
Returns the distance between two vertices.
Definition vertex_util_impl.hpp:54
#define UG_COND_THROW(cond, msg)
UG_COND_THROW(cond, msg) : performs a UG_THROW(msg) if cond == true.
Definition error.h:61
double number
Definition types.h:124
the ug namespace
AnisotropyState is_anisotropic(Edge *elem, const TAAPos &aaPos, number thresholdRatio)
Definition anisotropy_util_impl.h:43
AnisotropyState
Definition anisotropy_util.h:48
@ HEX_SHORTY
Definition anisotropy_util.h:55
@ QUAD_SHORTX
Definition anisotropy_util.h:50
@ HEX_SHORTXZ
Definition anisotropy_util.h:58
@ QUAD_SHORTY
Definition anisotropy_util.h:51
@ HEX_SHORTX
Definition anisotropy_util.h:54
@ PRISM_LONG
Definition anisotropy_util.h:53
@ HEX_SHORTXY
Definition anisotropy_util.h:57
@ PRISM_FLAT
Definition anisotropy_util.h:52
@ HEX_SHORTYZ
Definition anisotropy_util.h:59
@ ISOTROPIC
Definition anisotropy_util.h:49
@ HEX_SHORTZ
Definition anisotropy_util.h:56
AnisotropyState long_edges_of_anisotropic_elem(Edge *elem, Grid &grid, const TAAPos &aaPos, number thresholdRatio, std::vector< Edge * > &longEdges)
Definition anisotropy_util_impl.h:381
AnisotropyState close_sides_of_anisotropic_elem(Edge *elem, Grid &grid, const TAAPos &aaPos, number thresholdRatio, std::vector< Vertex * > &sidesOut)
Definition anisotropy_util_impl.h:215