![]() |
Mesh Oriented datABase
(version 5.4.1)
Array-based unstructured mesh datastructure
|
#include "moab/verdict.h"
#include "VerdictVector.hpp"
#include "verdict_defines.hpp"
#include "V_GaussIntegration.hpp"
#include <memory.h>
Go to the source code of this file.
Defines | |
#define | VERDICT_EXPORTS |
Functions | |
int | get_weight (double &m11, double &m21, double &m12, double &m22) |
C_FUNC_DEF void | v_set_quad_size (double size) |
return the average area of a quad | |
VerdictBoolean | is_collapsed_quad (double coordinates[][3]) |
returns whether the quad is collapsed or not | |
void | make_quad_edges (VerdictVector edges[4], double coordinates[][3]) |
void | signed_corner_areas (double areas[4], double coordinates[][3]) |
void | localize_quad_coordinates (VerdictVector nodes[4]) |
void | localize_quad_for_ef (VerdictVector node_pos[4]) |
VerdictVector | quad_normal (double coordinates[][3]) |
C_FUNC_DEF double | v_quad_edge_ratio (int, double coordinates[][3]) |
Calculates quad edge ratio. | |
C_FUNC_DEF double | v_quad_max_edge_ratio (int, double coordinates[][3]) |
Calculates quad maximum of edge ratio. | |
C_FUNC_DEF double | v_quad_aspect_ratio (int, double coordinates[][3]) |
Calculates quad aspect ratio. | |
C_FUNC_DEF double | v_quad_radius_ratio (int, double coordinates[][3]) |
Calculates quad radius ratio. | |
C_FUNC_DEF double | v_quad_med_aspect_frobenius (int, double coordinates[][3]) |
Calculates quad average Frobenius aspect. | |
C_FUNC_DEF double | v_quad_max_aspect_frobenius (int, double coordinates[][3]) |
Calculates quad maximum Frobenius aspect. | |
C_FUNC_DEF double | v_quad_skew (int, double coordinates[][3]) |
Calculates quad skew metric. | |
C_FUNC_DEF double | v_quad_taper (int, double coordinates[][3]) |
Calculates quad taper metric. | |
C_FUNC_DEF double | v_quad_warpage (int, double coordinates[][3]) |
Calculates quad warpage metric. | |
C_FUNC_DEF double | v_quad_area (int, double coordinates[][3]) |
Calculates quad area. | |
C_FUNC_DEF double | v_quad_stretch (int, double coordinates[][3]) |
Calculates quad strech metric. | |
C_FUNC_DEF double | v_quad_maximum_angle (int, double coordinates[][3]) |
Calculates quad's largest angle. | |
C_FUNC_DEF double | v_quad_minimum_angle (int, double coordinates[][3]) |
Calculates quad's smallest angle. | |
C_FUNC_DEF double | v_quad_oddy (int, double coordinates[][3]) |
Calculates quad oddy metric. | |
C_FUNC_DEF double | v_quad_condition (int, double coordinates[][3]) |
Calculates quad condition number metric. | |
C_FUNC_DEF double | v_quad_jacobian (int, double coordinates[][3]) |
Calculates quad jacobian. | |
C_FUNC_DEF double | v_quad_scaled_jacobian (int, double coordinates[][3]) |
Calculates quad scaled jacobian. | |
C_FUNC_DEF double | v_quad_shear (int, double coordinates[][3]) |
Calculates quad shear metric. | |
C_FUNC_DEF double | v_quad_shape (int, double coordinates[][3]) |
Calculates quad shape metric. | |
C_FUNC_DEF double | v_quad_relative_size_squared (int, double coordinates[][3]) |
Calculates quad relative size metric. | |
C_FUNC_DEF double | v_quad_shape_and_size (int num_nodes, double coordinates[][3]) |
Calculates quad shape-size metric. | |
C_FUNC_DEF double | v_quad_shear_and_size (int num_nodes, double coordinates[][3]) |
Calculates quad shear-size metric. | |
C_FUNC_DEF double | v_quad_distortion (int num_nodes, double coordinates[][3]) |
Calculates quad distortion metric. | |
C_FUNC_DEF void | v_quad_quality (int num_nodes, double coordinates[][3], unsigned int metrics_request_flag, QuadMetricVals *metric_vals) |
Calculates quality metrics for quadrilateral elements. | |
Variables | |
double | verdict_quad_size = 0 |
the average area of a quad |
#define VERDICT_EXPORTS |
Definition at line 23 of file V_QuadMetric.cpp.
int get_weight | ( | double & | m11, |
double & | m21, | ||
double & | m12, | ||
double & | m22 | ||
) |
weights based on the average size of a quad
Definition at line 37 of file V_QuadMetric.cpp.
References verdict_quad_size.
Referenced by v_quad_quality(), v_quad_relative_size_squared(), v_tet_quality(), and v_tet_relative_size_squared().
{
m11 = 1;
m21 = 0;
m12 = 0;
m22 = 1;
double scale = sqrt( verdict_quad_size / ( m11 * m22 - m21 * m12 ) );
m11 *= scale;
m21 *= scale;
m12 *= scale;
m22 *= scale;
return 1;
}
VerdictBoolean is_collapsed_quad | ( | double | coordinates[][3] | ) |
returns whether the quad is collapsed or not
Definition at line 62 of file V_QuadMetric.cpp.
References VERDICT_FALSE, and VERDICT_TRUE.
Referenced by v_quad_condition(), v_quad_distortion(), v_quad_jacobian(), v_quad_maximum_angle(), v_quad_minimum_angle(), v_quad_quality(), and v_quad_scaled_jacobian().
{
if( coordinates[3][0] == coordinates[2][0] && coordinates[3][1] == coordinates[2][1] &&
coordinates[3][2] == coordinates[2][2] )
return VERDICT_TRUE;
else
return VERDICT_FALSE;
}
void localize_quad_coordinates | ( | VerdictVector | nodes[4] | ) |
localize the coordinates of a quad
localizing puts the centriod of the quad at the orgin and also rotates the quad such that edge (0,1) is aligned with the x axis and the quad normal lines up with the y axis.
Definition at line 121 of file V_QuadMetric.cpp.
References center(), VerdictVector::normalize(), VerdictVector::set(), VerdictVector::x(), VerdictVector::y(), and VerdictVector::z().
{
int i;
VerdictVector global[4] = { nodes[0], nodes[1], nodes[2], nodes[3] };
VerdictVector center = ( global[0] + global[1] + global[2] + global[3] ) / 4.0;
for( i = 0; i < 4; i++ )
global[i] -= center;
VerdictVector vector_diff;
VerdictVector vector_sum;
VerdictVector ref_point( 0.0, 0.0, 0.0 );
VerdictVector tmp_vector, normal( 0.0, 0.0, 0.0 );
VerdictVector vector1, vector2;
for( i = 0; i < 4; i++ )
{
vector1 = global[i];
vector2 = global[( i + 1 ) % 4];
vector_diff = vector2 - vector1;
ref_point += vector1;
vector_sum = vector1 + vector2;
tmp_vector.set( vector_diff.y() * vector_sum.z(), vector_diff.z() * vector_sum.x(),
vector_diff.x() * vector_sum.y() );
normal += tmp_vector;
}
normal.normalize();
normal *= -1.0;
VerdictVector local_x_axis = global[1] - global[0];
local_x_axis.normalize();
VerdictVector local_y_axis = normal * local_x_axis;
local_y_axis.normalize();
for( i = 0; i < 4; i++ )
{
nodes[i].x( global[i] % local_x_axis );
nodes[i].y( global[i] % local_y_axis );
nodes[i].z( global[i] % normal );
}
}
void localize_quad_for_ef | ( | VerdictVector | node_pos[4] | ) |
moves and rotates the quad such that it enables us to use components of ef's
Definition at line 170 of file V_QuadMetric.cpp.
References VerdictVector::normalize(), VerdictVector::x(), and VerdictVector::y().
{
VerdictVector centroid( node_pos[0] );
centroid += node_pos[1];
centroid += node_pos[2];
centroid += node_pos[3];
centroid /= 4.0;
node_pos[0] -= centroid;
node_pos[1] -= centroid;
node_pos[2] -= centroid;
node_pos[3] -= centroid;
VerdictVector rotate = node_pos[1] + node_pos[2] - node_pos[3] - node_pos[0];
rotate.normalize();
double cosine = rotate.x();
double sine = rotate.y();
double xnew;
for( int i = 0; i < 4; i++ )
{
xnew = cosine * node_pos[i].x() + sine * node_pos[i].y();
node_pos[i].y( -sine * node_pos[i].x() + cosine * node_pos[i].y() );
node_pos[i].x( xnew );
}
}
void make_quad_edges | ( | VerdictVector | edges[4], |
double | coordinates[][3] | ||
) |
Definition at line 72 of file V_QuadMetric.cpp.
References VerdictVector::set().
Referenced by signed_corner_areas(), v_quad_aspect_ratio(), v_quad_edge_ratio(), v_quad_max_aspect_frobenius(), v_quad_med_aspect_frobenius(), v_quad_quality(), v_quad_radius_ratio(), v_quad_scaled_jacobian(), v_quad_shape(), v_quad_stretch(), and v_quad_warpage().
{
edges[0].set( coordinates[1][0] - coordinates[0][0], coordinates[1][1] - coordinates[0][1],
coordinates[1][2] - coordinates[0][2] );
edges[1].set( coordinates[2][0] - coordinates[1][0], coordinates[2][1] - coordinates[1][1],
coordinates[2][2] - coordinates[1][2] );
edges[2].set( coordinates[3][0] - coordinates[2][0], coordinates[3][1] - coordinates[2][1],
coordinates[3][2] - coordinates[2][2] );
edges[3].set( coordinates[0][0] - coordinates[3][0], coordinates[0][1] - coordinates[3][1],
coordinates[0][2] - coordinates[3][2] );
}
VerdictVector quad_normal | ( | double | coordinates[][3] | ) |
returns the normal vector of a quad
Definition at line 204 of file V_QuadMetric.cpp.
References VerdictVector::normalize(), and VerdictVector::set().
Referenced by v_quad_distortion().
{
// get normal at node 0
VerdictVector edge0, edge1;
edge0.set( coordinates[1][0] - coordinates[0][0], coordinates[1][1] - coordinates[0][1],
coordinates[1][2] - coordinates[0][2] );
edge1.set( coordinates[3][0] - coordinates[0][0], coordinates[3][1] - coordinates[0][1],
coordinates[3][2] - coordinates[0][2] );
VerdictVector norm0 = edge0 * edge1;
norm0.normalize();
// because some faces may have obtuse angles, check another normal at
// node 2 for consistent sense
edge0.set( coordinates[2][0] - coordinates[3][0], coordinates[2][1] - coordinates[3][1],
coordinates[2][2] - coordinates[3][2] );
edge1.set( coordinates[2][0] - coordinates[1][0], coordinates[2][1] - coordinates[1][1],
coordinates[2][2] - coordinates[1][2] );
VerdictVector norm2 = edge0 * edge1;
norm2.normalize();
// if these two agree, we are done, else test a third to decide
if( ( norm0 % norm2 ) > 0.0 )
{
norm0 += norm2;
norm0 *= 0.5;
return norm0;
}
// test normal at node1
edge0.set( coordinates[1][0] - coordinates[2][0], coordinates[1][1] - coordinates[2][1],
coordinates[1][2] - coordinates[2][2] );
edge1.set( coordinates[1][0] - coordinates[0][0], coordinates[1][1] - coordinates[0][1],
coordinates[1][2] - coordinates[0][2] );
VerdictVector norm1 = edge0 * edge1;
norm1.normalize();
if( ( norm0 % norm1 ) > 0.0 )
{
norm0 += norm1;
norm0 *= 0.5;
return norm0;
}
else
{
norm2 += norm1;
norm2 *= 0.5;
return norm2;
}
}
void signed_corner_areas | ( | double | areas[4], |
double | coordinates[][3] | ||
) |
Definition at line 85 of file V_QuadMetric.cpp.
References make_quad_edges(), and VerdictVector::normalize().
Referenced by v_quad_area(), v_quad_condition(), v_quad_jacobian(), v_quad_maximum_angle(), v_quad_quality(), v_quad_scaled_jacobian(), and v_quad_shape().
{
VerdictVector edges[4];
make_quad_edges( edges, coordinates );
VerdictVector corner_normals[4];
corner_normals[0] = edges[3] * edges[0];
corner_normals[1] = edges[0] * edges[1];
corner_normals[2] = edges[1] * edges[2];
corner_normals[3] = edges[2] * edges[3];
// principal axes
VerdictVector principal_axes[2];
principal_axes[0] = edges[0] - edges[2];
principal_axes[1] = edges[1] - edges[3];
// quad center unit normal
VerdictVector unit_center_normal;
unit_center_normal = principal_axes[0] * principal_axes[1];
unit_center_normal.normalize();
areas[0] = unit_center_normal % corner_normals[0];
areas[1] = unit_center_normal % corner_normals[1];
areas[2] = unit_center_normal % corner_normals[2];
areas[3] = unit_center_normal % corner_normals[3];
}
C_FUNC_DEF double v_quad_area | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad area.
the area of a quad
jacobian at quad center
Definition at line 611 of file V_QuadMetric.cpp.
References signed_corner_areas(), VERDICT_DBL_MAX, VERDICT_MAX, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure(), v_quad_quality(), and v_quad_relative_size_squared().
{
double corner_areas[4];
signed_corner_areas( corner_areas, coordinates );
double area = 0.25 * ( corner_areas[0] + corner_areas[1] + corner_areas[2] + corner_areas[3] );
if( area > 0 ) return (double)VERDICT_MIN( area, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( area, -VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_aspect_ratio | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad aspect ratio.
the aspect ratio of a quad
NB (P. Pebay 01/20/07): this is a generalization of the triangle aspect ratio using Heron's formula.
Definition at line 351 of file V_QuadMetric.cpp.
References VerdictVector::length(), make_quad_edges(), mb, VERDICT_DBL_MAX, VERDICT_DBL_MIN, VERDICT_MAX, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure(), and v_quad_quality().
{
VerdictVector edges[4];
make_quad_edges( edges, coordinates );
double a1 = edges[0].length();
double b1 = edges[1].length();
double c1 = edges[2].length();
double d1 = edges[3].length();
double ma = a1 > b1 ? a1 : b1;
double mb = c1 > d1 ? c1 : d1;
double hm = ma > mb ? ma : mb;
VerdictVector ab = edges[0] * edges[1];
VerdictVector cd = edges[2] * edges[3];
double denominator = ab.length() + cd.length();
if( denominator < VERDICT_DBL_MIN ) return (double)VERDICT_DBL_MAX;
double aspect_ratio = .5 * hm * ( a1 + b1 + c1 + d1 ) / denominator;
if( aspect_ratio > 0 ) return (double)VERDICT_MIN( aspect_ratio, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( aspect_ratio, -VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_condition | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad condition number metric.
the condition of a quad
maximum condition number of the Jacobian matrix at 4 corners
Definition at line 828 of file V_QuadMetric.cpp.
References is_collapsed_quad(), VerdictVector::set(), signed_corner_areas(), v_tri_condition(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, VERDICT_MAX, VERDICT_MIN, and VERDICT_TRUE.
Referenced by moab::VerdictWrapper::quality_measure().
{
if( is_collapsed_quad( coordinates ) == VERDICT_TRUE ) return v_tri_condition( 3, coordinates );
double areas[4];
signed_corner_areas( areas, coordinates );
double max_condition = 0.;
VerdictVector xxi, xet;
double condition;
for( int i = 0; i < 4; i++ )
{
xxi.set( coordinates[i][0] - coordinates[( i + 1 ) % 4][0], coordinates[i][1] - coordinates[( i + 1 ) % 4][1],
coordinates[i][2] - coordinates[( i + 1 ) % 4][2] );
xet.set( coordinates[i][0] - coordinates[( i + 3 ) % 4][0], coordinates[i][1] - coordinates[( i + 3 ) % 4][1],
coordinates[i][2] - coordinates[( i + 3 ) % 4][2] );
if( areas[i] < VERDICT_DBL_MIN )
condition = VERDICT_DBL_MAX;
else
condition = ( xxi % xxi + xet % xet ) / areas[i];
max_condition = VERDICT_MAX( max_condition, condition );
}
max_condition /= 2;
if( max_condition > 0 ) return (double)VERDICT_MIN( max_condition, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( max_condition, -VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_distortion | ( | int | num_nodes, |
double | coordinates[][3] | ||
) |
Calculates quad distortion metric.
the distortion of a quad
Definition at line 1052 of file V_QuadMetric.cpp.
References GaussIntegration::calculate_derivative_at_nodes(), GaussIntegration::calculate_shape_function_2d_quad(), dot_product(), moab::GeomUtil::first(), GaussIntegration::get_shape_func(), GaussIntegration::initialize(), is_collapsed_quad(), length(), VerdictVector::length(), maxNumberNodes, maxTotalNumberGaussPoints, VerdictVector::normalize(), quad_normal(), VerdictVector::set(), VERDICT_DBL_MAX, VERDICT_MIN, and VERDICT_TRUE.
Referenced by moab::VerdictWrapper::quality_measure(), and v_quad_quality().
{
// To calculate distortion for linear and 2nd order quads
// distortion = {min(|J|)/actual area}*{parent area}
// parent area = 4 for a quad.
// min |J| is the minimum over nodes and gaussian integration points
// created by Ling Pan, CAT on 4/30/01
double element_area = 0.0, distrt, thickness_gauss;
double cur_jacobian = 0., sign_jacobian, jacobian;
VerdictVector aa, bb, cc, normal_at_point, xin;
// use 2x2 gauss points for linear quads and 3x3 for 2nd order quads
int number_of_gauss_points = 0;
if( num_nodes == 4 )
{ // 2x2 quadrature rule
number_of_gauss_points = 2;
}
else if( num_nodes == 8 )
{ // 3x3 quadrature rule
number_of_gauss_points = 3;
}
int total_number_of_gauss_points = number_of_gauss_points * number_of_gauss_points;
VerdictVector face_normal = quad_normal( coordinates );
double distortion = VERDICT_DBL_MAX;
VerdictVector first, second;
int i;
// Will work out the case for collapsed quad later
if( is_collapsed_quad( coordinates ) == VERDICT_TRUE )
{
for( i = 0; i < 3; i++ )
{
first.set( coordinates[i][0] - coordinates[( i + 1 ) % 3][0],
coordinates[i][1] - coordinates[( i + 1 ) % 3][1],
coordinates[i][2] - coordinates[( i + 1 ) % 3][2] );
second.set( coordinates[i][0] - coordinates[( i + 2 ) % 3][0],
coordinates[i][1] - coordinates[( i + 2 ) % 3][1],
coordinates[i][2] - coordinates[( i + 2 ) % 3][2] );
sign_jacobian = ( face_normal % ( first * second ) ) > 0 ? 1. : -1.;
cur_jacobian = sign_jacobian * ( first * second ).length();
distortion = VERDICT_MIN( distortion, cur_jacobian );
}
element_area = ( first * second ).length() / 2.0;
distortion /= element_area;
}
else
{
double shape_function[maxTotalNumberGaussPoints][maxNumberNodes];
double dndy1[maxTotalNumberGaussPoints][maxNumberNodes];
double dndy2[maxTotalNumberGaussPoints][maxNumberNodes];
double weight[maxTotalNumberGaussPoints];
// create an object of GaussIntegration
GaussIntegration::initialize( number_of_gauss_points, num_nodes );
GaussIntegration::calculate_shape_function_2d_quad();
GaussIntegration::get_shape_func( shape_function[0], dndy1[0], dndy2[0], weight );
// calculate element area
int ife, ja;
for( ife = 0; ife < total_number_of_gauss_points; ife++ )
{
aa.set( 0.0, 0.0, 0.0 );
bb.set( 0.0, 0.0, 0.0 );
for( ja = 0; ja < num_nodes; ja++ )
{
xin.set( coordinates[ja][0], coordinates[ja][1], coordinates[ja][2] );
aa += dndy1[ife][ja] * xin;
bb += dndy2[ife][ja] * xin;
}
normal_at_point = aa * bb;
jacobian = normal_at_point.length();
element_area += weight[ife] * jacobian;
}
double dndy1_at_node[maxNumberNodes][maxNumberNodes];
double dndy2_at_node[maxNumberNodes][maxNumberNodes];
GaussIntegration::calculate_derivative_at_nodes( dndy1_at_node, dndy2_at_node );
VerdictVector normal_at_nodes[9];
// evaluate normal at nodes and distortion values at nodes
int jai;
for( ja = 0; ja < num_nodes; ja++ )
{
aa.set( 0.0, 0.0, 0.0 );
bb.set( 0.0, 0.0, 0.0 );
for( jai = 0; jai < num_nodes; jai++ )
{
xin.set( coordinates[jai][0], coordinates[jai][1], coordinates[jai][2] );
aa += dndy1_at_node[ja][jai] * xin;
bb += dndy2_at_node[ja][jai] * xin;
}
normal_at_nodes[ja] = aa * bb;
normal_at_nodes[ja].normalize();
}
// determine if element is flat
bool flat_element = true;
double dot_product;
for( ja = 0; ja < num_nodes; ja++ )
{
dot_product = normal_at_nodes[0] % normal_at_nodes[ja];
if( fabs( dot_product ) < 0.99 )
{
flat_element = false;
break;
}
}
// take into consideration of the thickness of the element
double thickness;
// get_quad_thickness(face, element_area, thickness );
thickness = 0.001 * sqrt( element_area );
// set thickness gauss point location
double zl = 0.5773502691896;
if( flat_element ) zl = 0.0;
int no_gauss_pts_z = ( flat_element ) ? 1 : 2;
double thickness_z;
int igz;
// loop on Gauss points
for( ife = 0; ife < total_number_of_gauss_points; ife++ )
{
// loop on the thickness direction gauss points
for( igz = 0; igz < no_gauss_pts_z; igz++ )
{
zl = -zl;
thickness_z = zl * thickness / 2.0;
aa.set( 0.0, 0.0, 0.0 );
bb.set( 0.0, 0.0, 0.0 );
cc.set( 0.0, 0.0, 0.0 );
for( ja = 0; ja < num_nodes; ja++ )
{
xin.set( coordinates[ja][0], coordinates[ja][1], coordinates[ja][2] );
xin += thickness_z * normal_at_nodes[ja];
aa += dndy1[ife][ja] * xin;
bb += dndy2[ife][ja] * xin;
thickness_gauss = shape_function[ife][ja] * thickness / 2.0;
cc += thickness_gauss * normal_at_nodes[ja];
}
normal_at_point = aa * bb;
// jacobian = normal_at_point.length();
distrt = cc % normal_at_point;
if( distrt < distortion ) distortion = distrt;
}
}
// loop through nodal points
for( ja = 0; ja < num_nodes; ja++ )
{
for( igz = 0; igz < no_gauss_pts_z; igz++ )
{
zl = -zl;
thickness_z = zl * thickness / 2.0;
aa.set( 0.0, 0.0, 0.0 );
bb.set( 0.0, 0.0, 0.0 );
cc.set( 0.0, 0.0, 0.0 );
for( jai = 0; jai < num_nodes; jai++ )
{
xin.set( coordinates[jai][0], coordinates[jai][1], coordinates[jai][2] );
xin += thickness_z * normal_at_nodes[ja];
aa += dndy1_at_node[ja][jai] * xin;
bb += dndy2_at_node[ja][jai] * xin;
if( jai == ja )
thickness_gauss = thickness / 2.0;
else
thickness_gauss = 0.;
cc += thickness_gauss * normal_at_nodes[jai];
}
}
normal_at_point = aa * bb;
sign_jacobian = ( face_normal % normal_at_point ) > 0 ? 1. : -1.;
distrt = sign_jacobian * ( cc % normal_at_point );
if( distrt < distortion ) distortion = distrt;
}
if( element_area * thickness != 0 )
distortion *= 8. / ( element_area * thickness );
else
distortion *= 8.;
}
return (double)distortion;
}
C_FUNC_DEF double v_quad_edge_ratio | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad edge ratio.
the edge ratio of a quad
NB (P. Pebay 01/19/07): Hmax / Hmin where Hmax and Hmin are respectively the maximum and the minimum edge lengths
Definition at line 271 of file V_QuadMetric.cpp.
References VerdictVector::length_squared(), make_quad_edges(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, VERDICT_MAX, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure(), and v_quad_quality().
{
VerdictVector edges[4];
make_quad_edges( edges, coordinates );
double a2 = edges[0].length_squared();
double b2 = edges[1].length_squared();
double c2 = edges[2].length_squared();
double d2 = edges[3].length_squared();
double mab, Mab, mcd, Mcd, m2, M2;
if( a2 < b2 )
{
mab = a2;
Mab = b2;
}
else // b2 <= a2
{
mab = b2;
Mab = a2;
}
if( c2 < d2 )
{
mcd = c2;
Mcd = d2;
}
else // d2 <= c2
{
mcd = d2;
Mcd = c2;
}
m2 = mab < mcd ? mab : mcd;
M2 = Mab > Mcd ? Mab : Mcd;
if( m2 < VERDICT_DBL_MIN )
return (double)VERDICT_DBL_MAX;
else
{
double edge_ratio = sqrt( M2 / m2 );
if( edge_ratio > 0 ) return (double)VERDICT_MIN( edge_ratio, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( edge_ratio, -VERDICT_DBL_MAX );
}
}
C_FUNC_DEF double v_quad_jacobian | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad jacobian.
the jacobian of a quad
minimum pointwise volume of local map at 4 corners and center of quad
Definition at line 870 of file V_QuadMetric.cpp.
References is_collapsed_quad(), signed_corner_areas(), v_tri_area(), VERDICT_DBL_MAX, VERDICT_MAX, VERDICT_MIN, and VERDICT_TRUE.
Referenced by moab::VerdictWrapper::quality_measure().
{
if( is_collapsed_quad( coordinates ) == VERDICT_TRUE ) return (double)( v_tri_area( 3, coordinates ) * 2.0 );
double areas[4];
signed_corner_areas( areas, coordinates );
double jacobian = VERDICT_MIN( VERDICT_MIN( areas[0], areas[1] ), VERDICT_MIN( areas[2], areas[3] ) );
if( jacobian > 0 ) return (double)VERDICT_MIN( jacobian, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( jacobian, -VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_max_aspect_frobenius | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad maximum Frobenius aspect.
the maximum Frobenius aspect of a quad
NB (P. Pebay 01/20/07): this metric is calculated by taking the maximum of the 4 Frobenius aspects at each corner of the quad, when the reference triangle is right isosceles.
Definition at line 484 of file V_QuadMetric.cpp.
References VerdictVector::length(), VerdictVector::length_squared(), make_quad_edges(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, VERDICT_MAX, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure(), and v_quad_quality().
{
VerdictVector edges[4];
make_quad_edges( edges, coordinates );
double a2 = edges[0].length_squared();
double b2 = edges[1].length_squared();
double c2 = edges[2].length_squared();
double d2 = edges[3].length_squared();
VerdictVector ab = edges[0] * edges[1];
VerdictVector bc = edges[1] * edges[2];
VerdictVector cd = edges[2] * edges[3];
VerdictVector da = edges[3] * edges[0];
double ab1 = ab.length();
double bc1 = bc.length();
double cd1 = cd.length();
double da1 = da.length();
if( ab1 < VERDICT_DBL_MIN || bc1 < VERDICT_DBL_MIN || cd1 < VERDICT_DBL_MIN || da1 < VERDICT_DBL_MIN )
return (double)VERDICT_DBL_MAX;
double qmax = ( a2 + b2 ) / ab1;
double qcur = ( b2 + c2 ) / bc1;
qmax = qmax > qcur ? qmax : qcur;
qcur = ( c2 + d2 ) / cd1;
qmax = qmax > qcur ? qmax : qcur;
qcur = ( d2 + a2 ) / da1;
qmax = qmax > qcur ? qmax : qcur;
double max_aspect_frobenius = .5 * qmax;
if( max_aspect_frobenius > 0 ) return (double)VERDICT_MIN( max_aspect_frobenius, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( max_aspect_frobenius, -VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_max_edge_ratio | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad maximum of edge ratio.
maximum of edge ratio of a quad
maximum edge length ratio at quad center
Definition at line 321 of file V_QuadMetric.cpp.
References VerdictVector::length(), VerdictVector::set(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, VERDICT_MAX, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure().
{
VerdictVector quad_nodes[4];
quad_nodes[0].set( coordinates[0][0], coordinates[0][1], coordinates[0][2] );
quad_nodes[1].set( coordinates[1][0], coordinates[1][1], coordinates[1][2] );
quad_nodes[2].set( coordinates[2][0], coordinates[2][1], coordinates[2][2] );
quad_nodes[3].set( coordinates[3][0], coordinates[3][1], coordinates[3][2] );
VerdictVector principal_axes[2];
principal_axes[0] = quad_nodes[1] + quad_nodes[2] - quad_nodes[0] - quad_nodes[3];
principal_axes[1] = quad_nodes[2] + quad_nodes[3] - quad_nodes[0] - quad_nodes[1];
double len1 = principal_axes[0].length();
double len2 = principal_axes[1].length();
if( len1 < VERDICT_DBL_MIN || len2 < VERDICT_DBL_MIN ) return (double)VERDICT_DBL_MAX;
double max_edge_ratio = VERDICT_MAX( len1 / len2, len2 / len1 );
if( max_edge_ratio > 0 ) return (double)VERDICT_MIN( max_edge_ratio, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( max_edge_ratio, -VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_maximum_angle | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad's largest angle.
the largest angle of a quad
largest included quad area (degrees)
Definition at line 670 of file V_QuadMetric.cpp.
References moab::angle(), is_collapsed_quad(), length(), VerdictVector::length(), VerdictVector::set(), signed_corner_areas(), v_tri_maximum_angle(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, VERDICT_MAX, VERDICT_MIN, VERDICT_PI, and VERDICT_TRUE.
Referenced by moab::VerdictWrapper::quality_measure().
{
// if this is a collapsed quad, just pass it on to
// the tri_largest_angle routine
if( is_collapsed_quad( coordinates ) == VERDICT_TRUE ) return v_tri_maximum_angle( 3, coordinates );
double angle;
double max_angle = 0.0;
VerdictVector edges[4];
edges[0].set( coordinates[1][0] - coordinates[0][0], coordinates[1][1] - coordinates[0][1],
coordinates[1][2] - coordinates[0][2] );
edges[1].set( coordinates[2][0] - coordinates[1][0], coordinates[2][1] - coordinates[1][1],
coordinates[2][2] - coordinates[1][2] );
edges[2].set( coordinates[3][0] - coordinates[2][0], coordinates[3][1] - coordinates[2][1],
coordinates[3][2] - coordinates[2][2] );
edges[3].set( coordinates[0][0] - coordinates[3][0], coordinates[0][1] - coordinates[3][1],
coordinates[0][2] - coordinates[3][2] );
// go around each node and calculate the angle
// at each node
double length[4];
length[0] = edges[0].length();
length[1] = edges[1].length();
length[2] = edges[2].length();
length[3] = edges[3].length();
if( length[0] <= VERDICT_DBL_MIN || length[1] <= VERDICT_DBL_MIN || length[2] <= VERDICT_DBL_MIN ||
length[3] <= VERDICT_DBL_MIN )
return 0.0;
angle = acos( -( edges[0] % edges[1] ) / ( length[0] * length[1] ) );
max_angle = VERDICT_MAX( angle, max_angle );
angle = acos( -( edges[1] % edges[2] ) / ( length[1] * length[2] ) );
max_angle = VERDICT_MAX( angle, max_angle );
angle = acos( -( edges[2] % edges[3] ) / ( length[2] * length[3] ) );
max_angle = VERDICT_MAX( angle, max_angle );
angle = acos( -( edges[3] % edges[0] ) / ( length[3] * length[0] ) );
max_angle = VERDICT_MAX( angle, max_angle );
max_angle = max_angle * 180.0 / VERDICT_PI;
// if any signed areas are < 0, then you are getting the wrong angle
double areas[4];
signed_corner_areas( areas, coordinates );
if( areas[0] < 0 || areas[1] < 0 || areas[2] < 0 || areas[3] < 0 )
{
max_angle = 360 - max_angle;
}
if( max_angle > 0 ) return (double)VERDICT_MIN( max_angle, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( max_angle, -VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_med_aspect_frobenius | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad average Frobenius aspect.
the average Frobenius aspect of a quad
NB (P. Pebay 01/20/07): this metric is calculated by averaging the 4 Frobenius aspects at each corner of the quad, when the reference triangle is right isosceles.
Definition at line 442 of file V_QuadMetric.cpp.
References VerdictVector::length(), VerdictVector::length_squared(), make_quad_edges(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, VERDICT_MAX, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure(), and v_quad_quality().
{
VerdictVector edges[4];
make_quad_edges( edges, coordinates );
double a2 = edges[0].length_squared();
double b2 = edges[1].length_squared();
double c2 = edges[2].length_squared();
double d2 = edges[3].length_squared();
VerdictVector ab = edges[0] * edges[1];
VerdictVector bc = edges[1] * edges[2];
VerdictVector cd = edges[2] * edges[3];
VerdictVector da = edges[3] * edges[0];
double ab1 = ab.length();
double bc1 = bc.length();
double cd1 = cd.length();
double da1 = da.length();
if( ab1 < VERDICT_DBL_MIN || bc1 < VERDICT_DBL_MIN || cd1 < VERDICT_DBL_MIN || da1 < VERDICT_DBL_MIN )
return (double)VERDICT_DBL_MAX;
double qsum = ( a2 + b2 ) / ab1;
qsum += ( b2 + c2 ) / bc1;
qsum += ( c2 + d2 ) / cd1;
qsum += ( d2 + a2 ) / da1;
double med_aspect_frobenius = .125 * qsum;
if( med_aspect_frobenius > 0 ) return (double)VERDICT_MIN( med_aspect_frobenius, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( med_aspect_frobenius, -VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_minimum_angle | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad's smallest angle.
the smallest angle of a quad
smallest included quad angle (degrees)
Definition at line 734 of file V_QuadMetric.cpp.
References moab::angle(), is_collapsed_quad(), length(), VerdictVector::length(), VerdictVector::set(), v_tri_minimum_angle(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, VERDICT_MAX, VERDICT_MIN, VERDICT_PI, and VERDICT_TRUE.
Referenced by moab::VerdictWrapper::quality_measure().
{
// if this quad is a collapsed quad, then just
// send it to the tri_smallest_angle routine
if( is_collapsed_quad( coordinates ) == VERDICT_TRUE ) return v_tri_minimum_angle( 3, coordinates );
double angle;
double min_angle = 360.0;
VerdictVector edges[4];
edges[0].set( coordinates[1][0] - coordinates[0][0], coordinates[1][1] - coordinates[0][1],
coordinates[1][2] - coordinates[0][2] );
edges[1].set( coordinates[2][0] - coordinates[1][0], coordinates[2][1] - coordinates[1][1],
coordinates[2][2] - coordinates[1][2] );
edges[2].set( coordinates[3][0] - coordinates[2][0], coordinates[3][1] - coordinates[2][1],
coordinates[3][2] - coordinates[2][2] );
edges[3].set( coordinates[0][0] - coordinates[3][0], coordinates[0][1] - coordinates[3][1],
coordinates[0][2] - coordinates[3][2] );
// go around each node and calculate the angle
// at each node
double length[4];
length[0] = edges[0].length();
length[1] = edges[1].length();
length[2] = edges[2].length();
length[3] = edges[3].length();
if( length[0] <= VERDICT_DBL_MIN || length[1] <= VERDICT_DBL_MIN || length[2] <= VERDICT_DBL_MIN ||
length[3] <= VERDICT_DBL_MIN )
return 360.0;
angle = acos( -( edges[0] % edges[1] ) / ( length[0] * length[1] ) );
min_angle = VERDICT_MIN( angle, min_angle );
angle = acos( -( edges[1] % edges[2] ) / ( length[1] * length[2] ) );
min_angle = VERDICT_MIN( angle, min_angle );
angle = acos( -( edges[2] % edges[3] ) / ( length[2] * length[3] ) );
min_angle = VERDICT_MIN( angle, min_angle );
angle = acos( -( edges[3] % edges[0] ) / ( length[3] * length[0] ) );
min_angle = VERDICT_MIN( angle, min_angle );
min_angle = min_angle * 180.0 / VERDICT_PI;
if( min_angle > 0 ) return (double)VERDICT_MIN( min_angle, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( min_angle, -VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_oddy | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad oddy metric.
the oddy of a quad
general distortion measure based on left Cauchy-Green Tensor
Definition at line 788 of file V_QuadMetric.cpp.
References moab::GeomUtil::first(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, VERDICT_MAX, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure().
{
double max_oddy = 0.;
VerdictVector first, second, node_pos[4];
double g, g11, g12, g22, cur_oddy;
int i;
for( i = 0; i < 4; i++ )
node_pos[i].set( coordinates[i][0], coordinates[i][1], coordinates[i][2] );
for( i = 0; i < 4; i++ )
{
first = node_pos[i] - node_pos[( i + 1 ) % 4];
second = node_pos[i] - node_pos[( i + 3 ) % 4];
g11 = first % first;
g12 = first % second;
g22 = second % second;
g = g11 * g22 - g12 * g12;
if( g < VERDICT_DBL_MIN )
cur_oddy = VERDICT_DBL_MAX;
else
cur_oddy = ( ( g11 - g22 ) * ( g11 - g22 ) + 4. * g12 * g12 ) / 2. / g;
max_oddy = VERDICT_MAX( max_oddy, cur_oddy );
}
if( max_oddy > 0 ) return (double)VERDICT_MIN( max_oddy, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( max_oddy, -VERDICT_DBL_MAX );
}
C_FUNC_DEF void v_quad_quality | ( | int | num_nodes, |
double | coordinates[][3], | ||
unsigned int | metrics_request_flag, | ||
QuadMetricVals * | metric_vals | ||
) |
Calculates quality metrics for quadrilateral elements.
multiple quality measures of a quad
Definition at line 1258 of file V_QuadMetric.cpp.
References QuadMetricVals::area, QuadMetricVals::aspect_ratio, QuadMetricVals::condition, determinant(), QuadMetricVals::distortion, QuadMetricVals::edge_ratio, get_weight(), is_collapsed_quad(), QuadMetricVals::jacobian, VerdictVector::length(), length_squared(), VerdictVector::length_squared(), make_quad_edges(), QuadMetricVals::max_aspect_frobenius, QuadMetricVals::max_edge_ratio, QuadMetricVals::maximum_angle, QuadMetricVals::med_aspect_frobenius, QuadMetricVals::minimum_angle, normalize(), QuadMetricVals::oddy, QuadMetricVals::radius_ratio, QuadMetricVals::relative_size_squared, QuadMetricVals::scaled_jacobian, VerdictVector::set(), QuadMetricVals::shape, QuadMetricVals::shape_and_size, QuadMetricVals::shear, QuadMetricVals::shear_and_size, signed_corner_areas(), QuadMetricVals::skew, QuadMetricVals::stretch, QuadMetricVals::taper, V_QUAD_AREA, v_quad_area(), V_QUAD_ASPECT_RATIO, v_quad_aspect_ratio(), V_QUAD_CONDITION, V_QUAD_DISTORTION, v_quad_distortion(), V_QUAD_EDGE_RATIO, v_quad_edge_ratio(), V_QUAD_JACOBIAN, V_QUAD_MAX_ASPECT_FROBENIUS, v_quad_max_aspect_frobenius(), V_QUAD_MAX_EDGE_RATIO, V_QUAD_MAXIMUM_ANGLE, V_QUAD_MED_ASPECT_FROBENIUS, v_quad_med_aspect_frobenius(), V_QUAD_MINIMUM_ANGLE, V_QUAD_ODDY, V_QUAD_RADIUS_RATIO, v_quad_radius_ratio(), V_QUAD_RELATIVE_SIZE_SQUARED, V_QUAD_SCALED_JACOBIAN, V_QUAD_SHAPE, V_QUAD_SHAPE_AND_SIZE, V_QUAD_SHEAR, V_QUAD_SHEAR_AND_SIZE, V_QUAD_SKEW, V_QUAD_STRETCH, V_QUAD_TAPER, V_QUAD_WARPAGE, v_set_quad_size(), v_tri_area(), v_tri_maximum_angle(), v_tri_minimum_angle(), v_tri_scaled_jacobian(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, VERDICT_FALSE, VERDICT_MAX, VERDICT_MIN, VERDICT_PI, VERDICT_TRUE, and QuadMetricVals::warpage.
Referenced by moab::VerdictWrapper::all_quality_measures().
{
memset( metric_vals, 0, sizeof( QuadMetricVals ) );
// for starts, lets set up some basic and common information
/* node numbers and side numbers used below
2
3 +--------- 2
/ +
/ |
3 / | 1
/ |
+ |
0 -------------+ 1
0
*/
// vectors for each side
VerdictVector edges[4];
make_quad_edges( edges, coordinates );
double areas[4];
signed_corner_areas( areas, coordinates );
double lengths[4];
lengths[0] = edges[0].length();
lengths[1] = edges[1].length();
lengths[2] = edges[2].length();
lengths[3] = edges[3].length();
VerdictBoolean is_collapsed = is_collapsed_quad( coordinates );
// handle collapsed quads metrics here
if( is_collapsed == VERDICT_TRUE && metrics_request_flag & ( V_QUAD_MINIMUM_ANGLE | V_QUAD_MAXIMUM_ANGLE |
V_QUAD_JACOBIAN | V_QUAD_SCALED_JACOBIAN ) )
{
if( metrics_request_flag & V_QUAD_MINIMUM_ANGLE )
metric_vals->minimum_angle = v_tri_minimum_angle( 3, coordinates );
if( metrics_request_flag & V_QUAD_MAXIMUM_ANGLE )
metric_vals->maximum_angle = v_tri_maximum_angle( 3, coordinates );
if( metrics_request_flag & V_QUAD_JACOBIAN )
metric_vals->jacobian = (double)( v_tri_area( 3, coordinates ) * 2.0 );
if( metrics_request_flag & V_QUAD_SCALED_JACOBIAN )
metric_vals->jacobian = (double)( v_tri_scaled_jacobian( 3, coordinates ) * 2.0 );
}
// calculate both largest and smallest angles
if( metrics_request_flag & ( V_QUAD_MINIMUM_ANGLE | V_QUAD_MAXIMUM_ANGLE ) && is_collapsed == VERDICT_FALSE )
{
// gather the angles
double angles[4];
angles[0] = acos( -( edges[0] % edges[1] ) / ( lengths[0] * lengths[1] ) );
angles[1] = acos( -( edges[1] % edges[2] ) / ( lengths[1] * lengths[2] ) );
angles[2] = acos( -( edges[2] % edges[3] ) / ( lengths[2] * lengths[3] ) );
angles[3] = acos( -( edges[3] % edges[0] ) / ( lengths[3] * lengths[0] ) );
if( lengths[0] <= VERDICT_DBL_MIN || lengths[1] <= VERDICT_DBL_MIN || lengths[2] <= VERDICT_DBL_MIN ||
lengths[3] <= VERDICT_DBL_MIN )
{
metric_vals->minimum_angle = 360.0;
metric_vals->maximum_angle = 0.0;
}
else
{
// if smallest angle, find the smallest angle
if( metrics_request_flag & V_QUAD_MINIMUM_ANGLE )
{
metric_vals->minimum_angle = VERDICT_DBL_MAX;
for( int i = 0; i < 4; i++ )
metric_vals->minimum_angle = VERDICT_MIN( angles[i], metric_vals->minimum_angle );
metric_vals->minimum_angle *= 180.0 / VERDICT_PI;
}
// if largest angle, find the largest angle
if( metrics_request_flag & V_QUAD_MAXIMUM_ANGLE )
{
metric_vals->maximum_angle = 0.0;
for( int i = 0; i < 4; i++ )
metric_vals->maximum_angle = VERDICT_MAX( angles[i], metric_vals->maximum_angle );
metric_vals->maximum_angle *= 180.0 / VERDICT_PI;
if( areas[0] < 0 || areas[1] < 0 || areas[2] < 0 || areas[3] < 0 )
metric_vals->maximum_angle = 360 - metric_vals->maximum_angle;
}
}
}
// handle max_edge_ratio, skew, taper, and area together
if( metrics_request_flag & ( V_QUAD_MAX_EDGE_RATIO | V_QUAD_SKEW | V_QUAD_TAPER ) )
{
// get principle axes
VerdictVector principal_axes[2];
principal_axes[0] = edges[0] - edges[2];
principal_axes[1] = edges[1] - edges[3];
if( metrics_request_flag & ( V_QUAD_MAX_EDGE_RATIO | V_QUAD_SKEW | V_QUAD_TAPER ) )
{
double len1 = principal_axes[0].length();
double len2 = principal_axes[1].length();
// calculate the max_edge_ratio ratio
if( metrics_request_flag & V_QUAD_MAX_EDGE_RATIO )
{
if( len1 < VERDICT_DBL_MIN || len2 < VERDICT_DBL_MIN )
metric_vals->max_edge_ratio = VERDICT_DBL_MAX;
else
metric_vals->max_edge_ratio = VERDICT_MAX( len1 / len2, len2 / len1 );
}
// calculate the taper
if( metrics_request_flag & V_QUAD_TAPER )
{
double min_length = VERDICT_MIN( len1, len2 );
VerdictVector cross_derivative = edges[1] + edges[3];
if( min_length < VERDICT_DBL_MIN )
metric_vals->taper = VERDICT_DBL_MAX;
else
metric_vals->taper = cross_derivative.length() / min_length;
}
// calculate the skew
if( metrics_request_flag & V_QUAD_SKEW )
{
if( principal_axes[0].normalize() < VERDICT_DBL_MIN || principal_axes[1].normalize() < VERDICT_DBL_MIN )
metric_vals->skew = 0.0;
else
metric_vals->skew = fabs( principal_axes[0] % principal_axes[1] );
}
}
}
// calculate the area
if( metrics_request_flag & ( V_QUAD_AREA | V_QUAD_RELATIVE_SIZE_SQUARED ) )
{
metric_vals->area = 0.25 * ( areas[0] + areas[1] + areas[2] + areas[3] );
}
// calculate the relative size
if( metrics_request_flag & ( V_QUAD_RELATIVE_SIZE_SQUARED | V_QUAD_SHAPE_AND_SIZE | V_QUAD_SHEAR_AND_SIZE ) )
{
double quad_area = v_quad_area( 4, coordinates );
v_set_quad_size( quad_area );
double w11, w21, w12, w22;
get_weight( w11, w21, w12, w22 );
double avg_area = determinant( w11, w21, w12, w22 );
if( avg_area < VERDICT_DBL_MIN )
metric_vals->relative_size_squared = 0.0;
else
metric_vals->relative_size_squared =
pow( VERDICT_MIN( metric_vals->area / avg_area, avg_area / metric_vals->area ), 2 );
}
// calculate the jacobian
if( metrics_request_flag & V_QUAD_JACOBIAN )
{
metric_vals->jacobian = VERDICT_MIN( VERDICT_MIN( areas[0], areas[1] ), VERDICT_MIN( areas[2], areas[3] ) );
}
if( metrics_request_flag & ( V_QUAD_SCALED_JACOBIAN | V_QUAD_SHEAR | V_QUAD_SHEAR_AND_SIZE ) )
{
double scaled_jac, min_scaled_jac = VERDICT_DBL_MAX;
if( lengths[0] < VERDICT_DBL_MIN || lengths[1] < VERDICT_DBL_MIN || lengths[2] < VERDICT_DBL_MIN ||
lengths[3] < VERDICT_DBL_MIN )
{
metric_vals->scaled_jacobian = 0.0;
metric_vals->shear = 0.0;
}
else
{
scaled_jac = areas[0] / ( lengths[0] * lengths[3] );
min_scaled_jac = VERDICT_MIN( scaled_jac, min_scaled_jac );
scaled_jac = areas[1] / ( lengths[1] * lengths[0] );
min_scaled_jac = VERDICT_MIN( scaled_jac, min_scaled_jac );
scaled_jac = areas[2] / ( lengths[2] * lengths[1] );
min_scaled_jac = VERDICT_MIN( scaled_jac, min_scaled_jac );
scaled_jac = areas[3] / ( lengths[3] * lengths[2] );
min_scaled_jac = VERDICT_MIN( scaled_jac, min_scaled_jac );
metric_vals->scaled_jacobian = min_scaled_jac;
// what the heck...set shear as well
if( min_scaled_jac <= VERDICT_DBL_MIN )
metric_vals->shear = 0.0;
else
metric_vals->shear = min_scaled_jac;
}
}
if( metrics_request_flag & ( V_QUAD_WARPAGE | V_QUAD_ODDY ) )
{
VerdictVector corner_normals[4];
corner_normals[0] = edges[3] * edges[0];
corner_normals[1] = edges[0] * edges[1];
corner_normals[2] = edges[1] * edges[2];
corner_normals[3] = edges[2] * edges[3];
if( metrics_request_flag & V_QUAD_ODDY )
{
double oddy, max_oddy = 0.0;
double diff, dot_prod;
double length_squared[4];
length_squared[0] = corner_normals[0].length_squared();
length_squared[1] = corner_normals[1].length_squared();
length_squared[2] = corner_normals[2].length_squared();
length_squared[3] = corner_normals[3].length_squared();
if( length_squared[0] < VERDICT_DBL_MIN || length_squared[1] < VERDICT_DBL_MIN ||
length_squared[2] < VERDICT_DBL_MIN || length_squared[3] < VERDICT_DBL_MIN )
metric_vals->oddy = VERDICT_DBL_MAX;
else
{
diff = ( lengths[0] * lengths[0] ) - ( lengths[1] * lengths[1] );
dot_prod = edges[0] % edges[1];
oddy = ( ( diff * diff ) + 4 * dot_prod * dot_prod ) / ( 2 * length_squared[1] );
max_oddy = VERDICT_MAX( oddy, max_oddy );
diff = ( lengths[1] * lengths[1] ) - ( lengths[2] * lengths[2] );
dot_prod = edges[1] % edges[2];
oddy = ( ( diff * diff ) + 4 * dot_prod * dot_prod ) / ( 2 * length_squared[2] );
max_oddy = VERDICT_MAX( oddy, max_oddy );
diff = ( lengths[2] * lengths[2] ) - ( lengths[3] * lengths[3] );
dot_prod = edges[2] % edges[3];
oddy = ( ( diff * diff ) + 4 * dot_prod * dot_prod ) / ( 2 * length_squared[3] );
max_oddy = VERDICT_MAX( oddy, max_oddy );
diff = ( lengths[3] * lengths[3] ) - ( lengths[0] * lengths[0] );
dot_prod = edges[3] % edges[0];
oddy = ( ( diff * diff ) + 4 * dot_prod * dot_prod ) / ( 2 * length_squared[0] );
max_oddy = VERDICT_MAX( oddy, max_oddy );
metric_vals->oddy = max_oddy;
}
}
if( metrics_request_flag & V_QUAD_WARPAGE )
{
if( corner_normals[0].normalize() < VERDICT_DBL_MIN || corner_normals[1].normalize() < VERDICT_DBL_MIN ||
corner_normals[2].normalize() < VERDICT_DBL_MIN || corner_normals[3].normalize() < VERDICT_DBL_MIN )
metric_vals->warpage = VERDICT_DBL_MAX;
else
{
metric_vals->warpage =
pow( VERDICT_MIN( corner_normals[0] % corner_normals[2], corner_normals[1] % corner_normals[3] ),
3 );
}
}
}
if( metrics_request_flag & V_QUAD_STRETCH )
{
VerdictVector temp;
temp.set( coordinates[2][0] - coordinates[0][0], coordinates[2][1] - coordinates[0][1],
coordinates[2][2] - coordinates[0][2] );
double diag02 = temp.length_squared();
temp.set( coordinates[3][0] - coordinates[1][0], coordinates[3][1] - coordinates[1][1],
coordinates[3][2] - coordinates[1][2] );
double diag13 = temp.length_squared();
static const double QUAD_STRETCH_FACTOR = sqrt( 2.0 );
// 'diag02' is now the max diagonal of the quad
diag02 = VERDICT_MAX( diag02, diag13 );
if( diag02 < VERDICT_DBL_MIN )
metric_vals->stretch = VERDICT_DBL_MAX;
else
metric_vals->stretch =
QUAD_STRETCH_FACTOR *
VERDICT_MIN( VERDICT_MIN( lengths[0], lengths[1] ), VERDICT_MIN( lengths[2], lengths[3] ) ) /
sqrt( diag02 );
}
if( metrics_request_flag & ( V_QUAD_CONDITION | V_QUAD_SHAPE | V_QUAD_SHAPE_AND_SIZE ) )
{
double lengths_squared[4];
lengths_squared[0] = edges[0].length_squared();
lengths_squared[1] = edges[1].length_squared();
lengths_squared[2] = edges[2].length_squared();
lengths_squared[3] = edges[3].length_squared();
if( areas[0] < VERDICT_DBL_MIN || areas[1] < VERDICT_DBL_MIN || areas[2] < VERDICT_DBL_MIN ||
areas[3] < VERDICT_DBL_MIN )
{
metric_vals->condition = VERDICT_DBL_MAX;
metric_vals->shape = 0.0;
}
else
{
double max_condition = 0.0, condition;
condition = ( lengths_squared[0] + lengths_squared[3] ) / areas[0];
max_condition = VERDICT_MAX( max_condition, condition );
condition = ( lengths_squared[1] + lengths_squared[0] ) / areas[1];
max_condition = VERDICT_MAX( max_condition, condition );
condition = ( lengths_squared[2] + lengths_squared[1] ) / areas[2];
max_condition = VERDICT_MAX( max_condition, condition );
condition = ( lengths_squared[3] + lengths_squared[2] ) / areas[3];
max_condition = VERDICT_MAX( max_condition, condition );
metric_vals->condition = 0.5 * max_condition;
metric_vals->shape = 2 / max_condition;
}
}
if( metrics_request_flag & V_QUAD_AREA )
{
if( metric_vals->area > 0 ) metric_vals->area = (double)VERDICT_MIN( metric_vals->area, VERDICT_DBL_MAX );
metric_vals->area = (double)VERDICT_MAX( metric_vals->area, -VERDICT_DBL_MAX );
}
if( metrics_request_flag & V_QUAD_MAX_EDGE_RATIO )
{
if( metric_vals->max_edge_ratio > 0 )
metric_vals->max_edge_ratio = (double)VERDICT_MIN( metric_vals->max_edge_ratio, VERDICT_DBL_MAX );
metric_vals->max_edge_ratio = (double)VERDICT_MAX( metric_vals->max_edge_ratio, -VERDICT_DBL_MAX );
}
if( metrics_request_flag & V_QUAD_CONDITION )
{
if( metric_vals->condition > 0 )
metric_vals->condition = (double)VERDICT_MIN( metric_vals->condition, VERDICT_DBL_MAX );
metric_vals->condition = (double)VERDICT_MAX( metric_vals->condition, -VERDICT_DBL_MAX );
}
// calculate distortion
if( metrics_request_flag & V_QUAD_DISTORTION )
{
metric_vals->distortion = v_quad_distortion( num_nodes, coordinates );
if( metric_vals->distortion > 0 )
metric_vals->distortion = (double)VERDICT_MIN( metric_vals->distortion, VERDICT_DBL_MAX );
metric_vals->distortion = (double)VERDICT_MAX( metric_vals->distortion, -VERDICT_DBL_MAX );
}
if( metrics_request_flag & V_QUAD_JACOBIAN )
{
if( metric_vals->jacobian > 0 )
metric_vals->jacobian = (double)VERDICT_MIN( metric_vals->jacobian, VERDICT_DBL_MAX );
metric_vals->jacobian = (double)VERDICT_MAX( metric_vals->jacobian, -VERDICT_DBL_MAX );
}
if( metrics_request_flag & V_QUAD_MAXIMUM_ANGLE )
{
if( metric_vals->maximum_angle > 0 )
metric_vals->maximum_angle = (double)VERDICT_MIN( metric_vals->maximum_angle, VERDICT_DBL_MAX );
metric_vals->maximum_angle = (double)VERDICT_MAX( metric_vals->maximum_angle, -VERDICT_DBL_MAX );
}
if( metrics_request_flag & V_QUAD_MINIMUM_ANGLE )
{
if( metric_vals->minimum_angle > 0 )
metric_vals->minimum_angle = (double)VERDICT_MIN( metric_vals->minimum_angle, VERDICT_DBL_MAX );
metric_vals->minimum_angle = (double)VERDICT_MAX( metric_vals->minimum_angle, -VERDICT_DBL_MAX );
}
if( metrics_request_flag & V_QUAD_ODDY )
{
if( metric_vals->oddy > 0 ) metric_vals->oddy = (double)VERDICT_MIN( metric_vals->oddy, VERDICT_DBL_MAX );
metric_vals->oddy = (double)VERDICT_MAX( metric_vals->oddy, -VERDICT_DBL_MAX );
}
if( metrics_request_flag & V_QUAD_RELATIVE_SIZE_SQUARED )
{
if( metric_vals->relative_size_squared > 0 )
metric_vals->relative_size_squared =
(double)VERDICT_MIN( metric_vals->relative_size_squared, VERDICT_DBL_MAX );
metric_vals->relative_size_squared =
(double)VERDICT_MAX( metric_vals->relative_size_squared, -VERDICT_DBL_MAX );
}
if( metrics_request_flag & V_QUAD_SCALED_JACOBIAN )
{
if( metric_vals->scaled_jacobian > 0 )
metric_vals->scaled_jacobian = (double)VERDICT_MIN( metric_vals->scaled_jacobian, VERDICT_DBL_MAX );
metric_vals->scaled_jacobian = (double)VERDICT_MAX( metric_vals->scaled_jacobian, -VERDICT_DBL_MAX );
}
if( metrics_request_flag & V_QUAD_SHEAR )
{
if( metric_vals->shear > 0 ) metric_vals->shear = (double)VERDICT_MIN( metric_vals->shear, VERDICT_DBL_MAX );
metric_vals->shear = (double)VERDICT_MAX( metric_vals->shear, -VERDICT_DBL_MAX );
}
// calculate shear and size
// reuse values from above
if( metrics_request_flag & V_QUAD_SHEAR_AND_SIZE )
{
metric_vals->shear_and_size = metric_vals->shear * metric_vals->relative_size_squared;
if( metric_vals->shear_and_size > 0 )
metric_vals->shear_and_size = (double)VERDICT_MIN( metric_vals->shear_and_size, VERDICT_DBL_MAX );
metric_vals->shear_and_size = (double)VERDICT_MAX( metric_vals->shear_and_size, -VERDICT_DBL_MAX );
}
if( metrics_request_flag & V_QUAD_SHAPE )
{
if( metric_vals->shape > 0 ) metric_vals->shape = (double)VERDICT_MIN( metric_vals->shape, VERDICT_DBL_MAX );
metric_vals->shape = (double)VERDICT_MAX( metric_vals->shape, -VERDICT_DBL_MAX );
}
// calculate shape and size
// reuse values from above
if( metrics_request_flag & V_QUAD_SHAPE_AND_SIZE )
{
metric_vals->shape_and_size = metric_vals->shape * metric_vals->relative_size_squared;
if( metric_vals->shape_and_size > 0 )
metric_vals->shape_and_size = (double)VERDICT_MIN( metric_vals->shape_and_size, VERDICT_DBL_MAX );
metric_vals->shape_and_size = (double)VERDICT_MAX( metric_vals->shape_and_size, -VERDICT_DBL_MAX );
}
if( metrics_request_flag & V_QUAD_SKEW )
{
if( metric_vals->skew > 0 ) metric_vals->skew = (double)VERDICT_MIN( metric_vals->skew, VERDICT_DBL_MAX );
metric_vals->skew = (double)VERDICT_MAX( metric_vals->skew, -VERDICT_DBL_MAX );
}
if( metrics_request_flag & V_QUAD_STRETCH )
{
if( metric_vals->stretch > 0 )
metric_vals->stretch = (double)VERDICT_MIN( metric_vals->stretch, VERDICT_DBL_MAX );
metric_vals->stretch = (double)VERDICT_MAX( metric_vals->stretch, -VERDICT_DBL_MAX );
}
if( metrics_request_flag & V_QUAD_TAPER )
{
if( metric_vals->taper > 0 ) metric_vals->taper = (double)VERDICT_MIN( metric_vals->taper, VERDICT_DBL_MAX );
metric_vals->taper = (double)VERDICT_MAX( metric_vals->taper, -VERDICT_DBL_MAX );
}
if( metrics_request_flag & V_QUAD_WARPAGE )
{
if( metric_vals->warpage > 0 )
metric_vals->warpage = (double)VERDICT_MIN( metric_vals->warpage, VERDICT_DBL_MAX );
metric_vals->warpage = (double)VERDICT_MAX( metric_vals->warpage, -VERDICT_DBL_MAX );
}
if( metrics_request_flag & V_QUAD_MED_ASPECT_FROBENIUS )
metric_vals->med_aspect_frobenius = v_quad_med_aspect_frobenius( 4, coordinates );
if( metrics_request_flag & V_QUAD_MAX_ASPECT_FROBENIUS )
metric_vals->max_aspect_frobenius = v_quad_max_aspect_frobenius( 4, coordinates );
if( metrics_request_flag & V_QUAD_EDGE_RATIO ) metric_vals->edge_ratio = v_quad_edge_ratio( 4, coordinates );
if( metrics_request_flag & V_QUAD_ASPECT_RATIO ) metric_vals->aspect_ratio = v_quad_aspect_ratio( 4, coordinates );
if( metrics_request_flag & V_QUAD_RADIUS_RATIO ) metric_vals->radius_ratio = v_quad_radius_ratio( 4, coordinates );
}
C_FUNC_DEF double v_quad_radius_ratio | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad radius ratio.
the radius ratio of a quad
NB (P. Pebay 01/19/07): this metric is called "radius ratio" by extension of a concept that does not exist in general with quads -- although a different name should probably be used in the future.
Definition at line 386 of file V_QuadMetric.cpp.
References VerdictVector::length(), VerdictVector::length_squared(), make_quad_edges(), VerdictVector::set(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, VERDICT_MAX, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure(), and v_quad_quality().
{
static const double normal_coeff = 1. / ( 2. * sqrt( 2. ) );
VerdictVector edges[4];
make_quad_edges( edges, coordinates );
double a2 = edges[0].length_squared();
double b2 = edges[1].length_squared();
double c2 = edges[2].length_squared();
double d2 = edges[3].length_squared();
VerdictVector diag;
diag.set( coordinates[2][0] - coordinates[0][0], coordinates[2][1] - coordinates[0][1],
coordinates[2][2] - coordinates[0][2] );
double m2 = diag.length_squared();
diag.set( coordinates[3][0] - coordinates[1][0], coordinates[3][1] - coordinates[1][1],
coordinates[3][2] - coordinates[1][2] );
double n2 = diag.length_squared();
double t0 = a2 > b2 ? a2 : b2;
double t1 = c2 > d2 ? c2 : d2;
double t2 = m2 > n2 ? m2 : n2;
double h2 = t0 > t1 ? t0 : t1;
h2 = h2 > t2 ? h2 : t2;
VerdictVector ab = edges[0] * edges[1];
VerdictVector bc = edges[1] * edges[2];
VerdictVector cd = edges[2] * edges[3];
VerdictVector da = edges[3] * edges[0];
t0 = da.length();
t1 = ab.length();
t2 = bc.length();
double t3 = cd.length();
t0 = t0 < t1 ? t0 : t1;
t2 = t2 < t3 ? t2 : t3;
t0 = t0 < t2 ? t0 : t2;
if( t0 < VERDICT_DBL_MIN ) return (double)VERDICT_DBL_MAX;
double radius_ratio = normal_coeff * sqrt( ( a2 + b2 + c2 + d2 ) * h2 ) / t0;
if( radius_ratio > 0 ) return (double)VERDICT_MIN( radius_ratio, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( radius_ratio, -VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_relative_size_squared | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad relative size metric.
the relative size of a quad
Min( J, 1/J ), where J is determinant of weighted Jacobian matrix
Definition at line 988 of file V_QuadMetric.cpp.
References determinant(), get_weight(), v_quad_area(), v_set_quad_size(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, VERDICT_MAX, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure(), v_quad_shape_and_size(), and v_quad_shear_and_size().
{
double quad_area = v_quad_area( 4, coordinates );
double rel_size = 0;
v_set_quad_size( quad_area );
double w11, w21, w12, w22;
get_weight( w11, w21, w12, w22 );
double avg_area = determinant( w11, w21, w12, w22 );
if( avg_area > VERDICT_DBL_MIN )
{
w11 = quad_area / avg_area;
if( w11 > VERDICT_DBL_MIN )
{
rel_size = VERDICT_MIN( w11, 1 / w11 );
rel_size *= rel_size;
}
}
if( rel_size > 0 ) return (double)VERDICT_MIN( rel_size, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( rel_size, -VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_scaled_jacobian | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad scaled jacobian.
scaled jacobian of a quad
Minimum Jacobian divided by the lengths of the 2 edge vector
Definition at line 888 of file V_QuadMetric.cpp.
References is_collapsed_quad(), length(), VerdictVector::length(), make_quad_edges(), signed_corner_areas(), v_tri_scaled_jacobian(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, VERDICT_MAX, VERDICT_MIN, and VERDICT_TRUE.
Referenced by moab::VerdictWrapper::quality_measure(), and v_quad_shear().
{
if( is_collapsed_quad( coordinates ) == VERDICT_TRUE ) return v_tri_scaled_jacobian( 3, coordinates );
double corner_areas[4], min_scaled_jac = VERDICT_DBL_MAX, scaled_jac;
signed_corner_areas( corner_areas, coordinates );
VerdictVector edges[4];
make_quad_edges( edges, coordinates );
double length[4];
length[0] = edges[0].length();
length[1] = edges[1].length();
length[2] = edges[2].length();
length[3] = edges[3].length();
if( length[0] < VERDICT_DBL_MIN || length[1] < VERDICT_DBL_MIN || length[2] < VERDICT_DBL_MIN ||
length[3] < VERDICT_DBL_MIN )
return 0.0;
scaled_jac = corner_areas[0] / ( length[0] * length[3] );
min_scaled_jac = VERDICT_MIN( scaled_jac, min_scaled_jac );
scaled_jac = corner_areas[1] / ( length[1] * length[0] );
min_scaled_jac = VERDICT_MIN( scaled_jac, min_scaled_jac );
scaled_jac = corner_areas[2] / ( length[2] * length[1] );
min_scaled_jac = VERDICT_MIN( scaled_jac, min_scaled_jac );
scaled_jac = corner_areas[3] / ( length[3] * length[2] );
min_scaled_jac = VERDICT_MIN( scaled_jac, min_scaled_jac );
if( min_scaled_jac > 0 ) return (double)VERDICT_MIN( min_scaled_jac, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( min_scaled_jac, -VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_shape | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad shape metric.
the shape of a quad
2/Condition number of weighted Jacobian matrix
Definition at line 944 of file V_QuadMetric.cpp.
References length_squared(), VerdictVector::length_squared(), make_quad_edges(), signed_corner_areas(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, VERDICT_MAX, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure(), and v_quad_shape_and_size().
{
double corner_areas[4], min_shape = VERDICT_DBL_MAX, shape;
signed_corner_areas( corner_areas, coordinates );
VerdictVector edges[4];
make_quad_edges( edges, coordinates );
double length_squared[4];
length_squared[0] = edges[0].length_squared();
length_squared[1] = edges[1].length_squared();
length_squared[2] = edges[2].length_squared();
length_squared[3] = edges[3].length_squared();
if( length_squared[0] <= VERDICT_DBL_MIN || length_squared[1] <= VERDICT_DBL_MIN ||
length_squared[2] <= VERDICT_DBL_MIN || length_squared[3] <= VERDICT_DBL_MIN )
return 0.0;
shape = corner_areas[0] / ( length_squared[0] + length_squared[3] );
min_shape = VERDICT_MIN( shape, min_shape );
shape = corner_areas[1] / ( length_squared[1] + length_squared[0] );
min_shape = VERDICT_MIN( shape, min_shape );
shape = corner_areas[2] / ( length_squared[2] + length_squared[1] );
min_shape = VERDICT_MIN( shape, min_shape );
shape = corner_areas[3] / ( length_squared[3] + length_squared[2] );
min_shape = VERDICT_MIN( shape, min_shape );
min_shape *= 2;
if( min_shape < VERDICT_DBL_MIN ) min_shape = 0;
if( min_shape > 0 ) return (double)VERDICT_MIN( min_shape, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( min_shape, -VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_shape_and_size | ( | int | num_nodes, |
double | coordinates[][3] | ||
) |
Calculates quad shape-size metric.
the relative shape and size of a quad
Product of Shape and Relative Size
Definition at line 1020 of file V_QuadMetric.cpp.
References size, v_quad_relative_size_squared(), v_quad_shape(), VERDICT_DBL_MAX, VERDICT_MAX, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure().
{
double shape, size;
size = v_quad_relative_size_squared( num_nodes, coordinates );
shape = v_quad_shape( num_nodes, coordinates );
double shape_and_size = shape * size;
if( shape_and_size > 0 ) return (double)VERDICT_MIN( shape_and_size, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( shape_and_size, -VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_shear | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad shear metric.
the shear of a quad
2/Condition number of Jacobian Skew matrix
Definition at line 929 of file V_QuadMetric.cpp.
References v_quad_scaled_jacobian(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure(), and v_quad_shear_and_size().
{
double scaled_jacobian = v_quad_scaled_jacobian( 4, coordinates );
if( scaled_jacobian <= VERDICT_DBL_MIN )
return 0.0;
else
return (double)VERDICT_MIN( scaled_jacobian, VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_shear_and_size | ( | int | num_nodes, |
double | coordinates[][3] | ||
) |
Calculates quad shear-size metric.
the shear and size of a quad
product of shear and relative size
Definition at line 1037 of file V_QuadMetric.cpp.
References size, v_quad_relative_size_squared(), v_quad_shear(), VERDICT_DBL_MAX, VERDICT_MAX, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure().
{
double shear, size;
shear = v_quad_shear( num_nodes, coordinates );
size = v_quad_relative_size_squared( num_nodes, coordinates );
double shear_and_size = shear * size;
if( shear_and_size > 0 ) return (double)VERDICT_MIN( shear_and_size, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( shear_and_size, -VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_skew | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad skew metric.
skew of a quad
maximum ||cos A|| where A is the angle between edges at quad center
Definition at line 530 of file V_QuadMetric.cpp.
References normalize(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure().
{
VerdictVector node_pos[4];
for( int i = 0; i < 4; i++ )
node_pos[i].set( coordinates[i][0], coordinates[i][1], coordinates[i][2] );
VerdictVector principle_axes[2];
principle_axes[0] = node_pos[1] + node_pos[2] - node_pos[3] - node_pos[0];
principle_axes[1] = node_pos[2] + node_pos[3] - node_pos[0] - node_pos[1];
if( principle_axes[0].normalize() < VERDICT_DBL_MIN ) return 0.0;
if( principle_axes[1].normalize() < VERDICT_DBL_MIN ) return 0.0;
double skew = fabs( principle_axes[0] % principle_axes[1] );
return (double)VERDICT_MIN( skew, VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_stretch | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad strech metric.
the stretch of a quad
sqrt(2) * minimum edge length / maximum diagonal length
Definition at line 628 of file V_QuadMetric.cpp.
References VerdictVector::length_squared(), make_quad_edges(), VerdictVector::set(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, VERDICT_MAX, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure().
{
VerdictVector edges[4], temp;
make_quad_edges( edges, coordinates );
double lengths_squared[4];
lengths_squared[0] = edges[0].length_squared();
lengths_squared[1] = edges[1].length_squared();
lengths_squared[2] = edges[2].length_squared();
lengths_squared[3] = edges[3].length_squared();
temp.set( coordinates[2][0] - coordinates[0][0], coordinates[2][1] - coordinates[0][1],
coordinates[2][2] - coordinates[0][2] );
double diag02 = temp.length_squared();
temp.set( coordinates[3][0] - coordinates[1][0], coordinates[3][1] - coordinates[1][1],
coordinates[3][2] - coordinates[1][2] );
double diag13 = temp.length_squared();
static const double QUAD_STRETCH_FACTOR = sqrt( 2.0 );
// 'diag02' is now the max diagonal of the quad
diag02 = VERDICT_MAX( diag02, diag13 );
if( diag02 < VERDICT_DBL_MIN )
return (double)VERDICT_DBL_MAX;
else
{
double stretch =
(double)( QUAD_STRETCH_FACTOR * sqrt( VERDICT_MIN( VERDICT_MIN( lengths_squared[0], lengths_squared[1] ),
VERDICT_MIN( lengths_squared[2], lengths_squared[3] ) ) /
diag02 ) );
return (double)VERDICT_MIN( stretch, VERDICT_DBL_MAX );
}
}
C_FUNC_DEF double v_quad_taper | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad taper metric.
taper of a quad
maximum ratio of lengths derived from opposite edges
Definition at line 553 of file V_QuadMetric.cpp.
References VerdictVector::length(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure().
{
VerdictVector node_pos[4];
for( int i = 0; i < 4; i++ )
node_pos[i].set( coordinates[i][0], coordinates[i][1], coordinates[i][2] );
VerdictVector principle_axes[2];
principle_axes[0] = node_pos[1] + node_pos[2] - node_pos[3] - node_pos[0];
principle_axes[1] = node_pos[2] + node_pos[3] - node_pos[0] - node_pos[1];
VerdictVector cross_derivative = node_pos[0] + node_pos[2] - node_pos[1] - node_pos[3];
double lengths[2];
lengths[0] = principle_axes[0].length();
lengths[1] = principle_axes[1].length();
// get min length
lengths[0] = VERDICT_MIN( lengths[0], lengths[1] );
if( lengths[0] < VERDICT_DBL_MIN ) return VERDICT_DBL_MAX;
double taper = cross_derivative.length() / lengths[0];
return (double)VERDICT_MIN( taper, VERDICT_DBL_MAX );
}
C_FUNC_DEF double v_quad_warpage | ( | int | , |
double | coordinates[][3] | ||
) |
Calculates quad warpage metric.
warpage of a quad
deviation of element from planarity
Definition at line 583 of file V_QuadMetric.cpp.
References make_quad_edges(), normalize(), VERDICT_DBL_MAX, VERDICT_DBL_MIN, VERDICT_MAX, and VERDICT_MIN.
Referenced by moab::VerdictWrapper::quality_measure().
{
VerdictVector edges[4];
make_quad_edges( edges, coordinates );
VerdictVector corner_normals[4];
corner_normals[0] = edges[3] * edges[0];
corner_normals[1] = edges[0] * edges[1];
corner_normals[2] = edges[1] * edges[2];
corner_normals[3] = edges[2] * edges[3];
if( corner_normals[0].normalize() < VERDICT_DBL_MIN || corner_normals[1].normalize() < VERDICT_DBL_MIN ||
corner_normals[2].normalize() < VERDICT_DBL_MIN || corner_normals[3].normalize() < VERDICT_DBL_MIN )
return (double)VERDICT_DBL_MIN;
double warpage =
pow( VERDICT_MIN( corner_normals[0] % corner_normals[2], corner_normals[1] % corner_normals[3] ), 3 );
if( warpage > 0 ) return (double)VERDICT_MIN( warpage, VERDICT_DBL_MAX );
return (double)VERDICT_MAX( warpage, -VERDICT_DBL_MAX );
}
C_FUNC_DEF void v_set_quad_size | ( | double | size | ) |
return the average area of a quad
Sets average size (area) of quad, needed for v_quad_relative_size(...)
Definition at line 56 of file V_QuadMetric.cpp.
References size, and verdict_quad_size.
Referenced by moab::VerdictWrapper::set_size(), v_quad_quality(), and v_quad_relative_size_squared().
{
verdict_quad_size = size;
}
double verdict_quad_size = 0 |
the average area of a quad
Definition at line 32 of file V_QuadMetric.cpp.
Referenced by get_weight(), and v_set_quad_size().