MOAB: Mesh Oriented datABase
(version 5.4.1)
|
Computes the condition number of given element. More...
#include <ConditionNumberQualityMetric.hpp>
Public Member Functions | |
MESQUITE_EXPORT | ConditionNumberQualityMetric () |
virtual MESQUITE_EXPORT | ~ConditionNumberQualityMetric () |
virtual destructor ensures use of polymorphism during destruction | |
virtual MESQUITE_EXPORT std::string | get_name () const |
virtual MESQUITE_EXPORT int | get_negate_flag () const |
1 if metric should be minimized, -1 if metric should be maximized. | |
virtual MESQUITE_EXPORT bool | evaluate (PatchData &pd, size_t handle, double &value, MsqError &err) |
Get metric value at a logical location in the patch. |
Computes the condition number of given element.
The metric does not use the sample point functionality or the compute_weighted_jacobian. It evaluates the metric at the element vertices, and uses the isotropic ideal element. It does require a feasible region, and the metric needs to be minimized.
Definition at line 56 of file ConditionNumberQualityMetric.hpp.
Definition at line 42 of file ConditionNumberQualityMetric.cpp.
: AveragingQM( QualityMetric::LINEAR ) {}
virtual MESQUITE_EXPORT MBMesquite::ConditionNumberQualityMetric::~ConditionNumberQualityMetric | ( | ) | [inline, virtual] |
virtual destructor ensures use of polymorphism during destruction
Definition at line 62 of file ConditionNumberQualityMetric.hpp.
{}
bool ConditionNumberQualityMetric::evaluate | ( | PatchData & | pd, |
size_t | handle, | ||
double & | value, | ||
MsqError & | err | ||
) | [virtual] |
Get metric value at a logical location in the patch.
Evaluate the metric at one location in the PatchData.
pd | The patch. |
handle | The location in the patch (as passed back from get_evaluations). |
value | The output metric value. |
Implements MBMesquite::QualityMetric.
Definition at line 54 of file ConditionNumberQualityMetric.cpp.
References MBMesquite::AveragingQM::average_metrics(), MBMesquite::condition_number_2d(), MBMesquite::condition_number_3d(), MBMesquite::PatchData::element_by_index(), MBMesquite::MsqMeshEntity::get_element_type(), MBMesquite::PatchData::get_vertex_array(), MBMesquite::MsqMeshEntity::get_vertex_index_array(), MBMesquite::HEXAHEDRON, MSQ_ERRZERO, MBMesquite::MSQ_MAX_CAP, MBMesquite::MSQ_MAX_NUM_VERT_PER_ENT, MSQ_SETERR, MBMesquite::MSQ_SQRT_THREE, MBMesquite::MSQ_SQRT_THREE_INV, MBMesquite::MSQ_SQRT_TWO, MBMesquite::PRISM, MBMesquite::PYRAMID, MBMesquite::QUADRILATERAL, MBMesquite::TETRAHEDRON, MBMesquite::TRIANGLE, and MBMesquite::MsqError::UNSUPPORTED_ELEMENT.
{ const MsqMeshEntity* const element = &pd.element_by_index( handle ); bool return_flag; double met_vals[MSQ_MAX_NUM_VERT_PER_ENT]; fval = MSQ_MAX_CAP; const size_t* v_i = element->get_vertex_index_array(); // only 3 temp_vec will be sent to cond-num calculator, but the // additional vector3Ds may be needed during the calculations Vector3D temp_vec[6]; const MsqVertex* vertices = pd.get_vertex_array( err ); EntityTopology type = element->get_element_type(); switch( type ) { case TRIANGLE: temp_vec[0] = vertices[v_i[1]] - vertices[v_i[0]]; temp_vec[2] = vertices[v_i[2]] - vertices[v_i[0]]; // make relative to equilateral temp_vec[1] = ( ( 2 * temp_vec[2] ) - temp_vec[0] ) * MSQ_SQRT_THREE_INV; return_flag = condition_number_2d( temp_vec, handle, pd, fval, err ); MSQ_ERRZERO( err ); return return_flag; case QUADRILATERAL: temp_vec[0] = vertices[v_i[1]] - vertices[v_i[0]]; temp_vec[1] = vertices[v_i[3]] - vertices[v_i[0]]; return_flag = condition_number_2d( temp_vec, handle, pd, met_vals[0], err ); MSQ_ERRZERO( err ); if( !return_flag ) return return_flag; temp_vec[0] = vertices[v_i[2]] - vertices[v_i[1]]; temp_vec[1] = vertices[v_i[0]] - vertices[v_i[1]]; return_flag = condition_number_2d( temp_vec, handle, pd, met_vals[1], err ); MSQ_ERRZERO( err ); if( !return_flag ) return return_flag; temp_vec[0] = vertices[v_i[3]] - vertices[v_i[2]]; temp_vec[1] = vertices[v_i[1]] - vertices[v_i[2]]; return_flag = condition_number_2d( temp_vec, handle, pd, met_vals[2], err ); MSQ_ERRZERO( err ); if( !return_flag ) return return_flag; temp_vec[0] = vertices[v_i[0]] - vertices[v_i[3]]; temp_vec[1] = vertices[v_i[2]] - vertices[v_i[3]]; return_flag = condition_number_2d( temp_vec, handle, pd, met_vals[3], err ); MSQ_ERRZERO( err ); fval = average_metrics( met_vals, 4, err ); return return_flag; case TETRAHEDRON: temp_vec[0] = vertices[v_i[1]] - vertices[v_i[0]]; temp_vec[3] = vertices[v_i[2]] - vertices[v_i[0]]; temp_vec[4] = vertices[v_i[3]] - vertices[v_i[0]]; // transform to equilateral tet temp_vec[1] = ( ( 2 * temp_vec[3] ) - temp_vec[0] ) / MSQ_SQRT_THREE; temp_vec[2] = ( ( 3 * temp_vec[4] ) - temp_vec[0] - temp_vec[3] ) / ( MSQ_SQRT_THREE * MSQ_SQRT_TWO ); return_flag = condition_number_3d( temp_vec, pd, fval, err ); MSQ_ERRZERO( err ); return return_flag; case HEXAHEDRON: // transform to v_i[0] temp_vec[0] = vertices[v_i[1]] - vertices[v_i[0]]; temp_vec[1] = vertices[v_i[3]] - vertices[v_i[0]]; temp_vec[2] = vertices[v_i[4]] - vertices[v_i[0]]; return_flag = condition_number_3d( temp_vec, pd, met_vals[0], err ); MSQ_ERRZERO( err ); if( !return_flag ) return return_flag; temp_vec[0] = vertices[v_i[2]] - vertices[v_i[1]]; temp_vec[1] = vertices[v_i[0]] - vertices[v_i[1]]; temp_vec[2] = vertices[v_i[5]] - vertices[v_i[1]]; return_flag = condition_number_3d( temp_vec, pd, met_vals[1], err ); MSQ_ERRZERO( err ); if( !return_flag ) return return_flag; temp_vec[0] = vertices[v_i[3]] - vertices[v_i[2]]; temp_vec[1] = vertices[v_i[1]] - vertices[v_i[2]]; temp_vec[2] = vertices[v_i[6]] - vertices[v_i[2]]; return_flag = condition_number_3d( temp_vec, pd, met_vals[2], err ); MSQ_ERRZERO( err ); if( !return_flag ) return return_flag; temp_vec[0] = vertices[v_i[0]] - vertices[v_i[3]]; temp_vec[1] = vertices[v_i[2]] - vertices[v_i[3]]; temp_vec[2] = vertices[v_i[7]] - vertices[v_i[3]]; return_flag = condition_number_3d( temp_vec, pd, met_vals[3], err ); MSQ_ERRZERO( err ); if( !return_flag ) return return_flag; temp_vec[0] = vertices[v_i[7]] - vertices[v_i[4]]; temp_vec[1] = vertices[v_i[5]] - vertices[v_i[4]]; temp_vec[2] = vertices[v_i[0]] - vertices[v_i[4]]; return_flag = condition_number_3d( temp_vec, pd, met_vals[4], err ); MSQ_ERRZERO( err ); if( !return_flag ) return return_flag; temp_vec[0] = vertices[v_i[4]] - vertices[v_i[5]]; temp_vec[1] = vertices[v_i[6]] - vertices[v_i[5]]; temp_vec[2] = vertices[v_i[1]] - vertices[v_i[5]]; return_flag = condition_number_3d( temp_vec, pd, met_vals[5], err ); MSQ_ERRZERO( err ); if( !return_flag ) return return_flag; temp_vec[0] = vertices[v_i[5]] - vertices[v_i[6]]; temp_vec[1] = vertices[v_i[7]] - vertices[v_i[6]]; temp_vec[2] = vertices[v_i[2]] - vertices[v_i[6]]; return_flag = condition_number_3d( temp_vec, pd, met_vals[6], err ); MSQ_ERRZERO( err ); if( !return_flag ) return return_flag; temp_vec[0] = vertices[v_i[6]] - vertices[v_i[7]]; temp_vec[1] = vertices[v_i[4]] - vertices[v_i[7]]; temp_vec[2] = vertices[v_i[3]] - vertices[v_i[7]]; return_flag = condition_number_3d( temp_vec, pd, met_vals[7], err ); MSQ_ERRZERO( err ); fval = average_metrics( met_vals, 8, err ); MSQ_ERRZERO( err ); return return_flag; case PYRAMID: { unsigned num_adj; const unsigned* adj_idx; return_flag = true; for( size_t j = 0; j < 4; ++j ) // skip fifth vertex (apex) { adj_idx = TopologyInfo::adjacent_vertices( type, j, num_adj ); assert( num_adj == 3 ); temp_vec[0] = vertices[v_i[adj_idx[0]]] - vertices[v_i[j]]; temp_vec[1] = vertices[v_i[adj_idx[1]]] - vertices[v_i[j]]; // calculate last vect map to right tetrahedron temp_vec[3] = vertices[v_i[adj_idx[2]]] - vertices[v_i[adj_idx[0]]]; temp_vec[4] = vertices[v_i[adj_idx[2]]] - vertices[v_i[adj_idx[1]]]; temp_vec[2] = 0.5 * ( temp_vec[3] + temp_vec[4] ); return_flag = return_flag && condition_number_3d( temp_vec, pd, met_vals[j], err ); } fval = average_metrics( met_vals, 4, err ); MSQ_ERRZERO( err ); return return_flag; } case PRISM: { unsigned num_adj; const unsigned* adj_idx; return_flag = true; for( size_t j = 0; j < 6; ++j ) { adj_idx = TopologyInfo::adjacent_vertices( type, j, num_adj ); assert( num_adj == 3 ); temp_vec[0] = vertices[v_i[adj_idx[0]]] - vertices[v_i[j]]; temp_vec[1] = vertices[v_i[adj_idx[1]]] - vertices[v_i[j]]; temp_vec[2] = vertices[v_i[adj_idx[2]]] - vertices[v_i[j]]; // map to right tetrahedron temp_vec[1] += vertices[v_i[adj_idx[1]]]; temp_vec[1] -= vertices[v_i[adj_idx[0]]]; temp_vec[1] *= MSQ_SQRT_THREE_INV; return_flag = return_flag && condition_number_3d( temp_vec, pd, met_vals[j], err ); } fval = average_metrics( met_vals, 6, err ); MSQ_ERRZERO( err ); return return_flag; } default: MSQ_SETERR( err ) ( MsqError::UNSUPPORTED_ELEMENT, "Unsupported cell type (%d) for Condition Number quality metric.", type ); fval = MSQ_MAX_CAP; return false; } // end switch over element type return false; }
std::string ConditionNumberQualityMetric::get_name | ( | ) | const [virtual] |
Implements MBMesquite::QualityMetric.
Definition at line 44 of file ConditionNumberQualityMetric.cpp.
Referenced by QualityAssessorTest::test_print_stats().
{ return "Condition Number"; }
int ConditionNumberQualityMetric::get_negate_flag | ( | ) | const [virtual] |
1 if metric should be minimized, -1 if metric should be maximized.
Implements MBMesquite::QualityMetric.
Definition at line 49 of file ConditionNumberQualityMetric.cpp.
{
return 1;
}