Branch data Line data Source code
1 : : /* *****************************************************************
2 : : MESQUITE -- The Mesh Quality Improvement Toolkit
3 : :
4 : : Copyright 2004 Sandia Corporation and Argonne National
5 : : Laboratory. Under the terms of Contract DE-AC04-94AL85000
6 : : with Sandia Corporation, the U.S. Government retains certain
7 : : rights in this software.
8 : :
9 : : This library is free software; you can redistribute it and/or
10 : : modify it under the terms of the GNU Lesser General Public
11 : : License as published by the Free Software Foundation; either
12 : : version 2.1 of the License, or (at your option) any later version.
13 : :
14 : : This library is distributed in the hope that it will be useful,
15 : : but WITHOUT ANY WARRANTY; without even the implied warranty of
16 : : MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 : : Lesser General Public License for more details.
18 : :
19 : : You should have received a copy of the GNU Lesser General Public License
20 : : (lgpl.txt) along with this library; if not, write to the Free Software
21 : : Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
22 : :
23 : : [email protected], [email protected], [email protected],
24 : : [email protected], [email protected], [email protected]
25 : :
26 : : ***************************************************************** */
27 : : /*! \file UntangleBetaQualityMetric.cpp
28 : : UntangleBeta is an untangle quality metric which can be used to evaluate
29 : : the quality of two- or three-dimensional elements.
30 : :
31 : : \author Michael Brewer
32 : : \date 2002-10-10
33 : : */
34 : :
35 : : #include "UntangleBetaQualityMetric.hpp"
36 : : #include "Vector3D.hpp"
37 : : #include "PatchData.hpp"
38 : : #include "MsqError.hpp"
39 : :
40 : : using namespace MBMesquite;
41 : :
42 : 135216 : static inline void untangle_function_2d( double beta, const Vector3D temp_vec[], size_t e_ind, PatchData& pd,
43 : : double& fval, MsqError& err )
44 : : {
45 [ + - ]: 135216 : Vector3D surface_normal;
46 [ + - ][ + - ]: 270432 : pd.get_domain_normal_at_element( e_ind, surface_normal, err );MSQ_ERRRTN( err );
[ - + ][ # # ]
[ # # ][ - + ]
47 [ + - ]: 135216 : Vector3D cross_vec = temp_vec[0] * temp_vec[1];
48 : : // cout<<"\nsurface_normal "<<surface_normal;
49 : : // cout<<"\cross_vec "<<cross_vec;
50 [ + - ]: 135216 : double temp_var = cross_vec.length();
51 [ + - ][ + + ]: 135216 : if( cross_vec % surface_normal < 0.0 ) { temp_var *= -1; }
52 : 135216 : temp_var -= beta;
53 : : // cout<<"temp_var == "<<temp_var;
54 : 135216 : fval = 0.0;
55 [ + + ]: 135216 : if( temp_var < 0.0 ) { fval = fabs( temp_var ) - temp_var; }
56 : : // cout<<"\nfval == "<<fval<<" e_ind "<<e_ind;
57 : : }
58 : :
59 : 8079689 : static inline void untangle_function_3d( double beta, const Vector3D temp_vec[], double& fval )
60 : : {
61 [ + - ]: 8079689 : double temp_var = temp_vec[0] % ( temp_vec[1] * temp_vec[2] );
62 : 8079689 : temp_var -= beta;
63 : 8079689 : fval = 0.0;
64 [ + + ]: 8079689 : if( temp_var < 0.0 ) { fval = fabs( temp_var ) - temp_var; }
65 : 8079689 : }
66 : :
67 : : /*! \fn UntangleBetaQualityMetric::UntangleBetaQualityMetric(double bet)
68 : : \brief For untangle beta, the constructor defaults to the SUM
69 : : averaging method, and to the ELEMENT_VERTICES evaluation mode.
70 : : */
71 [ + - ]: 4 : UntangleBetaQualityMetric::UntangleBetaQualityMetric( double bet ) : AveragingQM( RMS ), mBeta( bet ) {}
72 : :
73 : 2 : std::string UntangleBetaQualityMetric::get_name() const
74 : : {
75 [ + - ]: 2 : return "Untangle Beta";
76 : : }
77 : :
78 : 168127 : int UntangleBetaQualityMetric::get_negate_flag() const
79 : : {
80 : 168127 : return 1;
81 : : }
82 : :
83 : 2001808 : bool UntangleBetaQualityMetric::evaluate( PatchData& pd, size_t e_ind, double& fval, MsqError& err )
84 : : {
85 : :
86 : : double met_vals[MSQ_MAX_NUM_VERT_PER_ENT];
87 : 2001808 : fval = MSQ_MAX_CAP;
88 [ + - ]: 2001808 : const MsqMeshEntity* element = &pd.element_by_index( e_ind );
89 [ + - ]: 2001808 : const size_t* v_i = element->get_vertex_index_array();
90 : : // only 3 temp_vec will be sent to untangle calculator, but the
91 : : // additional vector3Ds may be needed during the calculations
92 [ + - ][ + + ]: 12010848 : Vector3D temp_vec[5];
93 [ + - ]: 2001808 : const MsqVertex* vertices = pd.get_vertex_array( err );
94 [ + - ][ - + ]: 2001808 : MSQ_ERRZERO( err );
[ # # ][ # # ]
[ - + ]
95 [ + - ]: 2001808 : EntityTopology type = element->get_element_type();
96 [ - + + + : 2001808 : switch( type )
+ - - ]
97 : : {
98 : : case TRIANGLE:
99 [ # # ][ # # ]: 0 : temp_vec[0] = vertices[v_i[1]] - vertices[v_i[0]];
100 [ # # ][ # # ]: 0 : temp_vec[2] = vertices[v_i[2]] - vertices[v_i[0]];
101 : : // make relative to equilateral
102 [ # # ][ # # ]: 0 : temp_vec[1] = ( ( 2 * temp_vec[2] ) - temp_vec[0] ) * MSQ_SQRT_THREE_INV;
[ # # ][ # # ]
103 [ # # ]: 0 : untangle_function_2d( mBeta, temp_vec, e_ind, pd, fval, err );
104 : 0 : break;
105 : : case QUADRILATERAL:
106 [ + - ][ + - ]: 33804 : temp_vec[0] = vertices[v_i[1]] - vertices[v_i[0]];
107 [ + - ][ + - ]: 33804 : temp_vec[1] = vertices[v_i[3]] - vertices[v_i[0]];
108 [ + - ]: 33804 : untangle_function_2d( mBeta, temp_vec, e_ind, pd, met_vals[0], err );
109 [ + - ][ - + ]: 33804 : MSQ_ERRZERO( err );
[ # # ][ # # ]
[ - + ]
110 : :
111 [ + - ][ + - ]: 33804 : temp_vec[0] = vertices[v_i[2]] - vertices[v_i[1]];
112 [ + - ][ + - ]: 33804 : temp_vec[1] = vertices[v_i[0]] - vertices[v_i[1]];
113 [ + - ]: 33804 : untangle_function_2d( mBeta, temp_vec, e_ind, pd, met_vals[1], err );
114 [ + - ][ - + ]: 33804 : MSQ_ERRZERO( err );
[ # # ][ # # ]
[ - + ]
115 : :
116 [ + - ][ + - ]: 33804 : temp_vec[0] = vertices[v_i[3]] - vertices[v_i[2]];
117 [ + - ][ + - ]: 33804 : temp_vec[1] = vertices[v_i[1]] - vertices[v_i[2]];
118 [ + - ]: 33804 : untangle_function_2d( mBeta, temp_vec, e_ind, pd, met_vals[2], err );
119 [ + - ][ - + ]: 33804 : MSQ_ERRZERO( err );
[ # # ][ # # ]
[ - + ]
120 : :
121 [ + - ][ + - ]: 33804 : temp_vec[0] = vertices[v_i[0]] - vertices[v_i[3]];
122 [ + - ][ + - ]: 33804 : temp_vec[1] = vertices[v_i[2]] - vertices[v_i[3]];
123 [ + - ]: 33804 : untangle_function_2d( mBeta, temp_vec, e_ind, pd, met_vals[3], err );
124 [ + - ][ - + ]: 33804 : MSQ_ERRZERO( err );
[ # # ][ # # ]
[ - + ]
125 [ + - ]: 33804 : fval = average_metrics( met_vals, 4, err );
126 [ + - ][ - + ]: 33804 : MSQ_ERRZERO( err );
[ # # ][ # # ]
[ - + ]
127 : 33804 : break;
128 : : case TETRAHEDRON:
129 [ + - ][ + - ]: 871889 : temp_vec[0] = vertices[v_i[1]] - vertices[v_i[0]];
130 [ + - ][ + - ]: 871889 : temp_vec[3] = vertices[v_i[2]] - vertices[v_i[0]];
131 [ + - ][ + - ]: 871889 : temp_vec[4] = vertices[v_i[3]] - vertices[v_i[0]];
132 : : // transform to equilateral tet
133 [ + - ][ + - ]: 871889 : temp_vec[1] = ( ( 2 * temp_vec[3] ) - temp_vec[0] ) / MSQ_SQRT_THREE;
[ + - ][ + - ]
134 [ + - ][ + - ]: 871889 : temp_vec[2] = ( ( 3 * temp_vec[4] ) - temp_vec[0] - temp_vec[3] ) / ( MSQ_SQRT_THREE * MSQ_SQRT_TWO );
[ + - ][ + - ]
[ + - ]
135 [ + - ]: 871889 : untangle_function_3d( mBeta, temp_vec, fval );
136 : 871889 : break;
137 : : case HEXAHEDRON:
138 : : // transform to v_i[0]
139 [ + - ][ + - ]: 705835 : temp_vec[0] = vertices[v_i[1]] - vertices[v_i[0]];
140 [ + - ][ + - ]: 705835 : temp_vec[1] = vertices[v_i[3]] - vertices[v_i[0]];
141 [ + - ][ + - ]: 705835 : temp_vec[2] = vertices[v_i[4]] - vertices[v_i[0]];
142 [ + - ]: 705835 : untangle_function_3d( mBeta, temp_vec, met_vals[0] );
143 : :
144 [ + - ][ + - ]: 705835 : temp_vec[0] = vertices[v_i[2]] - vertices[v_i[1]];
145 [ + - ][ + - ]: 705835 : temp_vec[1] = vertices[v_i[0]] - vertices[v_i[1]];
146 [ + - ][ + - ]: 705835 : temp_vec[2] = vertices[v_i[5]] - vertices[v_i[1]];
147 [ + - ]: 705835 : untangle_function_3d( mBeta, temp_vec, met_vals[1] );
148 : :
149 [ + - ][ + - ]: 705835 : temp_vec[0] = vertices[v_i[3]] - vertices[v_i[2]];
150 [ + - ][ + - ]: 705835 : temp_vec[1] = vertices[v_i[1]] - vertices[v_i[2]];
151 [ + - ][ + - ]: 705835 : temp_vec[2] = vertices[v_i[6]] - vertices[v_i[2]];
152 [ + - ]: 705835 : untangle_function_3d( mBeta, temp_vec, met_vals[2] );
153 : :
154 [ + - ][ + - ]: 705835 : temp_vec[0] = vertices[v_i[0]] - vertices[v_i[3]];
155 [ + - ][ + - ]: 705835 : temp_vec[1] = vertices[v_i[2]] - vertices[v_i[3]];
156 [ + - ][ + - ]: 705835 : temp_vec[2] = vertices[v_i[7]] - vertices[v_i[3]];
157 [ + - ]: 705835 : untangle_function_3d( mBeta, temp_vec, met_vals[3] );
158 : :
159 [ + - ][ + - ]: 705835 : temp_vec[0] = vertices[v_i[7]] - vertices[v_i[4]];
160 [ + - ][ + - ]: 705835 : temp_vec[1] = vertices[v_i[5]] - vertices[v_i[4]];
161 [ + - ][ + - ]: 705835 : temp_vec[2] = vertices[v_i[0]] - vertices[v_i[4]];
162 [ + - ]: 705835 : untangle_function_3d( mBeta, temp_vec, met_vals[4] );
163 : :
164 [ + - ][ + - ]: 705835 : temp_vec[0] = vertices[v_i[4]] - vertices[v_i[5]];
165 [ + - ][ + - ]: 705835 : temp_vec[1] = vertices[v_i[6]] - vertices[v_i[5]];
166 [ + - ][ + - ]: 705835 : temp_vec[2] = vertices[v_i[1]] - vertices[v_i[5]];
167 [ + - ]: 705835 : untangle_function_3d( mBeta, temp_vec, met_vals[5] );
168 : :
169 [ + - ][ + - ]: 705835 : temp_vec[0] = vertices[v_i[5]] - vertices[v_i[6]];
170 [ + - ][ + - ]: 705835 : temp_vec[1] = vertices[v_i[7]] - vertices[v_i[6]];
171 [ + - ][ + - ]: 705835 : temp_vec[2] = vertices[v_i[2]] - vertices[v_i[6]];
172 [ + - ]: 705835 : untangle_function_3d( mBeta, temp_vec, met_vals[6] );
173 : :
174 [ + - ][ + - ]: 705835 : temp_vec[0] = vertices[v_i[6]] - vertices[v_i[7]];
175 [ + - ][ + - ]: 705835 : temp_vec[1] = vertices[v_i[4]] - vertices[v_i[7]];
176 [ + - ][ + - ]: 705835 : temp_vec[2] = vertices[v_i[3]] - vertices[v_i[7]];
177 [ + - ]: 705835 : untangle_function_3d( mBeta, temp_vec, met_vals[7] );
178 [ + - ]: 705835 : fval = average_metrics( met_vals, 8, err );
179 [ + - ][ - + ]: 705835 : MSQ_ERRZERO( err );
[ # # ][ # # ]
[ - + ]
180 : 705835 : break;
181 : : case PYRAMID:
182 [ + + ]: 1951400 : for( unsigned i = 0; i < 4; ++i )
183 : : {
184 [ + - ][ + - ]: 1561120 : temp_vec[0] = vertices[v_i[( i + 1 ) % 4]] - vertices[v_i[i]];
185 [ + - ][ + - ]: 1561120 : temp_vec[1] = vertices[v_i[( i + 3 ) % 4]] - vertices[v_i[i]];
186 [ + - ][ + - ]: 1561120 : temp_vec[3] = vertices[v_i[4]] - vertices[v_i[i]];
187 [ + - ][ + - ]: 1561120 : temp_vec[2] = ( 4 * temp_vec[3] - temp_vec[0] - temp_vec[1] ) / ( 2.0 - MSQ_SQRT_TWO );
[ + - ][ + - ]
[ + - ]
188 [ + - ]: 1561120 : untangle_function_3d( mBeta, temp_vec, met_vals[i] );
189 : : }
190 [ + - ]: 390280 : fval = average_metrics( met_vals, 4, err );
191 [ + - ][ - + ]: 390280 : MSQ_ERRZERO( err );
[ # # ][ # # ]
[ - + ]
192 : 390280 : break;
193 : : case PRISM:
194 [ # # ][ # # ]: 0 : temp_vec[0] = vertices[v_i[1]] - vertices[v_i[0]];
195 [ # # ][ # # ]: 0 : temp_vec[3] = vertices[v_i[2]] - vertices[v_i[0]];
196 [ # # ][ # # ]: 0 : temp_vec[1] = ( ( 2 * temp_vec[3] ) - temp_vec[0] ) * MSQ_SQRT_THREE_INV;
[ # # ][ # # ]
197 [ # # ][ # # ]: 0 : temp_vec[2] = vertices[v_i[3]] - vertices[v_i[0]];
198 [ # # ]: 0 : untangle_function_3d( mBeta, temp_vec, met_vals[0] );
199 : :
200 [ # # ][ # # ]: 0 : temp_vec[0] = vertices[v_i[2]] - vertices[v_i[1]];
201 [ # # ][ # # ]: 0 : temp_vec[3] = vertices[v_i[0]] - vertices[v_i[1]];
202 [ # # ][ # # ]: 0 : temp_vec[1] = ( ( 2 * temp_vec[3] ) - temp_vec[0] ) * MSQ_SQRT_THREE_INV;
[ # # ][ # # ]
203 [ # # ][ # # ]: 0 : temp_vec[2] = vertices[v_i[4]] - vertices[v_i[1]];
204 [ # # ]: 0 : untangle_function_3d( mBeta, temp_vec, met_vals[1] );
205 : :
206 [ # # ][ # # ]: 0 : temp_vec[0] = vertices[v_i[0]] - vertices[v_i[2]];
207 [ # # ][ # # ]: 0 : temp_vec[3] = vertices[v_i[1]] - vertices[v_i[2]];
208 [ # # ][ # # ]: 0 : temp_vec[1] = ( ( 2 * temp_vec[3] ) - temp_vec[0] ) * MSQ_SQRT_THREE_INV;
[ # # ][ # # ]
209 [ # # ][ # # ]: 0 : temp_vec[2] = vertices[v_i[5]] - vertices[v_i[2]];
210 [ # # ]: 0 : untangle_function_3d( mBeta, temp_vec, met_vals[2] );
211 : :
212 [ # # ][ # # ]: 0 : temp_vec[0] = vertices[v_i[5]] - vertices[v_i[3]];
213 [ # # ][ # # ]: 0 : temp_vec[3] = vertices[v_i[4]] - vertices[v_i[3]];
214 [ # # ][ # # ]: 0 : temp_vec[1] = ( ( 2 * temp_vec[3] ) - temp_vec[0] ) * MSQ_SQRT_THREE_INV;
[ # # ][ # # ]
215 [ # # ][ # # ]: 0 : temp_vec[2] = vertices[v_i[0]] - vertices[v_i[3]];
216 [ # # ]: 0 : untangle_function_3d( mBeta, temp_vec, met_vals[3] );
217 : :
218 [ # # ][ # # ]: 0 : temp_vec[0] = vertices[v_i[3]] - vertices[v_i[4]];
219 [ # # ][ # # ]: 0 : temp_vec[3] = vertices[v_i[5]] - vertices[v_i[4]];
220 [ # # ][ # # ]: 0 : temp_vec[1] = ( ( 2 * temp_vec[3] ) - temp_vec[0] ) * MSQ_SQRT_THREE_INV;
[ # # ][ # # ]
221 [ # # ][ # # ]: 0 : temp_vec[2] = vertices[v_i[1]] - vertices[v_i[4]];
222 [ # # ]: 0 : untangle_function_3d( mBeta, temp_vec, met_vals[4] );
223 : :
224 [ # # ][ # # ]: 0 : temp_vec[0] = vertices[v_i[4]] - vertices[v_i[5]];
225 [ # # ][ # # ]: 0 : temp_vec[3] = vertices[v_i[3]] - vertices[v_i[5]];
226 [ # # ][ # # ]: 0 : temp_vec[1] = ( ( 2 * temp_vec[3] ) - temp_vec[0] ) * MSQ_SQRT_THREE_INV;
[ # # ][ # # ]
227 [ # # ][ # # ]: 0 : temp_vec[2] = vertices[v_i[2]] - vertices[v_i[5]];
228 [ # # ]: 0 : untangle_function_3d( mBeta, temp_vec, met_vals[5] );
229 : :
230 [ # # ]: 0 : fval = average_metrics( met_vals, 6, err );
231 [ # # ][ # # ]: 0 : MSQ_ERRZERO( err );
[ # # ][ # # ]
[ # # ]
232 : 0 : break;
233 : : default:
234 : : MSQ_SETERR( err )
235 [ # # ][ # # ]: 0 : ( MsqError::UNSUPPORTED_ELEMENT, "Unsupported cell type (%d) for Untangle quality metric.", type );
236 : :
237 : 0 : fval = MSQ_MAX_CAP;
238 : 0 : return false;
239 : : } // end switch over element type
240 : 2001808 : return true;
241 [ + - ][ + - ]: 16 : }
|