libMesh::Quality Namespace Reference

Utility functions for computing element quality indicators. More...

Functions

std::string name (const ElemQuality q)
 
std::string describe (const ElemQuality q)
 
std::vector< ElemQualityvalid (const ElemType t)
 

Variables

const unsigned int num_quals = 16
 

Detailed Description

Utility functions for computing element quality indicators.

A namespace for quality utility functions.

Author
John W. Peterson
Date
2002

Function Documentation

◆ describe()

std::string libMesh::Quality::describe ( const ElemQuality  q)
Returns
A description for a ElemQuality enum

This function returns a string containing a short description of q. Useful for asking the enum what it computes.

Definition at line 130 of file elem_quality.C.

References libMesh::ASPECT_RATIO, libMesh::ASPECT_RATIO_BETA, libMesh::ASPECT_RATIO_GAMMA, libMesh::CONDITION, libMesh::DIAGONAL, libMesh::DISTORTION, libMesh::JACOBIAN, libMesh::MAX_ANGLE, libMesh::MIN_ANGLE, libMesh::SHAPE, libMesh::SHEAR, libMesh::SIZE, libMesh::SKEW, libMesh::STRETCH, libMesh::TAPER, and libMesh::WARP.

131 {
132 
133  std::ostringstream desc;
134 
135  switch (q)
136  {
137 
138  case ASPECT_RATIO:
139  desc << "Max edge length ratio\n"
140  << "at element center.\n"
141  << '\n'
142  << "Suggested ranges:\n"
143  << "Hexes: (1 -> 4)\n"
144  << "Quads: (1 -> 4)";
145  break;
146 
147  case SKEW:
148  desc << "Maximum |cos A|, where A\n"
149  << "is the angle between edges\n"
150  << "at element center.\n"
151  << '\n'
152  << "Suggested ranges:\n"
153  << "Hexes: (0 -> 0.5)\n"
154  << "Quads: (0 -> 0.5)";
155  break;
156 
157  case SHEAR:
158  desc << "LIBMESH_DIM / K(Js)\n"
159  << '\n'
160  << "LIBMESH_DIM = element dimension.\n"
161  << "K(Js) = Condition number of \n"
162  << " Jacobian skew matrix.\n"
163  << '\n'
164  << "Suggested ranges:\n"
165  << "Hexes(LIBMESH_DIM=3): (0.3 -> 1)\n"
166  << "Quads(LIBMESH_DIM=2): (0.3 -> 1)";
167  break;
168 
169  case SHAPE:
170  desc << "LIBMESH_DIM / K(Jw)\n"
171  << '\n'
172  << "LIBMESH_DIM = element dimension.\n"
173  << "K(Jw) = Condition number of \n"
174  << " weighted Jacobian\n"
175  << " matrix.\n"
176  << '\n'
177  << "Suggested ranges:\n"
178  << "Hexes(LIBMESH_DIM=3): (0.3 -> 1)\n"
179  << "Tets(LIBMESH_DIM=3): (0.2 -> 1)\n"
180  << "Quads(LIBMESH_DIM=2): (0.3 -> 1).";
181  break;
182 
183  case MAX_ANGLE:
184  desc << "Largest included angle.\n"
185  << '\n'
186  << "Suggested ranges:\n"
187  << "Quads: (90 -> 135)\n"
188  << "Triangles: (60 -> 90)";
189  break;
190 
191  case MIN_ANGLE:
192  desc << "Smallest included angle.\n"
193  << '\n'
194  << "Suggested ranges:\n"
195  << "Quads: (45 -> 90)\n"
196  << "Triangles: (30 -> 60)";
197  break;
198 
199  case CONDITION:
200  desc << "Condition number of the\n"
201  << "Jacobian matrix.\n"
202  << '\n'
203  << "Suggested ranges:\n"
204  << "Quads: (1 -> 4)\n"
205  << "Hexes: (1 -> 8)\n"
206  << "Tris: (1 -> 1.3)\n"
207  << "Tets: (1 -> 3)";
208  break;
209 
210  case DISTORTION:
211  desc << "min |J| * A / <A>\n"
212  << '\n'
213  << "|J| = norm of Jacobian matrix\n"
214  << " A = actual area\n"
215  << "<A> = reference area\n"
216  << '\n'
217  << "Suggested ranges:\n"
218  << "Quads: (0.6 -> 1), <A>=4\n"
219  << "Hexes: (0.6 -> 1), <A>=8\n"
220  << "Tris: (0.6 -> 1), <A>=1/2\n"
221  << "Tets: (0.6 -> 1), <A>=1/6";
222  break;
223 
224  case TAPER:
225  desc << "Maximum ratio of lengths\n"
226  << "derived from opposite edges.\n"
227  << '\n'
228  << "Suggested ranges:\n"
229  << "Quads: (0.7 -> 1)\n"
230  << "Hexes: (0.4 -> 1)";
231  break;
232 
233  case WARP:
234  desc << "cos D\n"
235  << '\n'
236  << "D = minimum dihedral angle\n"
237  << " formed by diagonals.\n"
238  << '\n'
239  << "Suggested ranges:\n"
240  << "Quads: (0.9 -> 1)";
241  break;
242 
243  case STRETCH:
244  desc << "Sqrt(3) * L_min / L_max\n"
245  << '\n'
246  << "L_min = minimum edge length.\n"
247  << "L_max = maximum edge length.\n"
248  << '\n'
249  << "Suggested ranges:\n"
250  << "Quads: (0.25 -> 1)\n"
251  << "Hexes: (0.25 -> 1)";
252  break;
253 
254  case DIAGONAL:
255  desc << "D_min / D_max\n"
256  << '\n'
257  << "D_min = minimum diagonal.\n"
258  << "D_max = maximum diagonal.\n"
259  << '\n'
260  << "Suggested ranges:\n"
261  << "Hexes: (0.65 -> 1)";
262  break;
263 
264  case ASPECT_RATIO_BETA:
265  desc << "CR / (3 * IR)\n"
266  << '\n'
267  << "CR = circumsphere radius\n"
268  << "IR = inscribed sphere radius\n"
269  << '\n'
270  << "Suggested ranges:\n"
271  << "Tets: (1 -> 3)";
272  break;
273 
274  case ASPECT_RATIO_GAMMA:
275  desc << "S^(3/2) / 8.479670 * V\n"
276  << '\n'
277  << "S = sum(si*si/6)\n"
278  << "si = edge length\n"
279  << "V = volume\n"
280  << '\n'
281  << "Suggested ranges:\n"
282  << "Tets: (1 -> 3)";
283  break;
284 
285  case SIZE:
286  desc << "min (|J|, |1/J|)\n"
287  << '\n'
288  << "|J| = norm of Jacobian matrix.\n"
289  << '\n'
290  << "Suggested ranges:\n"
291  << "Quads: (0.3 -> 1)\n"
292  << "Hexes: (0.5 -> 1)\n"
293  << "Tris: (0.25 -> 1)\n"
294  << "Tets: (0.2 -> 1)";
295  break;
296 
297  case JACOBIAN:
298  desc << "Minimum Jacobian divided by\n"
299  << "the lengths of the LIBMESH_DIM\n"
300  << "largest edge vectors.\n"
301  << '\n'
302  << "LIBMESH_DIM = element dimension.\n"
303  << '\n'
304  << "Suggested ranges:\n"
305  << "Quads: (0.5 -> 1)\n"
306  << "Hexes: (0.5 -> 1)\n"
307  << "Tris: (0.5 -> 1.155)\n"
308  << "Tets: (0.5 -> 1.414)";
309  break;
310 
311  default:
312  desc << "Unknown";
313  break;
314  }
315 
316  return desc.str();
317 }

◆ name()

std::string libMesh::Quality::name ( const ElemQuality  q)
Returns
A descriptive name for a ElemQuality enum

This function returns a string containing some name for q. Useful for asking the enum what its name is. I added this since you may want a simple way to attach a name or description to the ElemQuality enums. It can be removed if it is found to be useless.

Definition at line 42 of file elem_quality.C.

References libMesh::ASPECT_RATIO, libMesh::ASPECT_RATIO_BETA, libMesh::ASPECT_RATIO_GAMMA, libMesh::CONDITION, libMesh::DIAGONAL, libMesh::DISTORTION, libMesh::JACOBIAN, libMesh::MAX_ANGLE, libMesh::MIN_ANGLE, libMesh::SHAPE, libMesh::SHEAR, libMesh::SIZE, libMesh::SKEW, libMesh::STRETCH, libMesh::TAPER, and libMesh::WARP.

Referenced by GETPOT_NAMESPACE::GetPot::_convert_to_type_no_default(), libMesh::EquationSystems::add_system(), libMesh::Factory< Base >::build(), libMesh::cast_ptr(), libMesh::cast_ref(), libMesh::ExodusII_IO_Helper::check_existing_vars(), libMesh::command_line_next(), libMesh::command_line_value(), libMesh::command_line_vector(), libMesh::Utility::complex_filename(), libMesh::ElemCutter::cut_3D(), libMesh::EquationSystems::delete_system(), libMesh::demangle(), DMlibMeshSetUpName_Private(), DMView_libMesh(), libMesh::Factory< Base >::Factory(), libMesh::Parameters::get(), libMesh::BoundaryInfo::get_id_by_name(), libMesh::MeshBase::get_id_by_name(), libMesh::ReferenceCounter::get_info(), libMesh::EquationSystems::get_system(), libMesh::EquationSystems::has_system(), libMesh::Parameters::have_parameter(), libMesh::ReferenceCounter::increment_constructor_count(), libMesh::ReferenceCounter::increment_destructor_count(), libMesh::PetscNonlinearSolver< Number >::init(), libMesh::PetscLinearSolver< T >::init(), libMesh::Parameters::insert(), libMesh::NameBasedIO::is_parallel_file_format(), libMesh::Xdr::open(), GETPOT_NAMESPACE::GetPot::variable::operator=(), libMesh::PetscMatrix< T >::print_matlab(), libMesh::PetscVector< T >::print_matlab(), libMesh::ExodusII_IO_Helper::NamesData::push_back_entry(), libMesh::OFFIO::read(), libMesh::NameBasedIO::read(), libMesh::TetGenIO::read(), libMesh::PltLoader::read(), libMesh::GmshIO::read(), libMesh::GMVIO::read(), libMesh::UnstructuredMesh::read(), libMesh::MatlabIO::read(), libMesh::VTKIO::read(), libMesh::EquationSystems::read(), libMesh::CheckpointIO::read_header(), libMesh::PltLoader::read_header(), libMesh::VariationalMeshSmoother::readmetr(), libMesh::ReferenceCountedObject< RBParametrized >::ReferenceCountedObject(), libMesh::Parameters::remove(), libMesh::Parameters::set(), libMesh::Parameters::Parameter< T >::type(), libMesh::NameBasedIO::write(), libMesh::GmshIO::write(), libMesh::UnstructuredMesh::write(), libMesh::CheckpointIO::write(), libMesh::EnsightIO::write(), libMesh::EquationSystems::write(), libMesh::TecplotIO::write_binary(), libMesh::PltLoader::write_dat(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::NameBasedIO::write_nodal_data(), libMesh::UCDIO::write_soln(), and libMesh::ReferenceCountedObject< RBParametrized >::~ReferenceCountedObject().

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 CONDITION:
74  its_name = "Condition Number";
75  break;
76 
77  case DISTORTION:
78  its_name = "Distortion";
79  break;
80 
81  case TAPER:
82  its_name = "Taper";
83  break;
84 
85  case WARP:
86  its_name = "Warp";
87  break;
88 
89  case STRETCH:
90  its_name = "Stretch";
91  break;
92 
93  case DIAGONAL:
94  its_name = "Diagonal";
95  break;
96 
97  case ASPECT_RATIO_BETA:
98  its_name = "AR Beta";
99  break;
100 
101  case ASPECT_RATIO_GAMMA:
102  its_name = "AR Gamma";
103  break;
104 
105  case SIZE:
106  its_name = "Size";
107  break;
108 
109  case JACOBIAN:
110  its_name = "Jacobian";
111  break;
112 
113  default:
114  its_name = "Unknown";
115  break;
116  }
117 
118  return its_name;
119 }

◆ valid()

std::vector< ElemQuality > libMesh::Quality::valid ( const ElemType  t)
Returns
The valid ElemQuality metrics for a given ElemType element type.

Returns all valid quality metrics for element type t.

Definition at line 324 of file elem_quality.C.

References libMesh::ASPECT_RATIO, libMesh::ASPECT_RATIO_BETA, libMesh::ASPECT_RATIO_GAMMA, libMesh::CONDITION, libMesh::DIAGONAL, libMesh::DISTORTION, libMesh::EDGE2, libMesh::EDGE3, libMesh::EDGE4, libMesh::HEX20, libMesh::HEX27, libMesh::HEX8, libMesh::INFEDGE2, libMesh::INFHEX16, libMesh::INFHEX18, libMesh::INFHEX8, libMesh::INFPRISM12, libMesh::INFPRISM6, libMesh::INFQUAD4, libMesh::INFQUAD6, libMesh::JACOBIAN, libMesh::MAX_ANGLE, libMesh::MIN_ANGLE, libMesh::PRISM18, libMesh::PRISM6, libMesh::PYRAMID13, libMesh::PYRAMID14, libMesh::PYRAMID5, libMesh::QUAD4, libMesh::QUAD8, libMesh::QUAD9, libMesh::QUADSHELL4, libMesh::QUADSHELL8, libMesh::SHAPE, libMesh::SHEAR, libMesh::SIZE, libMesh::SKEW, libMesh::STRETCH, libMesh::TAPER, libMesh::TET10, libMesh::TET4, libMesh::TRI3, libMesh::TRI6, libMesh::TRISHELL3, and libMesh::WARP.

325 {
326  std::vector<ElemQuality> v;
327 
328  switch (t)
329  {
330  case EDGE2:
331  case EDGE3:
332  case EDGE4:
333  {
334  // None yet
335  break;
336  }
337 
338  case TRI3:
339  case TRISHELL3:
340  case TRI6:
341  {
342  v.resize(7);
343  v[0] = MAX_ANGLE;
344  v[1] = MIN_ANGLE;
345  v[2] = CONDITION;
346  v[3] = JACOBIAN;
347  v[4] = SIZE;
348  v[5] = SHAPE;
349  v[6] = DISTORTION;
350  break;
351  }
352 
353  case QUAD4:
354  case QUADSHELL4:
355  case QUAD8:
356  case QUADSHELL8:
357  case QUAD9:
358  {
359  v.resize(13);
360  v[0] = ASPECT_RATIO;
361  v[1] = SKEW;
362  v[2] = TAPER;
363  v[3] = WARP;
364  v[4] = STRETCH;
365  v[5] = MIN_ANGLE;
366  v[6] = MAX_ANGLE;
367  v[7] = CONDITION;
368  v[8] = JACOBIAN;
369  v[9] = SHEAR;
370  v[10] = SHAPE;
371  v[11] = SIZE;
372  v[12] = DISTORTION;
373  break;
374  }
375 
376  case TET4:
377  case TET10:
378  {
379  v.resize(7);
380  v[0] = ASPECT_RATIO_BETA;
381  v[1] = ASPECT_RATIO_GAMMA;
382  v[2] = CONDITION;
383  v[3] = JACOBIAN;
384  v[4] = SHAPE;
385  v[5] = SIZE;
386  v[6] = DISTORTION;
387  break;
388  }
389 
390  case HEX8:
391  case HEX20:
392  case HEX27:
393  {
394  v.resize(11);
395  v[0] = ASPECT_RATIO;
396  v[1] = SKEW;
397  v[2] = SHEAR;
398  v[3] = SHAPE;
399  v[4] = CONDITION;
400  v[5] = JACOBIAN;
401  v[6] = DISTORTION;
402  v[7] = TAPER;
403  v[8] = STRETCH;
404  v[9] = DIAGONAL;
405  v[10] = SIZE;
406  break;
407  }
408 
409  case PRISM6:
410  case PRISM18:
411  {
412  // None yet
413  break;
414  }
415 
416  case PYRAMID5:
417  case PYRAMID13:
418  case PYRAMID14:
419  {
420  // None yet
421  break;
422  }
423 
424 
425 
426 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
427 
428  case INFEDGE2:
429  {
430  // None yet
431  break;
432  }
433 
434  case INFQUAD4:
435  case INFQUAD6:
436  {
437  // None yet
438  break;
439  }
440 
441  case INFHEX8:
442  case INFHEX16:
443  case INFHEX18:
444  {
445  // None yet
446  break;
447  }
448 
449  case INFPRISM6:
450  case INFPRISM12:
451  {
452  // None yet
453  break;
454  }
455 
456 #endif
457 
458 
459  default:
460  libmesh_error_msg("Undefined element type!");
461  }
462 
463  return v;
464 }

Variable Documentation

◆ num_quals

const unsigned int libMesh::Quality::num_quals = 16

The number of element quality types we have defined. This needs to be updated if you add one.

Definition at line 55 of file elem_quality.h.