libMesh::Quality Namespace Reference

Functions

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

Variables

const unsigned int num_quals = 16
 

Detailed Description

A namespace for quality utility functions.

Function Documentation

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 127 of file elem_quality.C.

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

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

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

Referenced by GETPOT_NAMESPACE::GetPot::_convert_to_type_no_default(), libMesh::EquationSystems::add_system(), libMesh::Factory< Base >::build(), libMesh::Utility::complex_filename(), libMesh::EquationSystems::delete_system(), libMesh::demangle(), DMLibMeshSetUpName_Private(), DMView_libMesh(), libMesh::Factory< Base >::Factory(), libMesh::Parameters::get(), libMesh::ReferenceCounter::get_info(), libMesh::ReferenceCounter::increment_constructor_count(), libMesh::ReferenceCounter::increment_destructor_count(), libMesh::XdrMGF::init(), libMesh::Parameters::insert(), libMesh::libmesh_cast_ptr(), libMesh::libmesh_cast_ref(), libMesh::Xdr::open(), GETPOT_NAMESPACE::GetPot::variable::operator=(), libMesh::SparseMatrix< Number >::print_matlab(), libMesh::NumericVector< Number >::print_matlab(), libMesh::process_trace(), libMesh::TetGenIO::read(), libMesh::CheckpointIO::read(), libMesh::UnstructuredMesh::read(), libMesh::EquationSystems::read(), libMesh::PltLoader::read_header(), libMesh::MeshData::read_tetgen(), libMesh::ReferenceCountedObject< RBParametrized >::ReferenceCountedObject(), libMesh::Parameters::set(), libMesh::Parameters::Parameter< T >::type(), libMesh::EnsightIO::write(), libMesh::CheckpointIO::write(), libMesh::UnstructuredMesh::write(), libMesh::TecplotIO::write_binary(), and libMesh::ReferenceCountedObject< RBParametrized >::~ReferenceCountedObject().

40 {
41  std::string its_name;
42 
43  switch (q)
44  {
45 
46  case ASPECT_RATIO:
47  its_name = "Aspect Ratio";
48  break;
49 
50  case SKEW:
51  its_name = "Skew";
52  break;
53 
54  case SHEAR:
55  its_name = "Shear";
56  break;
57 
58  case SHAPE:
59  its_name = "Shape";
60  break;
61 
62  case MAX_ANGLE:
63  its_name = "Maximum Angle";
64  break;
65 
66  case MIN_ANGLE:
67  its_name = "Minimum Angle";
68  break;
69 
70  case CONDITION:
71  its_name = "Condition Number";
72  break;
73 
74  case DISTORTION:
75  its_name = "Distortion";
76  break;
77 
78  case TAPER:
79  its_name = "Taper";
80  break;
81 
82  case WARP:
83  its_name = "Warp";
84  break;
85 
86  case STRETCH:
87  its_name = "Stretch";
88  break;
89 
90  case DIAGONAL:
91  its_name = "Diagonal";
92  break;
93 
94  case ASPECT_RATIO_BETA:
95  its_name = "AR Beta";
96  break;
97 
98  case ASPECT_RATIO_GAMMA:
99  its_name = "AR Gamma";
100  break;
101 
102  case SIZE:
103  its_name = "Size";
104  break;
105 
106  case JACOBIAN:
107  its_name = "Jacobian";
108  break;
109 
110  default:
111  its_name = "Unknown";
112  break;
113  }
114 
115  return its_name;
116 }
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 321 of file elem_quality.C.

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

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

Variable Documentation

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 45 of file elem_quality.h.


Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:33 UTC

Hosted By:
SourceForge.net Logo