Poster of Linux kernelThe best gift for a Linux geek
Quality

Quality

Section: C Library Functions (3) Updated: Thu Apr 7 2011
Local index Up
 

NAME

Quality -  

SYNOPSIS


 

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 Quality::describe (const ElemQualityq)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 126 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.

{
  
  std::ostringstream desc; 
  
  switch (q)
    {
      
    case ASPECT_RATIO:
      desc << 'Max edge length ratio
           << 'at element center.
           << '
           << 'Suggested ranges:
           << 'Hexes: (1 -> 4)
           << 'Quads: (1 -> 4)';
      break;
      
    case SKEW:
      desc << 'Maximum |cos A|, where A
           << 'is the angle between edges
           << 'at element center.
           << '
           << 'Suggested ranges:
           << 'Hexes: (0 -> 0.5)
           << 'Quads: (0 -> 0.5)';
      break;
      
    case SHEAR:
      desc << 'LIBMESH_DIM / K(Js)
           << '
           << 'LIBMESH_DIM   = element dimension.
           << 'K(Js) = Condition number of 
           << '        Jacobian skew matrix.
           << '
           << 'Suggested ranges:
           << 'Hexes(LIBMESH_DIM=3): (0.3 -> 1)
           << 'Quads(LIBMESH_DIM=2): (0.3 -> 1)';
      break;
      
    case SHAPE:
      desc << 'LIBMESH_DIM / K(Jw)
           << '
           << 'LIBMESH_DIM   = element dimension.
           << 'K(Jw) = Condition number of 
           << '        weighted Jacobian
           << '        matrix.
           << '
           << 'Suggested ranges:
           << 'Hexes(LIBMESH_DIM=3): (0.3 -> 1)
           << 'Tets(LIBMESH_DIM=3): (0.2 -> 1)
           << 'Quads(LIBMESH_DIM=2): (0.3 -> 1).';
      break;
      
    case MAX_ANGLE:
      desc << 'Largest included angle.
           << '
           << 'Suggested ranges:
           << 'Quads: (90 -> 135)
           << 'Triangles: (60 -> 90)';
      break;

    case MIN_ANGLE:
      desc << 'Smallest included angle.
           << '
           << 'Suggested ranges:
           << 'Quads: (45 -> 90)
           << 'Triangles: (30 -> 60)';
      break;
      
    case CONDITION:
      desc << 'Condition number of the
           << 'Jacobian matrix.
           << '
           << 'Suggested ranges:
           << 'Quads: (1 -> 4) 
           << 'Hexes: (1 -> 8)
           << 'Tris: (1 -> 1.3)
           << 'Tets: (1 -> 3)';
      break;
      
    case DISTORTION:
      desc << 'min |J| * A / <A>
           << '
           << '|J| = norm of Jacobian matrix
           << ' A  = actual area
           << '<A> = reference area
           << '
           << 'Suggested ranges:
           << 'Quads: (0.6 -> 1), <A>=4 
           << 'Hexes: (0.6 -> 1), <A>=8
           << 'Tris: (0.6 -> 1), <A>=1/2
           << 'Tets: (0.6 -> 1), <A>=1/6';
      break;
      
    case TAPER:
      desc << 'Maximum ratio of lengths
           << 'derived from opposited edges.
           << '
           << 'Suggested ranges:
           << 'Quads: (0.7 -> 1)
           << 'Hexes: (0.4 -> 1)';
      break;
      
    case WARP:
      desc << 'cos D
           << '
           << 'D = minimum dihedral angle
           << '    formed by diagonals.
           << '
           << 'Suggested ranges:
           << 'Quads: (0.9 -> 1)';
      break;
      
    case STRETCH:
      desc << 'Sqrt(3) * L_min / L_max
           << '
           << 'L_min = minimum edge length.
           << 'L_max = maximum edge length.
           << '
           << 'Suggested ranges:
           << 'Quads: (0.25 -> 1)
           << 'Hexes: (0.25 -> 1)';
      break;
      
    case DIAGONAL:
      desc << 'D_min / D_max
           << '
           << 'D_min = minimum diagonal.
           << 'D_max = maximum diagonal.
           << '
           << 'Suggested ranges:
           << 'Hexes: (0.65 -> 1)';
      break;

    case ASPECT_RATIO_BETA:
      desc << 'CR / (3 * IR)
           << '
           << 'CR = circumsphere radius
           << 'IR = inscribed sphere radius
           << '
           << 'Suggested ranges:
           << 'Tets: (1 -> 3)';
      break;
      
    case ASPECT_RATIO_GAMMA:
      desc << 'S^(3/2) / 8.479670 * V
           << '
           << 'S = sum(si*si/6)
           << 'si = edge length
           << 'V = volume
           << '
           << 'Suggested ranges:
           << 'Tets: (1 -> 3)';
      break;

    case SIZE:
      desc << 'min (|J|, |1/J|)
           << '
           << '|J| = norm of Jacobian matrix.
           << '
           << 'Suggested ranges:
           << 'Quads: (0.3 -> 1) 
           << 'Hexes: (0.5 -> 1)
           << 'Tris: (0.25 -> 1)
           << 'Tets: (0.2 -> 1)';
      break;
      
    case JACOBIAN:
      desc << 'Minimum Jacobian divided by
           << 'the lengths of the LIBMESH_DIM
           << 'largest edge vectors.
           << '
           << 'LIBMESH_DIM = element dimension.
           << '
           << 'Suggested ranges:
           << 'Quads: (0.5 -> 1) 
           << 'Hexes: (0.5 -> 1)
           << 'Tris: (0.5 -> 1.155)
           << 'Tets: (0.5 -> 1.414)';
      break;

    default:
      desc << 'Unknown';
      break;
    }
  
  return desc.str();
}
 

std::string Quality::name (const ElemQualityq)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 38 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::__find_variable(), abi_demangle(), EquationSystems::add_system(), Factory< Base >::build(), Utility::complex_filename(), EquationSystems::delete_system(), Factory< Base >::Factory(), Parameters::get(), ReferenceCounter::get_info(), ReferenceCounter::increment_constructor_count(), ReferenceCounter::increment_destructor_count(), XdrMGF::init(), Xdr::open(), GetPot::variable::operator=(), SparseMatrix< Number >::print_matlab(), NumericVector< Number >::print_matlab(), UnstructuredMesh::read(), TetGenIO::read(), PltLoader::read_header(), MeshData::read_tetgen(), ReferenceCountedObject< Value >::ReferenceCountedObject(), Parameters::set(), VTKIO::solution_to_vtk(), Parameters::Parameter< T >::type(), UnstructuredMesh::write(), EnsightIO::write(), and ReferenceCountedObject< Value >::~ReferenceCountedObject().

{
  std::string its_name;
  
  switch (q)
    {
      
    case ASPECT_RATIO:
      its_name = 'Aspect Ratio';
      break;
      
    case SKEW:
      its_name = 'Skew';
      break;
      
    case SHEAR:
      its_name = 'Shear';
      break;
      
    case SHAPE:
      its_name = 'Shape';
      break;
      
    case MAX_ANGLE:
      its_name = 'Maximum Angle';
      break;

    case MIN_ANGLE:
      its_name = 'Minimum Angle';
      break;
      
    case CONDITION:
      its_name = 'Condition Number';
      break;
      
    case DISTORTION:
      its_name = 'Distortion';
      break;
      
    case TAPER:
      its_name = 'Taper';
      break;
      
    case WARP:
      its_name = 'Warp';
      break;
      
    case STRETCH:
      its_name = 'Stretch';
      break;
      
    case DIAGONAL:
      its_name = 'Diagonal';
      break;

    case ASPECT_RATIO_BETA:
      its_name = 'AR Beta';
      break;
      
    case ASPECT_RATIO_GAMMA:
      its_name = 'AR Gamma';
      break;

    case SIZE:
      its_name = 'Size';
      break;
      
    case JACOBIAN:
      its_name = 'Jacobian';
      break;

    default:
      its_name = 'Unknown';
      break;
    }
  
  return its_name;
}
 

std::vector< ElemQuality > Quality::valid (const ElemTypet)Returns:

the valid ElemQuality metrics for a given ElemType element type.

Returns all valid quality metrics for element type t.

Definition at line 320 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, libMeshEnums::PRISM18, libMeshEnums::PRISM6, 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.

{
  std::vector<ElemQuality> v;
  
  switch (t)
    {
    case EDGE2:
    case EDGE3:
    case EDGE4:
      {
        // None yet
        break;
      }

    case TRI3:
    case TRI6:
      {
        v.resize(7);
        v[0] = MAX_ANGLE;
        v[1] = MIN_ANGLE;
        v[2] = CONDITION;
        v[3] = JACOBIAN;
        v[4] = SIZE;
        v[5] = SHAPE;
        v[6] = DISTORTION;
        break;
      }

    case QUAD4:
    case QUAD8:
    case QUAD9:
      {
        v.resize(13);
        v[0]  = ASPECT_RATIO;
        v[1]  = SKEW;
        v[2]  = TAPER;
        v[3]  = WARP;
        v[4]  = STRETCH;
        v[5]  = MIN_ANGLE;
        v[6]  = MAX_ANGLE;
        v[7]  = CONDITION;
        v[8]  = JACOBIAN;
        v[9]  = SHEAR;
        v[10] = SHAPE;
        v[11] = SIZE;
        v[12] = DISTORTION;
        break;
      }

    case TET4:
    case TET10:
      {
        v.resize(7);
        v[0]  = ASPECT_RATIO_BETA;
        v[1]  = ASPECT_RATIO_GAMMA;
        v[2]  = CONDITION;
        v[3]  = JACOBIAN;
        v[4]  = SHAPE;
        v[5]  = SIZE;
        v[6]  = DISTORTION;
        break;
      }

    case HEX8:
    case HEX20:
    case HEX27:
      {
        v.resize(11);
        v[0]  = ASPECT_RATIO;
        v[1]  = SKEW;
        v[2]  = SHEAR;
        v[3] = SHAPE;
        v[4]  = CONDITION;
        v[5]  = JACOBIAN;
        v[6]  = DISTORTION;
        v[7]  = TAPER;
        v[8]  = STRETCH;
        v[9]  = DIAGONAL;
        v[10]  = SIZE;
        break;
      }

    case PRISM6:
    case PRISM18:
      {
        // None yet
        break;
      }

    case PYRAMID5:
      {
        // None yet
        break;
      }



#ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS

    case INFEDGE2:
      {
        // None yet
        break;
      }

    case INFQUAD4:
    case INFQUAD6:
      {
        // None yet
        break;
      }

    case INFHEX8:
    case INFHEX16:
    case INFHEX18:
      {
        // None yet
        break;
      }

    case INFPRISM6:
    case INFPRISM12:
      {
        // None yet
        break;
      }

#endif

      
    default:
      {
        std::cout << 'Undefined element type!.' << std::endl;
        libmesh_error();
      }
    }
  
  return v;
}
 

Variable Documentation

 

const unsigned int Quality::num_quals = 16The number of element quality types we have defined. This needs to be updated if you add one.

Definition at line 44 of file elem_quality.h.  

Author

Generated automatically by Doxygen for libMesh from the source code.


 

Index

NAME
SYNOPSIS
Functions
Variables
Detailed Description
Function Documentation
std::string Quality::describe (const ElemQualityq)Returns:
std::string Quality::name (const ElemQualityq)Returns:
std::vector< ElemQuality > Quality::valid (const ElemTypet)Returns:
Variable Documentation
const unsigned int Quality::num_quals = 16The number of element quality types we have defined. This needs to be updated if you add one.
Author

This document was created by man2html, using the manual pages.
Time: 21:53:48 GMT, April 16, 2011