MOAB: Mesh Oriented datABase
(version 5.3.1)
|
00001 /* ***************************************************************** 00002 MESQUITE -- The Mesh Quality Improvement Toolkit 00003 00004 Copyright 2004 Sandia Corporation and Argonne National 00005 Laboratory. Under the terms of Contract DE-AC04-94AL85000 00006 with Sandia Corporation, the U.S. Government retains certain 00007 rights in this software. 00008 00009 This library is free software; you can redistribute it and/or 00010 modify it under the terms of the GNU Lesser General Public 00011 License as published by the Free Software Foundation; either 00012 version 2.1 of the License, or (at your option) any later version. 00013 00014 This library is distributed in the hope that it will be useful, 00015 but WITHOUT ANY WARRANTY; without even the implied warranty of 00016 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 00017 Lesser General Public License for more details. 00018 00019 You should have received a copy of the GNU Lesser General Public License 00020 (lgpl.txt) along with this library; if not, write to the Free Software 00021 Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA 00022 00023 diachin2@llnl.gov, djmelan@sandia.gov, mbrewer@sandia.gov, 00024 pknupp@sandia.gov, tleurent@mcs.anl.gov, tmunson@mcs.anl.gov 00025 00026 ***************************************************************** */ 00027 /*! 00028 \file IdealWeightMeanRatio.cpp 00029 \brief 00030 00031 \author Michael Brewer 00032 \date 2002-11-11 00033 */ 00034 #include "IdealWeightMeanRatio.hpp" 00035 #include "MeanRatioFunctions.hpp" 00036 #include "Vector3D.hpp" 00037 #include "MsqTimer.hpp" 00038 #include "MsqDebug.hpp" 00039 #include "MsqError.hpp" 00040 #include "PatchData.hpp" 00041 00042 #include <cmath> 00043 #include <vector> 00044 using std::vector; 00045 00046 using namespace MBMesquite; 00047 00048 std::string IdealWeightMeanRatio::get_name() const 00049 { 00050 return "Mean Ratio"; 00051 } 00052 00053 int IdealWeightMeanRatio::get_negate_flag() const 00054 { 00055 return -1; 00056 } 00057 00058 bool IdealWeightMeanRatio::evaluate( PatchData& pd, size_t handle, double& m, MsqError& err ) 00059 { 00060 const MsqMeshEntity* e = &pd.element_by_index( handle ); 00061 EntityTopology topo = e->get_element_type(); 00062 00063 const MsqVertex* vertices = pd.get_vertex_array( err ); 00064 const size_t* v_i = e->get_vertex_index_array(); 00065 00066 Vector3D n; // Surface normal for 2D objects 00067 00068 // Prism and Hex element descriptions 00069 static const int locs_pri[6][4] = { { 0, 1, 2, 3 }, { 1, 2, 0, 4 }, { 2, 0, 1, 5 }, 00070 { 3, 5, 4, 0 }, { 4, 3, 5, 1 }, { 5, 4, 3, 2 } }; 00071 static const int locs_hex[8][4] = { { 0, 1, 3, 4 }, { 1, 2, 0, 5 }, { 2, 3, 1, 6 }, { 3, 0, 2, 7 }, 00072 { 4, 7, 5, 0 }, { 5, 4, 6, 1 }, { 6, 5, 7, 2 }, { 7, 6, 4, 3 } }; 00073 00074 const Vector3D d_con( 1.0, 1.0, 1.0 ); 00075 00076 int i; 00077 00078 m = 0.0; 00079 bool metric_valid = false; 00080 switch( topo ) 00081 { 00082 case TRIANGLE: 00083 pd.get_domain_normal_at_element( e, n, err ); 00084 MSQ_ERRZERO( err ); 00085 n = n / n.length(); // Need unit normal 00086 mCoords[0] = vertices[v_i[0]]; 00087 mCoords[1] = vertices[v_i[1]]; 00088 mCoords[2] = vertices[v_i[2]]; 00089 metric_valid = m_fcn_2e( m, mCoords, n, a2Con, b2Con, c2Con ); 00090 if( !metric_valid ) return false; 00091 break; 00092 00093 case QUADRILATERAL: 00094 pd.get_domain_normal_at_element( e, n, err ); 00095 MSQ_ERRZERO( err ); 00096 for( i = 0; i < 4; ++i ) 00097 { 00098 n = n / n.length(); // Need unit normal 00099 mCoords[0] = vertices[v_i[locs_hex[i][0]]]; 00100 mCoords[1] = vertices[v_i[locs_hex[i][1]]]; 00101 mCoords[2] = vertices[v_i[locs_hex[i][2]]]; 00102 metric_valid = m_fcn_2i( mMetrics[i], mCoords, n, a2Con, b2Con, c2Con, d_con ); 00103 if( !metric_valid ) return false; 00104 } 00105 m = average_metrics( mMetrics, 4, err ); 00106 MSQ_ERRZERO( err ); 00107 break; 00108 00109 case TETRAHEDRON: 00110 mCoords[0] = vertices[v_i[0]]; 00111 mCoords[1] = vertices[v_i[1]]; 00112 mCoords[2] = vertices[v_i[2]]; 00113 mCoords[3] = vertices[v_i[3]]; 00114 metric_valid = m_fcn_3e( m, mCoords, a3Con, b3Con, c3Con ); 00115 if( !metric_valid ) return false; 00116 break; 00117 00118 case PYRAMID: 00119 for( i = 0; i < 4; ++i ) 00120 { 00121 mCoords[0] = vertices[v_i[i]]; 00122 mCoords[1] = vertices[v_i[( i + 1 ) % 4]]; 00123 mCoords[2] = vertices[v_i[( i + 3 ) % 4]]; 00124 mCoords[3] = vertices[v_i[4]]; 00125 metric_valid = m_fcn_3p( mMetrics[i], mCoords, a3Con, b3Con, c3Con ); 00126 if( !metric_valid ) return false; 00127 } 00128 m = average_metrics( mMetrics, 4, err ); 00129 MSQ_ERRZERO( err ); 00130 break; 00131 00132 case PRISM: 00133 for( i = 0; i < 6; ++i ) 00134 { 00135 mCoords[0] = vertices[v_i[locs_pri[i][0]]]; 00136 mCoords[1] = vertices[v_i[locs_pri[i][1]]]; 00137 mCoords[2] = vertices[v_i[locs_pri[i][2]]]; 00138 mCoords[3] = vertices[v_i[locs_pri[i][3]]]; 00139 metric_valid = m_fcn_3w( mMetrics[i], mCoords, a3Con, b3Con, c3Con ); 00140 if( !metric_valid ) return false; 00141 } 00142 m = average_metrics( mMetrics, 6, err ); 00143 MSQ_ERRZERO( err ); 00144 break; 00145 00146 case HEXAHEDRON: 00147 for( i = 0; i < 8; ++i ) 00148 { 00149 mCoords[0] = vertices[v_i[locs_hex[i][0]]]; 00150 mCoords[1] = vertices[v_i[locs_hex[i][1]]]; 00151 mCoords[2] = vertices[v_i[locs_hex[i][2]]]; 00152 mCoords[3] = vertices[v_i[locs_hex[i][3]]]; 00153 metric_valid = m_fcn_3i( mMetrics[i], mCoords, a3Con, b3Con, c3Con, d_con ); 00154 if( !metric_valid ) return false; 00155 } 00156 m = average_metrics( mMetrics, 8, err ); 00157 MSQ_ERRZERO( err ); 00158 break; 00159 00160 default: 00161 MSQ_SETERR( err ) 00162 ( MsqError::UNSUPPORTED_ELEMENT, "Element type (%d) not supported in IdealWeightMeanRatio", (int)topo ); 00163 return false; 00164 } // end switch over element type 00165 return true; 00166 } 00167 00168 bool IdealWeightMeanRatio::evaluate_with_gradient( PatchData& pd, size_t handle, double& m, 00169 std::vector< size_t >& indices, std::vector< Vector3D >& g, 00170 MsqError& err ) 00171 { 00172 // FUNCTION_TIMER_START(__FUNC__); 00173 const MsqMeshEntity* e = &pd.element_by_index( handle ); 00174 EntityTopology topo = e->get_element_type(); 00175 00176 if( !analytical_average_gradient() && topo != TRIANGLE && topo != TETRAHEDRON ) 00177 { 00178 static bool print = true; 00179 if( print ) 00180 { 00181 MSQ_DBGOUT( 1 ) << "Analyical gradient not available for selected averaging scheme. " 00182 << "Using (possibly much slower) numerical approximation of gradient" 00183 << " of quality metric. " << std::endl; 00184 print = false; 00185 } 00186 return QualityMetric::evaluate_with_gradient( pd, handle, m, indices, g, err ); 00187 } 00188 00189 const MsqVertex* vertices = pd.get_vertex_array( err ); 00190 const size_t* v_i = e->get_vertex_index_array(); 00191 00192 Vector3D n; // Surface normal for 2D objects 00193 00194 // double nm, t=0; 00195 00196 // Prism and Hex element descriptions 00197 static const int locs_pri[6][4] = { { 0, 1, 2, 3 }, { 1, 2, 0, 4 }, { 2, 0, 1, 5 }, 00198 { 3, 5, 4, 0 }, { 4, 3, 5, 1 }, { 5, 4, 3, 2 } }; 00199 static const int locs_hex[8][4] = { { 0, 1, 3, 4 }, { 1, 2, 0, 5 }, { 2, 3, 1, 6 }, { 3, 0, 2, 7 }, 00200 { 4, 7, 5, 0 }, { 5, 4, 6, 1 }, { 6, 5, 7, 2 }, { 7, 6, 4, 3 } }; 00201 00202 const Vector3D d_con( 1.0, 1.0, 1.0 ); 00203 00204 int i; 00205 00206 bool metric_valid = false; 00207 const uint32_t fm = fixed_vertex_bitmap( pd, e, indices ); 00208 00209 m = 0.0; 00210 00211 switch( topo ) 00212 { 00213 case TRIANGLE: 00214 pd.get_domain_normal_at_element( e, n, err ); 00215 MSQ_ERRZERO( err ); 00216 n = n / n.length(); // Need unit normal 00217 mCoords[0] = vertices[v_i[0]]; 00218 mCoords[1] = vertices[v_i[1]]; 00219 mCoords[2] = vertices[v_i[2]]; 00220 g.resize( 3 ); 00221 if( !g_fcn_2e( m, arrptr( g ), mCoords, n, a2Con, b2Con, c2Con ) ) return false; 00222 break; 00223 00224 case QUADRILATERAL: 00225 pd.get_domain_normal_at_element( e, n, err ); 00226 MSQ_ERRZERO( err ); 00227 n /= n.length(); // Need unit normal 00228 for( i = 0; i < 4; ++i ) 00229 { 00230 mCoords[0] = vertices[v_i[locs_hex[i][0]]]; 00231 mCoords[1] = vertices[v_i[locs_hex[i][1]]]; 00232 mCoords[2] = vertices[v_i[locs_hex[i][2]]]; 00233 if( !g_fcn_2i( mMetrics[i], mGradients + 3 * i, mCoords, n, a2Con, b2Con, c2Con, d_con ) ) return false; 00234 } 00235 00236 g.resize( 4 ); 00237 m = average_corner_gradients( QUADRILATERAL, fm, 4, mMetrics, mGradients, arrptr( g ), err ); 00238 MSQ_ERRZERO( err ); 00239 break; 00240 00241 case TETRAHEDRON: 00242 mCoords[0] = vertices[v_i[0]]; 00243 mCoords[1] = vertices[v_i[1]]; 00244 mCoords[2] = vertices[v_i[2]]; 00245 mCoords[3] = vertices[v_i[3]]; 00246 g.resize( 4 ); 00247 metric_valid = g_fcn_3e( m, arrptr( g ), mCoords, a3Con, b3Con, c3Con ); 00248 if( !metric_valid ) return false; 00249 break; 00250 00251 case PYRAMID: 00252 for( i = 0; i < 4; ++i ) 00253 { 00254 mCoords[0] = vertices[v_i[i]]; 00255 mCoords[1] = vertices[v_i[( i + 1 ) % 4]]; 00256 mCoords[2] = vertices[v_i[( i + 3 ) % 4]]; 00257 mCoords[3] = vertices[v_i[4]]; 00258 metric_valid = g_fcn_3p( mMetrics[i], mGradients + 4 * i, mCoords, a3Con, b3Con, c3Con ); 00259 if( !metric_valid ) return false; 00260 } 00261 g.resize( 5 ); 00262 m = average_corner_gradients( PYRAMID, fm, 4, mMetrics, mGradients, arrptr( g ), err ); 00263 MSQ_ERRZERO( err ); 00264 break; 00265 00266 case PRISM: 00267 for( i = 0; i < 6; ++i ) 00268 { 00269 mCoords[0] = vertices[v_i[locs_pri[i][0]]]; 00270 mCoords[1] = vertices[v_i[locs_pri[i][1]]]; 00271 mCoords[2] = vertices[v_i[locs_pri[i][2]]]; 00272 mCoords[3] = vertices[v_i[locs_pri[i][3]]]; 00273 if( !g_fcn_3w( mMetrics[i], mGradients + 4 * i, mCoords, a3Con, b3Con, c3Con ) ) return false; 00274 } 00275 g.resize( 6 ); 00276 m = average_corner_gradients( PRISM, fm, 6, mMetrics, mGradients, arrptr( g ), err ); 00277 MSQ_ERRZERO( err ); 00278 break; 00279 00280 case HEXAHEDRON: 00281 for( i = 0; i < 8; ++i ) 00282 { 00283 mCoords[0] = vertices[v_i[locs_hex[i][0]]]; 00284 mCoords[1] = vertices[v_i[locs_hex[i][1]]]; 00285 mCoords[2] = vertices[v_i[locs_hex[i][2]]]; 00286 mCoords[3] = vertices[v_i[locs_hex[i][3]]]; 00287 if( !g_fcn_3i( mMetrics[i], mGradients + 4 * i, mCoords, a3Con, b3Con, c3Con, d_con ) ) return false; 00288 } 00289 g.resize( 8 ); 00290 m = average_corner_gradients( HEXAHEDRON, fm, 8, mMetrics, mGradients, arrptr( g ), err ); 00291 MSQ_ERRZERO( err ); 00292 break; 00293 00294 default: 00295 MSQ_SETERR( err ) 00296 ( MsqError::UNSUPPORTED_ELEMENT, "Element type (%d) not supported in IdealWeightMeanRatio", (int)topo ); 00297 return false; 00298 } 00299 00300 remove_fixed_gradients( topo, fm, g ); 00301 return true; 00302 } 00303 00304 bool IdealWeightMeanRatio::evaluate_with_Hessian_diagonal( PatchData& pd, size_t handle, double& m, 00305 std::vector< size_t >& indices, std::vector< Vector3D >& g, 00306 std::vector< SymMatrix3D >& h, MsqError& err ) 00307 { 00308 const MsqMeshEntity* e = &pd.element_by_index( handle ); 00309 EntityTopology topo = e->get_element_type(); 00310 00311 if( !analytical_average_hessian() && topo != TRIANGLE && topo != TETRAHEDRON ) 00312 { 00313 static bool print = true; 00314 if( print ) 00315 { 00316 MSQ_DBGOUT( 1 ) << "Analyical gradient not available for selected averaging scheme. " 00317 << "Using (possibly much slower) numerical approximation of gradient" 00318 << " of quality metric. " << std::endl; 00319 print = false; 00320 } 00321 return QualityMetric::evaluate_with_Hessian_diagonal( pd, handle, m, indices, g, h, err ); 00322 } 00323 00324 const MsqVertex* vertices = pd.get_vertex_array( err ); 00325 const size_t* v_i = e->get_vertex_index_array(); 00326 00327 Vector3D n; // Surface normal for 2D objects 00328 00329 // Prism and Hex element descriptions 00330 static const int locs_pri[6][4] = { { 0, 1, 2, 3 }, { 1, 2, 0, 4 }, { 2, 0, 1, 5 }, 00331 { 3, 5, 4, 0 }, { 4, 3, 5, 1 }, { 5, 4, 3, 2 } }; 00332 static const int locs_hex[8][4] = { { 0, 1, 3, 4 }, { 1, 2, 0, 5 }, { 2, 3, 1, 6 }, { 3, 0, 2, 7 }, 00333 { 4, 7, 5, 0 }, { 5, 4, 6, 1 }, { 6, 5, 7, 2 }, { 7, 6, 4, 3 } }; 00334 00335 const Vector3D d_con( 1.0, 1.0, 1.0 ); 00336 00337 int i; 00338 00339 bool metric_valid = false; 00340 const uint32_t fm = fixed_vertex_bitmap( pd, e, indices ); 00341 00342 m = 0.0; 00343 00344 switch( topo ) 00345 { 00346 case TRIANGLE: 00347 pd.get_domain_normal_at_element( e, n, err ); 00348 MSQ_ERRZERO( err ); 00349 n = n / n.length(); // Need unit normal 00350 mCoords[0] = vertices[v_i[0]]; 00351 mCoords[1] = vertices[v_i[1]]; 00352 mCoords[2] = vertices[v_i[2]]; 00353 g.resize( 3 ), h.resize( 3 ); 00354 if( !h_fcn_2e( m, arrptr( g ), mHessians, mCoords, n, a2Con, b2Con, c2Con ) ) return false; 00355 h[0] = mHessians[0].upper(); 00356 h[1] = mHessians[3].upper(); 00357 h[2] = mHessians[5].upper(); 00358 break; 00359 00360 case QUADRILATERAL: 00361 pd.get_domain_normal_at_element( e, n, err ); 00362 MSQ_ERRZERO( err ); 00363 n = n / n.length(); // Need unit normal 00364 for( i = 0; i < 4; ++i ) 00365 { 00366 mCoords[0] = vertices[v_i[locs_hex[i][0]]]; 00367 mCoords[1] = vertices[v_i[locs_hex[i][1]]]; 00368 mCoords[2] = vertices[v_i[locs_hex[i][2]]]; 00369 if( !h_fcn_2i( mMetrics[i], mGradients + 3 * i, mHessians + 6 * i, mCoords, n, a2Con, b2Con, c2Con, 00370 d_con ) ) 00371 return false; 00372 } 00373 00374 g.resize( 4 ), h.resize( 4 ); 00375 m = average_corner_hessian_diagonals( QUADRILATERAL, fm, 4, mMetrics, mGradients, mHessians, arrptr( g ), 00376 arrptr( h ), err ); 00377 MSQ_ERRZERO( err ); 00378 break; 00379 00380 case TETRAHEDRON: 00381 mCoords[0] = vertices[v_i[0]]; 00382 mCoords[1] = vertices[v_i[1]]; 00383 mCoords[2] = vertices[v_i[2]]; 00384 mCoords[3] = vertices[v_i[3]]; 00385 g.resize( 4 ), h.resize( 4 ); 00386 metric_valid = h_fcn_3e( m, arrptr( g ), mHessians, mCoords, a3Con, b3Con, c3Con ); 00387 if( !metric_valid ) return false; 00388 h[0] = mHessians[0].upper(); 00389 h[1] = mHessians[4].upper(); 00390 h[2] = mHessians[7].upper(); 00391 h[3] = mHessians[9].upper(); 00392 break; 00393 00394 case PYRAMID: 00395 for( i = 0; i < 4; ++i ) 00396 { 00397 mCoords[0] = vertices[v_i[i]]; 00398 mCoords[1] = vertices[v_i[( i + 1 ) % 4]]; 00399 mCoords[2] = vertices[v_i[( i + 3 ) % 4]]; 00400 mCoords[3] = vertices[v_i[4]]; 00401 metric_valid = 00402 h_fcn_3p( mMetrics[i], mGradients + 4 * i, mHessians + 10 * i, mCoords, a3Con, b3Con, c3Con ); 00403 if( !metric_valid ) return false; 00404 } 00405 00406 g.resize( 5 ), h.resize( 5 ); 00407 m = average_corner_hessian_diagonals( PYRAMID, fm, 4, mMetrics, mGradients, mHessians, arrptr( g ), 00408 arrptr( h ), err ); 00409 MSQ_ERRZERO( err ); 00410 break; 00411 00412 case PRISM: 00413 for( i = 0; i < 6; ++i ) 00414 { 00415 mCoords[0] = vertices[v_i[locs_pri[i][0]]]; 00416 mCoords[1] = vertices[v_i[locs_pri[i][1]]]; 00417 mCoords[2] = vertices[v_i[locs_pri[i][2]]]; 00418 mCoords[3] = vertices[v_i[locs_pri[i][3]]]; 00419 if( !h_fcn_3w( mMetrics[i], mGradients + 4 * i, mHessians + 10 * i, mCoords, a3Con, b3Con, c3Con ) ) 00420 return false; 00421 } 00422 00423 g.resize( 6 ), h.resize( 6 ); 00424 m = average_corner_hessian_diagonals( PRISM, fm, 6, mMetrics, mGradients, mHessians, arrptr( g ), 00425 arrptr( h ), err ); 00426 MSQ_ERRZERO( err ); 00427 break; 00428 00429 case HEXAHEDRON: 00430 for( i = 0; i < 8; ++i ) 00431 { 00432 mCoords[0] = vertices[v_i[locs_hex[i][0]]]; 00433 mCoords[1] = vertices[v_i[locs_hex[i][1]]]; 00434 mCoords[2] = vertices[v_i[locs_hex[i][2]]]; 00435 mCoords[3] = vertices[v_i[locs_hex[i][3]]]; 00436 if( !h_fcn_3i( mMetrics[i], mGradients + 4 * i, mHessians + 10 * i, mCoords, a3Con, b3Con, c3Con, 00437 d_con ) ) 00438 return false; 00439 } 00440 00441 g.resize( 8 ), h.resize( 8 ); 00442 m = average_corner_hessian_diagonals( HEXAHEDRON, fm, 8, mMetrics, mGradients, mHessians, arrptr( g ), 00443 arrptr( h ), err ); 00444 MSQ_ERRZERO( err ); 00445 break; 00446 00447 default: 00448 MSQ_SETERR( err ) 00449 ( MsqError::UNSUPPORTED_ELEMENT, "Element type (%d) not supported in IdealWeightMeanRatio", (int)topo ); 00450 return false; 00451 } // end switch over element type 00452 00453 remove_fixed_diagonals( topo, fm, g, h ); 00454 return true; 00455 } 00456 00457 bool IdealWeightMeanRatio::evaluate_with_Hessian( PatchData& pd, size_t handle, double& m, 00458 std::vector< size_t >& indices, std::vector< Vector3D >& g, 00459 std::vector< Matrix3D >& h, MsqError& err ) 00460 { 00461 // FUNCTION_TIMER_START(__FUNC__); 00462 const MsqMeshEntity* e = &pd.element_by_index( handle ); 00463 EntityTopology topo = e->get_element_type(); 00464 00465 if( !analytical_average_hessian() && topo != TRIANGLE && topo != TETRAHEDRON ) 00466 { 00467 static bool print = true; 00468 if( print ) 00469 { 00470 MSQ_DBGOUT( 1 ) << "Analyical gradient not available for selected averaging scheme. " 00471 << "Using (possibly much slower) numerical approximation of gradient" 00472 << " of quality metric. " << std::endl; 00473 print = false; 00474 } 00475 return QualityMetric::evaluate_with_Hessian( pd, handle, m, indices, g, h, err ); 00476 } 00477 00478 const MsqVertex* vertices = pd.get_vertex_array( err ); 00479 const size_t* v_i = e->get_vertex_index_array(); 00480 00481 Vector3D n; // Surface normal for 2D objects 00482 00483 // Prism and Hex element descriptions 00484 static const int locs_pri[6][4] = { { 0, 1, 2, 3 }, { 1, 2, 0, 4 }, { 2, 0, 1, 5 }, 00485 { 3, 5, 4, 0 }, { 4, 3, 5, 1 }, { 5, 4, 3, 2 } }; 00486 static const int locs_hex[8][4] = { { 0, 1, 3, 4 }, { 1, 2, 0, 5 }, { 2, 3, 1, 6 }, { 3, 0, 2, 7 }, 00487 { 4, 7, 5, 0 }, { 5, 4, 6, 1 }, { 6, 5, 7, 2 }, { 7, 6, 4, 3 } }; 00488 00489 const Vector3D d_con( 1.0, 1.0, 1.0 ); 00490 00491 int i; 00492 00493 bool metric_valid = false; 00494 const uint32_t fm = fixed_vertex_bitmap( pd, e, indices ); 00495 00496 m = 0.0; 00497 00498 switch( topo ) 00499 { 00500 case TRIANGLE: 00501 pd.get_domain_normal_at_element( e, n, err ); 00502 MSQ_ERRZERO( err ); 00503 n = n / n.length(); // Need unit normal 00504 mCoords[0] = vertices[v_i[0]]; 00505 mCoords[1] = vertices[v_i[1]]; 00506 mCoords[2] = vertices[v_i[2]]; 00507 g.resize( 3 ), h.resize( 6 ); 00508 if( !h_fcn_2e( m, arrptr( g ), arrptr( h ), mCoords, n, a2Con, b2Con, c2Con ) ) return false; 00509 break; 00510 00511 case QUADRILATERAL: 00512 pd.get_domain_normal_at_element( e, n, err ); 00513 MSQ_ERRZERO( err ); 00514 n = n / n.length(); // Need unit normal 00515 for( i = 0; i < 4; ++i ) 00516 { 00517 mCoords[0] = vertices[v_i[locs_hex[i][0]]]; 00518 mCoords[1] = vertices[v_i[locs_hex[i][1]]]; 00519 mCoords[2] = vertices[v_i[locs_hex[i][2]]]; 00520 if( !h_fcn_2i( mMetrics[i], mGradients + 3 * i, mHessians + 6 * i, mCoords, n, a2Con, b2Con, c2Con, 00521 d_con ) ) 00522 return false; 00523 } 00524 00525 g.resize( 4 ), h.resize( 10 ); 00526 m = average_corner_hessians( QUADRILATERAL, fm, 4, mMetrics, mGradients, mHessians, arrptr( g ), 00527 arrptr( h ), err ); 00528 MSQ_ERRZERO( err ); 00529 break; 00530 00531 case TETRAHEDRON: 00532 mCoords[0] = vertices[v_i[0]]; 00533 mCoords[1] = vertices[v_i[1]]; 00534 mCoords[2] = vertices[v_i[2]]; 00535 mCoords[3] = vertices[v_i[3]]; 00536 g.resize( 4 ), h.resize( 10 ); 00537 metric_valid = h_fcn_3e( m, arrptr( g ), arrptr( h ), mCoords, a3Con, b3Con, c3Con ); 00538 if( !metric_valid ) return false; 00539 break; 00540 00541 case PYRAMID: 00542 for( i = 0; i < 4; ++i ) 00543 { 00544 mCoords[0] = vertices[v_i[i]]; 00545 mCoords[1] = vertices[v_i[( i + 1 ) % 4]]; 00546 mCoords[2] = vertices[v_i[( i + 3 ) % 4]]; 00547 mCoords[3] = vertices[v_i[4]]; 00548 metric_valid = 00549 h_fcn_3p( mMetrics[i], mGradients + 4 * i, mHessians + 10 * i, mCoords, a3Con, b3Con, c3Con ); 00550 if( !metric_valid ) return false; 00551 } 00552 00553 g.resize( 5 ), h.resize( 15 ); 00554 m = average_corner_hessians( PYRAMID, fm, 4, mMetrics, mGradients, mHessians, arrptr( g ), arrptr( h ), 00555 err ); 00556 MSQ_ERRZERO( err ); 00557 break; 00558 00559 case PRISM: 00560 for( i = 0; i < 6; ++i ) 00561 { 00562 mCoords[0] = vertices[v_i[locs_pri[i][0]]]; 00563 mCoords[1] = vertices[v_i[locs_pri[i][1]]]; 00564 mCoords[2] = vertices[v_i[locs_pri[i][2]]]; 00565 mCoords[3] = vertices[v_i[locs_pri[i][3]]]; 00566 if( !h_fcn_3w( mMetrics[i], mGradients + 4 * i, mHessians + 10 * i, mCoords, a3Con, b3Con, c3Con ) ) 00567 return false; 00568 } 00569 00570 g.resize( 6 ), h.resize( 21 ); 00571 m = average_corner_hessians( PRISM, fm, 6, mMetrics, mGradients, mHessians, arrptr( g ), arrptr( h ), err ); 00572 MSQ_ERRZERO( err ); 00573 break; 00574 00575 case HEXAHEDRON: 00576 for( i = 0; i < 8; ++i ) 00577 { 00578 mCoords[0] = vertices[v_i[locs_hex[i][0]]]; 00579 mCoords[1] = vertices[v_i[locs_hex[i][1]]]; 00580 mCoords[2] = vertices[v_i[locs_hex[i][2]]]; 00581 mCoords[3] = vertices[v_i[locs_hex[i][3]]]; 00582 if( !h_fcn_3i( mMetrics[i], mGradients + 4 * i, mHessians + 10 * i, mCoords, a3Con, b3Con, c3Con, 00583 d_con ) ) 00584 return false; 00585 } 00586 00587 g.resize( 8 ), h.resize( 36 ); 00588 m = average_corner_hessians( HEXAHEDRON, fm, 8, mMetrics, mGradients, mHessians, arrptr( g ), arrptr( h ), 00589 err ); 00590 MSQ_ERRZERO( err ); 00591 break; 00592 00593 default: 00594 MSQ_SETERR( err ) 00595 ( MsqError::UNSUPPORTED_ELEMENT, "Element type (%d) not supported in IdealWeightMeanRatio", (int)topo ); 00596 return false; 00597 } // end switch over element type 00598 00599 remove_fixed_gradients( topo, fm, g ); 00600 remove_fixed_hessians( topo, fm, h ); 00601 // FUNCTION_TIMER_END(); 00602 return true; 00603 }