libMesh
elem_quality.C
Go to the documentation of this file.
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 // C++ includes
19 #include <iostream>
20 #include <sstream>
21 
22 // Local includes
23 #include "libmesh/libmesh_common.h"
24 #include "libmesh/elem_quality.h"
25 #include "libmesh/enum_elem_type.h"
26 #include "libmesh/enum_elem_quality.h"
27 
28 
29 namespace libMesh
30 {
31 
32 // ------------------------------------------------------------
33 // Quality function definitions
34 
42 std::string Quality::name (const ElemQuality q)
43 {
44  std::string its_name;
45 
46  switch (q)
47  {
48 
49  case ASPECT_RATIO:
50  its_name = "Aspect Ratio";
51  break;
52 
53  case SKEW:
54  its_name = "Skew";
55  break;
56 
57  case SHEAR:
58  its_name = "Shear";
59  break;
60 
61  case SHAPE:
62  its_name = "Shape";
63  break;
64 
65  case MAX_ANGLE:
66  its_name = "Maximum Angle";
67  break;
68 
69  case MIN_ANGLE:
70  its_name = "Minimum Angle";
71  break;
72 
73  case MAX_DIHEDRAL_ANGLE:
74  its_name = "Maximum Dihedral Angle";
75  break;
76 
77  case MIN_DIHEDRAL_ANGLE:
78  its_name = "Minimum Dihedral Angle";
79  break;
80 
81  case CONDITION:
82  its_name = "Condition Number";
83  break;
84 
85  case DISTORTION:
86  its_name = "Distortion";
87  break;
88 
89  case TAPER:
90  its_name = "Taper";
91  break;
92 
93  case WARP:
94  its_name = "Warp";
95  break;
96 
97  case STRETCH:
98  its_name = "Stretch";
99  break;
100 
101  case DIAGONAL:
102  its_name = "Diagonal";
103  break;
104 
105  case ASPECT_RATIO_BETA:
106  its_name = "AR Beta";
107  break;
108 
109  case ASPECT_RATIO_GAMMA:
110  its_name = "AR Gamma";
111  break;
112 
113  case SIZE:
114  its_name = "Size";
115  break;
116 
117  case JACOBIAN:
118  its_name = "Jacobian";
119  break;
120 
121  case SCALED_JACOBIAN:
122  its_name = "Scaled Jacobian";
123  break;
124 
125  default:
126  its_name = "Unknown";
127  break;
128  }
129 
130  return its_name;
131 }
132 
133 
134 
135 
136 
142 std::string Quality::describe (const ElemQuality q)
143 {
144 
145  std::ostringstream desc;
146 
147  switch (q)
148  {
149 
150  case EDGE_LENGTH_RATIO:
151  case ASPECT_RATIO:
152  desc << "Max edge length ratio\n"
153  << "at element center.\n"
154  << '\n'
155  << "Suggested ranges:\n"
156  << "Hexes: (1 -> 4)\n"
157  << "Quads: (1 -> 4)";
158  break;
159 
160  case SKEW:
161  desc << "Maximum |cos A|, where A\n"
162  << "is the angle between edges\n"
163  << "at element center.\n"
164  << '\n'
165  << "Suggested ranges:\n"
166  << "Hexes: (0 -> 0.5)\n"
167  << "Quads: (0 -> 0.5)";
168  break;
169 
170  case SHEAR:
171  desc << "LIBMESH_DIM / K(Js)\n"
172  << '\n'
173  << "LIBMESH_DIM = element dimension.\n"
174  << "K(Js) = Condition number of \n"
175  << " Jacobian skew matrix.\n"
176  << '\n'
177  << "Suggested ranges:\n"
178  << "Hexes(LIBMESH_DIM=3): (0.3 -> 1)\n"
179  << "Quads(LIBMESH_DIM=2): (0.3 -> 1)";
180  break;
181 
182  case SHAPE:
183  desc << "LIBMESH_DIM / K(Jw)\n"
184  << '\n'
185  << "LIBMESH_DIM = element dimension.\n"
186  << "K(Jw) = Condition number of \n"
187  << " weighted Jacobian\n"
188  << " matrix.\n"
189  << '\n'
190  << "Suggested ranges:\n"
191  << "Hexes(LIBMESH_DIM=3): (0.3 -> 1)\n"
192  << "Tets(LIBMESH_DIM=3): (0.2 -> 1)\n"
193  << "Quads(LIBMESH_DIM=2): (0.3 -> 1).";
194  break;
195 
196  case MAX_ANGLE:
197  desc << "Largest angle between all adjacent pairs of edges (in 2D, sides).\n"
198  << '\n'
199  << "Suggested ranges:\n"
200  << "Quads: (90 -> 135)\n"
201  << "Triangles: (60 -> 90)";
202  break;
203 
204  case MIN_ANGLE:
205  desc << "Smallest angle between all adjacent pairs of edges (in 2D, sides).\n"
206  << '\n'
207  << "Suggested ranges:\n"
208  << "Quads: (45 -> 90)\n"
209  << "Triangles: (30 -> 60)";
210  break;
211 
212  case MAX_DIHEDRAL_ANGLE:
213  desc << "Largest angle between all adjacent pairs of sides (in 2D, equivalent to MAX_ANGLE).\n"
214  << '\n'
215  << "Suggested ranges:\n"
216  << "Quads: (90 -> 135)\n"
217  << "Triangles: (60 -> 90)";
218  break;
219 
220  case MIN_DIHEDRAL_ANGLE:
221  desc << "Smallest angle between all adjacent pairs of sides (in 2D, equivalent to MIN_ANGLE).\n"
222  << '\n'
223  << "Suggested ranges:\n"
224  << "Quads: (45 -> 90)\n"
225  << "Triangles: (30 -> 60)";
226  break;
227 
228  case CONDITION:
229  desc << "Condition number of the\n"
230  << "Jacobian matrix.\n"
231  << '\n'
232  << "Suggested ranges:\n"
233  << "Quads: (1 -> 4)\n"
234  << "Hexes: (1 -> 8)\n"
235  << "Tris: (1 -> 1.3)\n"
236  << "Tets: (1 -> 3)";
237  break;
238 
239  case DISTORTION:
240  desc << "min |J| * A / <A>\n"
241  << '\n'
242  << "|J| = norm of Jacobian matrix\n"
243  << " A = actual area\n"
244  << "<A> = reference area\n"
245  << '\n'
246  << "Suggested ranges:\n"
247  << "Quads: (0.6 -> 1), <A>=4\n"
248  << "Hexes: (0.6 -> 1), <A>=8\n"
249  << "Tris: (0.6 -> 1), <A>=1/2\n"
250  << "Tets: (0.6 -> 1), <A>=1/6";
251  break;
252 
253  case TAPER:
254  desc << "Maximum ratio of lengths\n"
255  << "derived from opposite edges.\n"
256  << '\n'
257  << "Suggested ranges:\n"
258  << "Quads: (0.7 -> 1)\n"
259  << "Hexes: (0.4 -> 1)";
260  break;
261 
262  case WARP:
263  desc << "cos D\n"
264  << '\n'
265  << "D = minimum dihedral angle\n"
266  << " formed by diagonals.\n"
267  << '\n'
268  << "Suggested ranges:\n"
269  << "Quads: (0.9 -> 1)";
270  break;
271 
272  case STRETCH:
273  desc << "Sqrt(3) * L_min / L_max\n"
274  << '\n'
275  << "L_min = minimum edge length.\n"
276  << "L_max = maximum edge length.\n"
277  << '\n'
278  << "Suggested ranges:\n"
279  << "Quads: (0.25 -> 1)\n"
280  << "Hexes: (0.25 -> 1)";
281  break;
282 
283  case DIAGONAL:
284  desc << "D_min / D_max\n"
285  << '\n'
286  << "D_min = minimum diagonal.\n"
287  << "D_max = maximum diagonal.\n"
288  << '\n'
289  << "Suggested ranges:\n"
290  << "Hexes: (0.65 -> 1)";
291  break;
292 
293  case ASPECT_RATIO_BETA:
294  desc << "CR / (3 * IR)\n"
295  << '\n'
296  << "CR = circumsphere radius\n"
297  << "IR = inscribed sphere radius\n"
298  << '\n'
299  << "Suggested ranges:\n"
300  << "Tets: (1 -> 3)";
301  break;
302 
303  case ASPECT_RATIO_GAMMA:
304  desc << "S^(3/2) / 8.479670 * V\n"
305  << '\n'
306  << "S = sum(si*si/6)\n"
307  << "si = edge length\n"
308  << "V = volume\n"
309  << '\n'
310  << "Suggested ranges:\n"
311  << "Tets: (1 -> 3)";
312  break;
313 
314  case SIZE:
315  desc << "min (|J|, |1/J|)\n"
316  << '\n'
317  << "|J| = norm of Jacobian matrix.\n"
318  << '\n'
319  << "Suggested ranges:\n"
320  << "Quads: (0.3 -> 1)\n"
321  << "Hexes: (0.5 -> 1)\n"
322  << "Tris: (0.25 -> 1)\n"
323  << "Tets: (0.2 -> 1)";
324  break;
325 
326  case JACOBIAN:
327  case SCALED_JACOBIAN:
328  desc << "Minimum nodal Jacobian.\n"
329  << "The nodal Jacobians are computed by taking the cross product (2D) or scalar product (3D) of the adjacent edges that meet at that node.\n"
330  << "In the SCALED_JACOBIAN case, we also then divide by the lengths of each of the associated edges.\n"
331  << "For Pyramid elements where four edges meet at the apex node, special handling is required.\n"
332  << '\n'
333  << "Suggested acceptable ranges (from Cubit documentation) for SCALED_JACOBIAN metric:\n"
334  << "Quads/Hexes: (0.5 -> 1)\n"
335  << "Tris/Tets: (0.2 -> 1.0)";
336  break;
337 
338  default:
339  desc << "Unknown";
340  break;
341  }
342 
343  return desc.str();
344 }
345 
346 
347 std::vector<ElemQuality> Quality::valid(const ElemType t)
348 {
349  std::vector<ElemQuality> v;
350 
351  switch (t)
352  {
353  case EDGE2:
354  case EDGE3:
355  case EDGE4:
356  {
357  // None yet
358  break;
359  }
360 
361  case TRI3:
362  case TRISHELL3:
363  case TRI6:
364  case TRI7:
365  {
366  v = {
367  CONDITION,
368  DISTORTION,
369  EDGE_LENGTH_RATIO,
370  JACOBIAN,
371  SCALED_JACOBIAN,
372  MAX_ANGLE,
373  MIN_ANGLE,
374  MAX_DIHEDRAL_ANGLE,
375  MIN_DIHEDRAL_ANGLE,
376  SHAPE,
377  SIZE
378  };
379 
380  break;
381  }
382 
383  case QUAD4:
384  case QUADSHELL4:
385  case QUAD8:
386  case QUADSHELL8:
387  case QUAD9:
388  case QUADSHELL9:
389  {
390  v = {
391  ASPECT_RATIO,
392  CONDITION,
393  DISTORTION,
394  EDGE_LENGTH_RATIO,
395  JACOBIAN,
396  SCALED_JACOBIAN,
397  MAX_ANGLE,
398  MIN_ANGLE,
399  MAX_DIHEDRAL_ANGLE,
400  MIN_DIHEDRAL_ANGLE,
401  SHAPE,
402  SHEAR,
403  SIZE,
404  SKEW,
405  STRETCH,
406  TAPER,
407  WARP
408  };
409 
410  break;
411  }
412 
413  case TET4:
414  case TET10:
415  case TET14:
416  {
417  v = {
418  ASPECT_RATIO_BETA,
419  ASPECT_RATIO_GAMMA,
420  CONDITION,
421  DISTORTION,
422  JACOBIAN,
423  SCALED_JACOBIAN,
424  MAX_ANGLE,
425  MIN_ANGLE,
426  MAX_DIHEDRAL_ANGLE,
427  MIN_DIHEDRAL_ANGLE,
428  SHAPE,
429  SIZE
430  };
431 
432  break;
433  }
434 
435  case HEX8:
436  case HEX20:
437  case HEX27:
438  {
439  v = {
440  ASPECT_RATIO,
441  CONDITION,
442  DIAGONAL,
443  DISTORTION,
444  JACOBIAN,
445  SCALED_JACOBIAN,
446  MAX_ANGLE,
447  MIN_ANGLE,
448  MAX_DIHEDRAL_ANGLE,
449  MIN_DIHEDRAL_ANGLE,
450  SHAPE,
451  SHEAR,
452  SIZE,
453  SKEW,
454  STRETCH,
455  TAPER
456  };
457 
458  break;
459  }
460 
461  case PRISM6:
462  case PRISM18:
463  case PRISM20:
464  case PRISM21:
465  {
466  v = {
467  EDGE_LENGTH_RATIO,
468  MAX_ANGLE,
469  MIN_ANGLE,
470  MAX_DIHEDRAL_ANGLE,
471  MIN_DIHEDRAL_ANGLE,
472  };
473 
474  break;
475  }
476 
477  case PYRAMID5:
478  case PYRAMID13:
479  case PYRAMID14:
480  case PYRAMID18:
481  {
482  v = {
483  EDGE_LENGTH_RATIO,
484  MAX_ANGLE,
485  MIN_ANGLE,
486  MAX_DIHEDRAL_ANGLE,
487  MIN_DIHEDRAL_ANGLE,
488  };
489 
490  break;
491  }
492 
493 
494 
495 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
496 
497  case INFEDGE2:
498  {
499  // None yet
500  break;
501  }
502 
503  case INFQUAD4:
504  case INFQUAD6:
505  case INFHEX8:
506  case INFHEX16:
507  case INFHEX18:
508  case INFPRISM6:
509  case INFPRISM12:
510  {
511  v = {
512  MAX_ANGLE,
513  MIN_ANGLE,
514  MAX_DIHEDRAL_ANGLE,
515  MIN_DIHEDRAL_ANGLE,
516  };
517 
518  break;
519  }
520 
521 #endif
522 
523 
524  default:
525  libmesh_error_msg("Undefined element type!");
526  }
527 
528  return v;
529 }
530 
531 } // namespace libMesh
libMesh::Quality::name
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
libMesh::ElemType
ElemType
Defines an enum for geometric element types.
Definition: enum_elem_type.h:33
libMesh::Quality::describe
std::string describe(const ElemQuality q)
This function returns a string containing a short description of q.
Definition: elem_quality.C:142
libMesh
The libMesh namespace provides an interface to certain functionality in the library.
libMesh::ElemQuality
ElemQuality
Defines an enum for element quality metrics.
Definition: enum_elem_quality.h:34
libMesh::Quality::valid
std::vector< ElemQuality > valid(const ElemType t)
Definition: elem_quality.C:347
Generated on Sat Aug 2 2025 20:01:50 for libMesh by   doxygen 1.8.14