1/*========================================================================= 2 3 Module: verdict.h.in 4 5 Copyright (c) 2006 Sandia Corporation. 6 All rights reserved. 7 See Copyright.txt or http://www.kitware.com/Copyright.htm for details. 8 9 This software is distributed WITHOUT ANY WARRANTY; without even 10 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR 11 PURPOSE. See the above copyright notice for more information. 12 13=========================================================================*/ 14 15 16/*! \file verdict.h 17 \brief Header file for verdict library that calculates metrics for finite elements. 18 Also see: \ref index "Main Page" 19 * 20 * verdict.h is the header file for applications/libraries to include 21 * to compute quality metrics. 22 * 23 * This file is part of VERDICT 24 * 25 */ 26 27// .SECTION Thanks 28// Prior to its inclusion within VTK, this code was developed by the CUBIT 29// project at Sandia National Laboratories. 30 31#ifndef __verdict_h 32#define __verdict_h 33 34#define VERDICT_VERSION @verdict_VERSION_FLAT@ 35 36#cmakedefine BUILD_SHARED_LIBS 37#ifdef BUILD_SHARED_LIBS 38# define VERDICT_SHARED_LIB 39#endif 40 41/*!\page VerdictAsASubProject Verdict as a Subproject 42 * 43 * Should a project wish to include a private version of verdict, it may 44 * turn on the CMake \c VERDICT_MANGLE option and define \c VERDICT_MANGLE_PREFIX. 45 * This will use the \c verdict_mangle.h file to define macros that rename all 46 * verdict library symbols using the given prefix. 47 * 48 * This allows programs using the parent project to link to the parent project 49 * as well as other versions of verdict than the one used by the parent project. 50 * If you are using mangling, you should consider setting \c VERDICT_NO_LIBRARY_VERSION 51 * to 1 and \c VERDICT_LIBRARY_PROPERTIES to include your parent project's shared 52 * library version number. 53 * 54 * Here is an example. Assume we have a library named XXX which we would like to use 55 * a private copy of the verdict library (instead of the version installed in a system 56 * directory). The XXX project source directory would contain a copy of the verdict 57 * source and the CMake file for XXX would include:\code 58 * set( VERDICT_MANGLE "ON" CACHE BOOL "XXX requires Verdict to be mangled" FORCE ) 59 * set( VERDICT_MANGLE_PREFIX "xxx" CACHE STRING "Verdict routines prefixed with xxx_" FORCE ) 60 * set( VERDICT_NO_LIBRARY_VERSION 1 ) 61 * set( VERDICT_LIBRARY_PROPERTIES 62 * ${VERDICT_LIBRARY_PROPERTIES} 63 * VERSION "xxx${XXX_VERSION}" 64 * ) 65 * add_subdirectory( verdict ) 66 * \endcode 67 * By setting these variables before descending into the verdict directory, 68 * we force all the symbol names in the library to be prefaced with \c xxx_ 69 * and the resulting library will be named \a xxxverdict instead of simply \a verdict. 70 */ 71#cmakedefine VERDICT_MANGLE 72#ifdef VERDICT_MANGLE 73# include "verdict_mangle.h" 74#endif 75 76#define VERDICT_DBL_MIN 1.0E-30 77#define VERDICT_DBL_MAX 1.0E+30 78#define VERDICT_PI 3.1415926535897932384626 79 80#if defined(_WIN32) || defined (__CYGWIN__) 81# define VERDICT_ABI_IMPORT __declspec(dllimport) 82# define VERDICT_ABI_EXPORT __declspec(dllexport) 83#elif __GNUC__ >= 4 84# define VERDICT_ABI_IMPORT __attribute__ ((visibility("default"))) 85# define VERDICT_ABI_EXPORT __attribute__ ((visibility("default"))) 86#else 87# define VERDICT_ABI_IMPORT 88# define VERDICT_ABI_EXPORT 89#endif 90 91#ifdef __cplusplus 92# define VERDICT_EXTERN_C extern "C" 93#endif 94 95#if defined(VERDICT_SHARED_LIB) 96# ifdef verdict_EXPORTS 97# define C_FUNC_DEF VERDICT_EXTERN_C VERDICT_ABI_EXPORT 98# else 99# define C_FUNC_DEF VERDICT_EXTERN_C VERDICT_ABI_IMPORT 100# endif 101#else 102# define C_FUNC_DEF VERDICT_EXTERN_C 103#endif 104 105 106/* typedef for the user if they want to use 107 * function pointers */ 108#ifdef __cplusplus 109extern "C" { 110#endif 111 typedef double(*VerdictFunction)(int, double[][3]); 112 typedef int(*ComputeNormal)(double point[3], double normal[3]); 113#ifdef __cplusplus 114} 115#endif 116 117 118/**\brief A <em>struct</em> to hold values computed by \ref v_hex_quality. 119 * 120 * HexMetricVals is used just as QuadMetricVals. 121 * For an example, see the \ref UsingQuadMetricVals. 122 */ 123struct HexMetricVals 124{ 125 /** \sa v_hex_edge_ratio */ 126 double edge_ratio ; 127 /** \sa v_hex_max_edge_ratio */ 128 double max_edge_ratio ; 129 /** \sa v_hex_skew */ 130 double skew ; 131 /** \sa v_hex_taper */ 132 double taper ; 133 /** \sa v_hex_volume */ 134 double volume ; 135 /** \sa v_hex_stretch */ 136 double stretch ; 137 /** \sa v_hex_diagonal */ 138 double diagonal ; 139 /** \sa v_hex_dimension */ 140 double dimension ; 141 /** \sa v_hex_oddy */ 142 double oddy ; 143 /** \sa v_hex_med_aspect_frobenius */ 144 double med_aspect_frobenius ; 145 /** \sa v_hex_max_aspect_frobenius */ 146 double max_aspect_frobenius ; 147 /** \sa v_hex_condition */ 148 double condition ; 149 /** \sa v_hex_jacobian */ 150 double jacobian ; 151 /** \sa v_hex_scaled_jacobian */ 152 double scaled_jacobian ; 153 /** \sa v_hex_shear */ 154 double shear ; 155 /** \sa v_hex_shape */ 156 double shape ; 157 /** \sa v_hex_relative_size */ 158 double relative_size_squared; 159 /** \sa v_hex_shape_and_size */ 160 double shape_and_size ; 161 /** \sa v_hex_shear_and_size */ 162 double shear_and_size ; 163 /** \sa v_hex_distortion */ 164 double distortion; 165}; 166 167/**\brief A <em>struct</em> to hold values computed by \ref v_edge_quality. 168 * 169 * EdgeMetricVals is used just as QuadMetricVals. 170 * For an example, see the \ref UsingQuadMetricVals. 171 */ 172struct EdgeMetricVals 173{ 174 double length ; 175}; 176 177 178/**\brief A <em>struct</em> to hold values computed by \ref v_knife_quality. 179 * 180 * KnifeMetricVals is used just as QuadMetricVals. 181 * For an example, see the \ref UsingQuadMetricVals. 182 */ 183struct KnifeMetricVals 184{ 185 double volume ; 186}; 187 188/**\brief A <em>struct</em> to hold values computed by \ref v_quad_quality. 189 * 190 * \section UsingQuadMetricVals Example Usage of the QuadMetricVals Structure 191 * 192 * The following is an example of how this struct is used with Verdict. 193 * 194 * Example: \code 195 * QuadMetricVals quad_metrics = {0}; 196 * unsigned long metrics_flag = 0; 197 * metrics_flag += V_QUAD_SHAPE; 198 * metrics_flag += V_QUAD_DISTORTION; 199 * metrics_flag += V_QUAD_AREA; 200 * double quad_nodes[4][3]; 201 * get_quad_nodes( quad_nodes ); //some user-defined function to load 202 * //xyz coordinate info. into array 203 * v_quad_quality( 4, quad_nodes, metrics_flag, quad_metrics ); 204 * double my_shape = quad_metrics.shape; 205 * double my_distortion = quad_metrics.distortion; 206 * double my_area = quad_metrics.area; \endcode 207 * 208 */ 209struct QuadMetricVals 210{ 211 /** \sa v_quad_edge_ratio function */ 212 double edge_ratio ; 213 /** \sa v_quad_max_edge_ratio function */ 214 double max_edge_ratio ; 215 /** \sa v_quad_aspect_ratio function */ 216 double aspect_ratio ; 217 /** \sa v_quad_radius_ratio function */ 218 double radius_ratio ; 219 /** \sa v_quad_med_aspect_frobenius function */ 220 double med_aspect_frobenius ; 221 /** \sa v_quad_max_aspect_frobenius function */ 222 double max_aspect_frobenius ; 223 /** \sa v_quad_skew function */ 224 double skew ; 225 /** \sa v_quad_taper function */ 226 double taper ; 227 /** \sa v_quad_warpage function */ 228 double warpage ; 229 /** \sa v_quad_area function */ 230 double area ; 231 /** \sa v_quad_stretch function */ 232 double stretch ; 233 /** \sa v_quad_smallest_angle function */ 234 double minimum_angle ; 235 /** \sa v_quad_largest_angle function */ 236 double maximum_angle ; 237 /** \sa v_quad_oddy function */ 238 double oddy ; 239 /** \sa v_quad_condition function */ 240 double condition ; 241 /** \sa v_quad_jacobian function */ 242 double jacobian ; 243 /** \sa v_quad_scaled_jacobian function */ 244 double scaled_jacobian ; 245 /** \sa v_quad_shear function */ 246 double shear ; 247 /** \sa v_quad_shape function */ 248 double shape ; 249 /** \sa v_quad_relative_size_squared function */ 250 double relative_size_squared ; 251 /** \sa v_quad_shape_and_size function */ 252 double shape_and_size ; 253 /** \sa v_quad_shear_and_size function */ 254 double shear_and_size ; 255 /** \sa v_quad_distortion function */ 256 double distortion; 257}; 258 259/**\brief A <em>struct</em> to hold values computed by \ref v_pyramid_quality. 260 * 261 * PyramidMetricVals is used just as QuadMetricVals. 262 * For an example, see the \ref UsingQuadMetricVals. 263 */ 264struct PyramidMetricVals 265{ 266 double volume ; 267}; 268 269/**\brief A <em>struct</em> to hold values computed by \ref v_wedge_quality. 270 * 271 * WedgeMetricVals is used just as QuadMetricVals. 272 * For an example, see the \ref UsingQuadMetricVals. 273 */ 274struct WedgeMetricVals 275{ 276 /** \sa v_wedge_volume */ 277 double volume ; 278 /** \sa v_wedge_edge_ratio */ 279 double edge_ratio ; 280 /** \sa v_wedge_max_aspect_frobenius */ 281 double max_aspect_frobenius; 282 /** \sa v_wedge_mean_aspect_frobenius */ 283 double mean_aspect_frobenius; 284 /** \sa v_wedge_jacobian */ 285 double jacobian; 286 /** \sa v_wedge_distortion */ 287 double distortion; 288 /** \sa v_wedge_max_stretch */ 289 double max_stretch; 290 /** \sa v_wedge_scaled_jacobian */ 291 double scaled_jacobian; 292 /** \sa v_wedge_shape */ 293 double shape; 294 /** \sa v_wedge_condition */ 295 double condition; 296}; 297 298/**\brief A <em>struct</em> to hold values computed by \ref v_tet_quality. 299 * 300 * TetMetricVals is used just as QuadMetricVals. 301 * For an example, see the \ref UsingQuadMetricVals. 302 */ 303struct TetMetricVals 304{ 305 /** \sa v_tet_edge_ratio*/ 306 double edge_ratio; 307 /** \sa v_tet_radius_ratio*/ 308 double radius_ratio; 309 /** \sa v_tet_aspect_beta*/ 310 double aspect_beta; 311 /** \sa v_tet_aspect_ratio */ 312 double aspect_ratio ; 313 /** \sa v_tet_aspect_gamma */ 314 double aspect_gamma ; 315 /** \sa v_tet_aspect_frobenius */ 316 double aspect_frobenius ; 317 /** \sa v_tet_minimum_angle */ 318 double minimum_angle ; 319 /** \sa v_tet_collapse_ratio*/ 320 double collapse_ratio; 321 /** \sa v_tet_volume */ 322 double volume ; 323 /** \sa v_tet_condition */ 324 double condition ; 325 /** \sa v_tet_jacobian */ 326 double jacobian ; 327 /** \sa v_tet_scaled_jacobian */ 328 double scaled_jacobian ; 329 /** \sa v_tet_shape */ 330 double shape ; 331 /** \sa v_tet_relative_size */ 332 double relative_size_squared ; 333 /** \sa v_tet_shape_and_size*/ 334 double shape_and_size ; 335 /** \sa v_tet_distortion */ 336 double distortion; 337 /** \sa v_tet_equivolume_skew */ 338 double equivolume_skew; 339 /** \sa v_tet_squish_index */ 340 double squish_index; 341 342}; 343 344/**\brief A <em>struct</em> to hold values computed by \ref v_tri_quality. 345 * 346 * TriMetricVals is used just as QuadMetricVals. 347 * For an example, see the \ref UsingQuadMetricVals. 348 */ 349struct TriMetricVals 350{ 351 /** \sa v_tri_edge_ratio */ 352 double edge_ratio ; 353 /** \sa v_tri_aspect_ratio */ 354 double aspect_ratio ; 355 /** \sa v_tri_radius_ratio */ 356 double radius_ratio ; 357 /** \sa v_tri_aspect_frobenius */ 358 double aspect_frobenius ; 359 /** \sa v_tri_area*/ 360 double area ; 361 /** \sa v_tri_minimum_angle*/ 362 double minimum_angle ; 363 /** \sa v_tri_maximum_angle */ 364 double maximum_angle ; 365 /** \sa v_tri_condition */ 366 double condition ; 367 /** \sa v_tri_scaled_jacobian */ 368 double scaled_jacobian ; 369 /** \sa v_tri_shape */ 370 double shape ; 371 /** \sa v_tri_relative_size_squared */ 372 double relative_size_squared ; 373 /** \sa v_tri_shape_and_size */ 374 double shape_and_size ; 375 /** \sa v_tri_distortion */ 376 double distortion; 377}; 378 379 380 381/* definition of bit fields to determine which metrics to calculate */ 382 383//! \name Atomic hex quality flags 384//! 385//@{ 386 387/*!\hideinitializer \brief Request that the maximum edge ratio be computed. \sa v_hex_max_edge_ratio */ 388#define V_HEX_MAX_EDGE_RATIO 1 389/*!\hideinitializer \brief Request that skew be computed. \sa v_hex_skew */ 390#define V_HEX_SKEW 2 391/*!\hideinitializer \brief Request that taper be computed. \sa v_hex_taper */ 392#define V_HEX_TAPER 4 393/*!\hideinitializer \brief Request that volume be computed. \sa v_hex_volume */ 394#define V_HEX_VOLUME 8 395/*!\hideinitializer \brief Request that stretch be computed. \sa v_hex_stretch */ 396#define V_HEX_STRETCH 16 397/*!\hideinitializer \brief Request that the diagonal be computed. \sa v_hex_diagonal */ 398#define V_HEX_DIAGONAL 32 399/*!\hideinitializer \brief Request that the dimension be computed. \sa v_hex_dimension */ 400#define V_HEX_DIMENSION 64 401/*!\hideinitializer \brief Request that the Oddy quality be computed. \sa v_hex_oddy */ 402#define V_HEX_ODDY 128 403/*!\hideinitializer \brief Request that maximum Frobenius aspect be computed. \sa v_hex_condition */ 404#define V_HEX_MAX_ASPECT_FROBENIUS 256 405/*!\hideinitializer \brief Request that the condition be computed. \sa v_hex_condition */ 406#define V_HEX_CONDITION 256 407/*!\hideinitializer \brief Request that the Jacobian be computed. \sa v_hex_jacobian */ 408#define V_HEX_JACOBIAN 512 409/*!\hideinitializer \brief Request that the scaled Jacobian be computed. \sa v_hex_scaled_jacobian */ 410#define V_HEX_SCALED_JACOBIAN 1024 411/*!\hideinitializer \brief Request that shear be computed. \sa v_hex_shear */ 412#define V_HEX_SHEAR 2048 413/*!\hideinitializer \brief Request that shape be computed. \sa v_hex_shape */ 414#define V_HEX_SHAPE 4096 415/*!\hideinitializer \brief Request that the square of the relative size be computed. \sa v_hex_relative_size_squared */ 416#define V_HEX_RELATIVE_SIZE_SQUARED 8192 417/*!\hideinitializer \brief Request that the product of shape and size be computed. \sa v_hex_shape_and_size */ 418#define V_HEX_SHAPE_AND_SIZE 16384 419/*!\hideinitializer \brief Request that the product of shear and size be computed. \sa v_hex_shear_and_size */ 420#define V_HEX_SHEAR_AND_SIZE 32768 421/*!\hideinitializer \brief Request that distortion be computed. \sa v_hex_distortion */ 422#define V_HEX_DISTORTION 65536 423/*!\hideinitializer \brief Request that the edge ratio be computed. \sa v_hex_edge_ratio */ 424#define V_HEX_EDGE_RATIO 131072 425/*!\hideinitializer \brief Request that average Frobenius aspect be computed. \sa v_hex_med_aspect_frobenius */ 426#define V_HEX_MED_ASPECT_FROBENIUS 262144 427//!\name Composite hex quality flags 428//! These flags are bitwise combinations of the atomic flags above. 429//@{ 430/*!\hideinitializer \brief Request that all hex metrics be computed. */ 431#define V_HEX_ALL 524287 432/*!\hideinitializer \brief Request the traditional hex qualities be computed. 433 * 434 * This includes 435 * \ref V_HEX_MAX_EDGE_RATIO, \ref V_HEX_SKEW, \ref V_HEX_TAPER 436 * \ref V_HEX_STRETCH, \ref V_HEX_DIAGONAL, \ref V_HEX_ODDY, 437 * \ref V_HEX_CONDITION, \ref V_HEX_JACOBIAN, \ref V_HEX_SCALED_JACOBIAN, and 438 * \ref V_HEX_DIMENSION. 439 */ 440#define V_HEX_TRADITIONAL \ 441 V_HEX_MAX_EDGE_RATIO + \ 442 V_HEX_SKEW + \ 443 V_HEX_TAPER + \ 444 V_HEX_STRETCH + \ 445 V_HEX_DIAGONAL + \ 446 V_HEX_ODDY + \ 447 V_HEX_CONDITION + \ 448 V_HEX_JACOBIAN + \ 449 V_HEX_SCALED_JACOBIAN + \ 450 V_HEX_DIMENSION 451/*!\hideinitializer \brief Request hex qualities used to diagnose mesh or solver problems. 452 * 453 * This is currently a synonym for \ref V_HEX_VOLUME but may include more qualities in the future. 454 */ 455#define V_HEX_DIAGNOSTIC V_HEX_VOLUME 456/*!\hideinitializer \brief Request hex qualities that have an algebraic expression be computed. 457 * 458 * This includes 459 * \ref V_HEX_SHAPE, \ref V_HEX_SHEAR, \ref V_HEX_RELATIVE_SIZE_SQUARED, 460 * \ref V_HEX_SHAPE_AND_SIZE, and \ref V_HEX_SHEAR_AND_SIZE. 461 */ 462#define V_HEX_ALGEBRAIC \ 463 V_HEX_SHAPE + \ 464 V_HEX_SHEAR + \ 465 V_HEX_RELATIVE_SIZE_SQUARED + \ 466 V_HEX_SHAPE_AND_SIZE + \ 467 V_HEX_SHEAR_AND_SIZE 468/*!\hideinitializer \brief Request hex qualities for Mr. Robinson. 469 */ 470#define V_HEX_ROBINSON \ 471 V_HEX_SKEW + \ 472 V_HEX_TAPER 473//@} 474//@} 475 476//! \name Atomic tetrahedral quality flags 477//! 478//@{ 479/*!\hideinitializer \brief Request that the radius ratio be computed. \sa v_tet_radius_ratio */ 480#define V_TET_RADIUS_RATIO 1 481/*!\hideinitializer \brief Request that the aspect ratio beta be computed. \sa v_tet_aspect_beta */ 482#define V_TET_ASPECT_BETA 2 483/*!\hideinitializer \brief Request that the aspect ratio gamma be computed. \sa v_tet_aspect_gamma */ 484#define V_TET_ASPECT_GAMMA 4 485/*!\hideinitializer \brief Request that the volume be computed. \sa v_tet_volume */ 486#define V_TET_VOLUME 8 487/*!\hideinitializer \brief Request that the condition be computed. \sa v_tet_condition */ 488#define V_TET_CONDITION 16 489/*!\hideinitializer \brief Request that the Jacobian be computed. \sa v_tet_jacobian */ 490#define V_TET_JACOBIAN 32 491/*!\hideinitializer \brief Request that the scaled Jacobian be computed. \sa v_tet_scaled_jacobian */ 492#define V_TET_SCALED_JACOBIAN 64 493/*!\hideinitializer \brief Request that the shape be computed. \sa v_tet_shape */ 494#define V_TET_SHAPE 128 495/*!\hideinitializer \brief Request that the square of the relative size be computed. \sa v_tet_relative_size_squared */ 496#define V_TET_RELATIVE_SIZE_SQUARED 256 497/*!\hideinitializer \brief Request that the product of shape and size be computed. \sa v_tet_shape_and_size */ 498#define V_TET_SHAPE_AND_SIZE 512 499/*!\hideinitializer \brief Request that the distortion be computed. \sa v_tet_distortion */ 500#define V_TET_DISTORTION 1024 501/*!\hideinitializer \brief Request that the ratio of edge lengths be computed. \sa v_tet_edge_ratio */ 502#define V_TET_EDGE_RATIO 2048 503/*!\hideinitializer \brief Request that the aspect ratio be computed. \sa v_tet_aspect_ratio */ 504#define V_TET_ASPECT_RATIO 4096 505/*!\hideinitializer \brief Request that the Frobenius aspect be computed. \sa v_tet_aspect_frobenius */ 506#define V_TET_ASPECT_FROBENIUS 8192 507/*!\hideinitializer \brief Request that the minimum dihedral angle be computed. \sa v_tet_minimum_angle */ 508#define V_TET_MINIMUM_ANGLE 16384 509/*!\hideinitializer \brief Request that the collapse ratio be computed. \sa v_tet_collapse_ratio */ 510#define V_TET_COLLAPSE_RATIO 32768 511/*!\hideinitializer \brief Request that the equivolume skew as defined by Fluent be computed. \sa v_tet_equivolume_skew */ 512#define V_TET_EQUIVOLUME_SKEW 65536 513/*!\hideinitializer \brief Request that the squish index as defined by Fluent be computed. \sa v_tet_squish_index */ 514#define V_TET_SQUISH_INDEX 131072 515 516//! Composite tetrahedral quality flags 517//! These flags are bitwise combinations of the atomic flags above. 518//@{ 519/*!\hideinitializer \brief Request that the all tetrahedral qualities be computed. */ 520#define V_TET_ALL 262143 521/*!\hideinitializer \brief Request that the tradtional tetrahedral qualities be computed. 522 * 523 * This includes 524 * \ref V_TET_ASPECT_BETA, \ref V_TET_ASPECT_GAMMA, \ref V_TET_CONDITION, 525 * \ref V_TET_JACOBIAN, and \ref V_TET_SCALED_JACOBIAN. 526 */ 527#define V_TET_TRADITIONAL V_TET_ASPECT_BETA + \ 528 V_TET_ASPECT_GAMMA + \ 529 V_TET_CONDITION + \ 530 V_TET_JACOBIAN + \ 531 V_TET_SCALED_JACOBIAN 532/*!\hideinitializer \brief Request hex qualities used to diagnose mesh or solver problems. 533 * 534 * This is currently a synonym for \ref V_TET_VOLUME but may include more qualities in the future. 535 */ 536#define V_TET_DIAGNOSTIC V_TET_VOLUME 537/*!\hideinitializer \brief Request hex qualities that have an algebraic expression be computed. 538 * 539 * This includes 540 * \ref V_TET_SHAPE, \ref V_TET_RELATIVE_SIZE_SQUARED, and \ref V_TET_SHAPE_AND_SIZE. 541 */ 542#define V_TET_ALGEBRAIC V_TET_SHAPE + \ 543 V_TET_RELATIVE_SIZE_SQUARED + \ 544 V_TET_SHAPE_AND_SIZE 545//@} 546//@} 547 548 549//! \name Pyramid bit fields 550//! 551//@{ 552/*!\hideinitializer \brief */ 553#define V_PYRAMID_VOLUME 1 554/*!\hideinitializer \brief */ 555#define V_PYRAMID_ALL \ 556 V_PYRAMID_VOLUME 557//@} 558 559//! \name Wedge bit fields 560//! 561//@{ 562/*!\hideinitializer \brief */ 563#define V_WEDGE_VOLUME 1 564#define V_WEDGE_EDGE_RATIO 2 565#define V_WEDGE_MAX_ASPECT_FROBENIUS 4 566#define V_WEDGE_MEAN_ASPECT_FROBENIUS 8 567#define V_WEDGE_JACOBIAN 16 568#define V_WEDGE_SCALED_JACOBIAN 32 569#define V_WEDGE_DISTORTION 64 570#define V_WEDGE_MAX_STRETCH 128 571#define V_WEDGE_SHAPE 256 572#define V_WEDGE_CONDITION 512 573/*!\hideinitializer \brief */ 574#define V_WEDGE_ALL \ 575 V_WEDGE_VOLUME 576//@} 577 578//! \name Knife bit fields 579//! 580//@{ 581/*!\hideinitializer \brief */ 582#define V_KNIFE_VOLUME 1 583/*!\hideinitializer \brief */ 584#define V_KNIFE_ALL \ 585 V_KNIFE_VOLUME 586//@} 587 588//! \name Quad bit fields 589//! 590//@{ 591/*!\hideinitializer \brief */ 592#define V_QUAD_MAX_EDGE_RATIO 1 593/*!\hideinitializer \brief */ 594#define V_QUAD_SKEW 2 595/*!\hideinitializer \brief */ 596#define V_QUAD_TAPER 4 597/*!\hideinitializer \brief */ 598#define V_QUAD_WARPAGE 8 599/*!\hideinitializer \brief */ 600#define V_QUAD_AREA 16 601/*!\hideinitializer \brief */ 602#define V_QUAD_STRETCH 32 603/*!\hideinitializer \brief */ 604#define V_QUAD_MINIMUM_ANGLE 64 605/*!\hideinitializer \brief */ 606#define V_QUAD_MAXIMUM_ANGLE 128 607/*!\hideinitializer \brief */ 608#define V_QUAD_ODDY 256 609/*!\hideinitializer \brief */ 610#define V_QUAD_CONDITION 512 611/*!\hideinitializer \brief */ 612#define V_QUAD_JACOBIAN 1024 613/*!\hideinitializer \brief */ 614#define V_QUAD_SCALED_JACOBIAN 2048 615/*!\hideinitializer \brief */ 616#define V_QUAD_SHEAR 4096 617/*!\hideinitializer \brief */ 618#define V_QUAD_SHAPE 8192 619/*!\hideinitializer \brief */ 620#define V_QUAD_RELATIVE_SIZE_SQUARED 16384 621/*!\hideinitializer \brief */ 622#define V_QUAD_SHAPE_AND_SIZE 32768 623/*!\hideinitializer \brief */ 624#define V_QUAD_SHEAR_AND_SIZE 65536 625/*!\hideinitializer \brief */ 626#define V_QUAD_DISTORTION 131072 627/*!\hideinitializer \brief */ 628#define V_QUAD_EDGE_RATIO 262144 629/*!\hideinitializer \brief */ 630#define V_QUAD_ASPECT_RATIO 524288 631/*!\hideinitializer \brief */ 632#define V_QUAD_RADIUS_RATIO 1048576 633/*!\hideinitializer \brief */ 634#define V_QUAD_MED_ASPECT_FROBENIUS 2097152 635/*!\hideinitializer \brief */ 636#define V_QUAD_MAX_ASPECT_FROBENIUS 4194304 637/*!\hideinitializer \brief */ 638#define V_QUAD_ALL 8388607 639/*!\hideinitializer \brief */ 640#define V_QUAD_TRADITIONAL \ 641 V_QUAD_MAX_EDGE_RATIO + \ 642 V_QUAD_SKEW + \ 643 V_QUAD_TAPER + \ 644 V_QUAD_WARPAGE + \ 645 V_QUAD_STRETCH + \ 646 V_QUAD_MINIMUM_ANGLE + \ 647 V_QUAD_MAXIMUM_ANGLE + \ 648 V_QUAD_ODDY + \ 649 V_QUAD_CONDITION + \ 650 V_QUAD_JACOBIAN + \ 651 V_QUAD_SCALED_JACOBIAN 652/*!\hideinitializer \brief */ 653#define V_QUAD_DIAGNOSTIC \ 654 V_QUAD_AREA 655/*!\hideinitializer \brief */ 656#define V_QUAD_ALGEBRAIC \ 657 V_QUAD_SHEAR + \ 658 V_QUAD_SHAPE + \ 659 V_QUAD_RELATIVE_SIZE_SQUARED + \ 660 V_QUAD_SHAPE_AND_SIZE 661/*!\hideinitializer \brief */ 662#define V_QUAD_ROBINSON \ 663 V_QUAD_MAX_EDGE_RATIO + \ 664 V_QUAD_SKEW + \ 665 V_QUAD_TAPER 666//@} 667 668 669//! \name Tri bit fields 670//! 671//@{ 672/*!\hideinitializer \brief */ 673#define V_TRI_ASPECT_FROBENIUS 1 674/*!\hideinitializer \brief */ 675#define V_TRI_AREA 2 676/*!\hideinitializer \brief */ 677#define V_TRI_MINIMUM_ANGLE 4 678/*!\hideinitializer \brief */ 679#define V_TRI_MAXIMUM_ANGLE 8 680/*!\hideinitializer \brief */ 681#define V_TRI_CONDITION 16 682/*!\hideinitializer \brief */ 683#define V_TRI_SCALED_JACOBIAN 32 684/*!\hideinitializer \brief */ 685#define V_TRI_SHAPE 64 686/*!\hideinitializer \brief */ 687#define V_TRI_RELATIVE_SIZE_SQUARED 128 688/*!\hideinitializer \brief */ 689#define V_TRI_SHAPE_AND_SIZE 256 690/*!\hideinitializer \brief */ 691#define V_TRI_DISTORTION 512 692/*!\hideinitializer \brief */ 693#define V_TRI_RADIUS_RATIO 1024 694/*!\hideinitializer \brief */ 695#define V_TRI_EDGE_RATIO 2048 696/*!\hideinitializer \brief */ 697#define V_TRI_ALL 4095 698/*!\hideinitializer \brief */ 699#define V_TRI_TRADITIONAL \ 700 V_TRI_ASPECT_FROBENIUS + \ 701 V_TRI_MINIMUM_ANGLE + \ 702 V_TRI_MAXIMUM_ANGLE + \ 703 V_TRI_CONDITION + \ 704 V_TRI_SCALED_JACOBIAN 705/*!\hideinitializer \brief */ 706#define V_TRI_DIAGNOSTIC \ 707 V_TRI_AREA 708/*!\hideinitializer \brief */ 709#define V_TRI_ALGEBRAIC \ 710 V_TRI_SHAPE + \ 711 V_TRI_SHAPE_AND_SIZE + \ 712 V_TRI_RELATIVE_SIZE_SQUARED 713//@} 714 715//! \name Tri bit fields 716//! 717//@{ 718/*!\hideinitializer \brief */ 719#define V_EDGE_LENGTH 1 720/*!\hideinitializer \brief */ 721#define V_EDGE_ALL \ 722 V_EDGE_LENGTH 723//@} 724 725 726/*! \mainpage 727 Verdict is a library used to calculate metrics on the following type of elements: 728 729 \li Hexahedra 730 \li Tetrahedra 731 \li Pyramid 732 \li Wedge 733 \li Knife 734 \li Quadrilateral 735 \li Triangle 736 \li Edge 737 738 Verdict calculates individual or multiple metrics on a single elment. 739 The v_*_quality(...) functions allow for efficient calculations of 740 multiple metrics on a single element. Individual metrics may be 741 calculated on a single element as well. 742 743 \section GetVerdict Obtaining, Configuring, Building, and Installing Verdict 744 745 Since your are reading reference documentation generated from the source code, 746 you ostensibly already have Verdict. 747 Instructions for this are on a separate page: \ref ObtainConfigureBuildInstall . 748 749 \section UsingVerdict Using Verdict 750 751 The v_*_quality functions take the following parameters: 752 753 \param num_nodes Number of nodes in the element. 754 \param coordinates 2D array containing x,y,z coordinate data of the nodes. 755 \param metrics_request_flag Bitfield to define which metrics to calculate. 756 \param metric_vals Struct to hold the metric calculated values. 757 758 All other functions take these parameters below and return the calculated 759 metric value. 760 761 \param num_nodes Number of nodes in the element. 762 \param coordinates 2D array containing x,y,z coordinate data of the nodes. 763 764 765 \par Setting the metrics_request_flag: 766 In order to use v_*_quality functions you must know how to set the bitfield argument 767 correctly. To calculate aspect ratio, condition number, shape and shear metrics on a triangle, set 768 the "metrics_request_flag" like so: 769 770 \code 771 unsigned int metrics_request_flag = 0; 772 metrics_request_flag += V_TRI_ASPECT_FROBENIUS; 773 metrics_request_flag += V_TRI_CONDITION; 774 metrics_request_flag += V_TRI_SHAPE; 775 metrics_request_flag += V_TRI_SHEAR; 776 \endcode 777 778 The bitwise field can also be set for many metrics at once using summary tags such as \ref V_HEX_ALL, 779 \ref V_HEX_DIAGNOSTIC, or \ref V_TRI_ALGEBRAIC. 780 781 Below is an example of how use Verdict's functions: 782 783 784 Example: \code 785 QuadMetricVals quad_metrics = {0}; 786 unsigned long metrics_flag = 0; 787 metrics_flag += V_QUAD_SHAPE; 788 metrics_flag += V_QUAD_DISTORTION; 789 metrics_flag += V_QUAD_AREA; 790 double quad_nodes[4][3]; 791 792 //node 1 793 quad_node[0][0] = 0; //x 794 quad_node[0][1] = 0; //y 795 quad_node[0][2] = 0; //z 796 797 //node 2 798 quad_node[1][0] = 1; 799 quad_node[1][1] = 0.1; 800 quad_node[1][2] = 0.1; 801 802 //node 3 803 quad_node[2][0] = 0.9; 804 quad_node[2][1] = 0.9; 805 quad_node[2][2] = -0.1; 806 807 //node 4 808 quad_node[3][0] = -0.05; 809 quad_node[3][1] = 1; 810 quad_node[3][2] = 0; 811 812 //calculate multiple metrics with one call 813 v_quad_quality( 4, quad_nodes, metrics_flag, quad_metrics ); 814 double my_shape = quad_metrics.shape; 815 double my_distortion = quad_metrics.distortion; 816 double my_area = quad_metrics.area; 817 818 //calculate an individual metric 819 double my_relative_size = v_quad_relative_size( 4, quad_nodes ); 820 \endcode 821 822*/ 823 824/*!\page ObtainConfigureBuildInstall Obtain, Configure, Build, and Install Verdict 825 826 \section ObtainVerdict Obtaining the Verdict Source Code 827 828 The verdict source code repository is now maintained by <a href="http://www.kitware.com/">Kitware</a>. 829 You can check out the source with \code 830 cvs -d :pserver:anonymous@public.kitware.com:/cvsroot/Verdict login 831 cvs -d :pserver:anonymous@public.kitware.com:/cvsroot/Verdict -z3 co Verdict 832 \endcode 833 When asked for a password, enter "verdict". 834 835 \section ConfiguringVerdict Configuring Verdict 836 837 Verdict uses <a href="http://www.cmake.org/">CMake</a> to create project files. 838 To build Verdict from source, you will need to obtain CMake. 839 840 We strongly suggest that you build Verdict in a separate directory tree from the source code. 841 Life is just too short to explain why; if you're curious, search the CMake mailing list. 842 To make the instructions clear, we'll assume that you have the Verdict source code in \c /tmp/Verdict. 843 Configuring Verdict is then as simple as \code 844 mkdir /tmp/VerdictBuild 845 cd /tmp/VerdictBuild 846 cmake ../Verdict 847 \endcode 848 At this point, you may wish to edit <tt>/tmp/VerdictBuild/CMakeCache.txt</tt> to change settings 849 from their default values. 850 851 If you are going to include Verdict inside another library, you might want to read 852 how to use \ref VerdictAsASubProject. 853 854 \section BuildingAndInstallingVerdict Building and Installing Verdict 855 856 Assuming you are using CMake's Makefile generator, building and installing Verdict 857 is accomplished with \code 858 make 859 make install 860 \endcode 861 By default, Verdict will be installed in \c /usr/local. 862 If you wish to change the location, edit your <tt>CMakeCache.txt</tt> file so that 863 \c CMAKE_INSTALL_PREFIX contains a different directory and re-run \c cmake. 864 */ 865 866 //! Calculates quality metrics for hexahedral elements. 867 C_FUNC_DEF void v_hex_quality( int num_nodes, double coordinates[][3], 868 unsigned int metrics_request_flag, struct HexMetricVals *metric_vals ); 869 870 //! Calculates quality metrics for tetrahedral elements. 871 C_FUNC_DEF void v_tet_quality( int num_nodes, double coordinates[][3], 872 unsigned int metrics_request_flag, struct TetMetricVals *metric_vals ); 873 874 //! Calculates quality metrics for pyramid elements. 875 C_FUNC_DEF void v_pyramid_quality( int num_nodes, double coordinates[][3], 876 unsigned int metrics_request_flag, struct PyramidMetricVals *metric_vals ); 877 878 //! Calculates quality metrics for wedge elements. 879 C_FUNC_DEF void v_wedge_quality( int num_nodes, double coordinates[][3], 880 unsigned int metrics_request_flag, struct WedgeMetricVals *metric_vals ); 881 882 //! Calculates quality metrics for knife elements. 883 C_FUNC_DEF void v_knife_quality( int num_nodes, double coordinates[][3], 884 unsigned int metrics_request_flag, struct KnifeMetricVals *metric_vals ); 885 886 //! Calculates quality metrics for quadralateral elements. 887 C_FUNC_DEF void v_quad_quality( int num_nodes, double coordinates[][3], 888 unsigned int metrics_request_flag, struct QuadMetricVals *metric_vals ); 889 890 //! Calculates quality metrics for triangle elements. 891 C_FUNC_DEF void v_tri_quality( int num_nodes, double coordinates[][3], 892 unsigned int metrics_request_flag, struct TriMetricVals *metric_vals ); 893 894 //! Calculates quality metrics for edge elements. 895 C_FUNC_DEF void v_edge_quality( int num_nodes, double coordinates[][3], 896 unsigned int metrics_request_flag, struct EdgeMetricVals *metric_vals ); 897 898 899 900/* individual quality functions for hex elements */ 901 902 //! Sets average size (volume) of hex, needed for v_hex_relative_size(...) 903 C_FUNC_DEF void v_set_hex_size( double size ); 904 905 //! Calculates hex edge ratio metric. 906 /** Hmax / Hmin where Hmax and Hmin are respectively the maximum and the 907 minimum edge lengths */ 908 C_FUNC_DEF double v_hex_edge_ratio( int num_nodes, double coordinates[][3] ); 909 910 //! Calculates hex maximum of edge ratio 911 /**Maximum edge length ratio at hex center. 912 Reference --- L.M. Taylor, and D.P. Flanagan, Pronto3D - A Three Dimensional Transient 913 Solid Dynamics Program, SAND87-1912, Sandia National Laboratories, 1989. */ 914 C_FUNC_DEF double v_hex_max_edge_ratio( int num_nodes, double coordinates[][3] ); 915 916 //! Calculates hex skew metric. 917 /** Maximum |cos A| where A is the angle between edges at hex center. 918 Reference --- L.M. Taylor, and D.P. Flanagan, Pronto3D - A Three Dimensional Transient 919 Solid Dynamics Program, SAND87-1912, Sandia National Laboratories, 1989. */ 920 C_FUNC_DEF double v_hex_skew( int num_nodes, double coordinates[][3] ); 921 922 //! Calculates hex taper metric 923 /** Maximum ratio of lengths derived from opposite edges. 924 Reference --- L.M. Taylor, and D.P. Flanagan, Pronto3D - A Three Dimensional Transient 925 Solid Dynamics Program, SAND87-1912, Sandia National Laboratories, 1989. */ 926 C_FUNC_DEF double v_hex_taper( int num_nodes, double coordinates[][3] ); 927 928 //! Calculates hex volume 929 /** Jacobian at hex center. 930 Reference --- L.M. Taylor, and D.P. Flanagan, Pronto3D - A Three Dimensional Transient 931 Solid Dynamics Program, SAND87-1912, Sandia National Laboratories, 1989. */ 932 C_FUNC_DEF double v_hex_volume( int num_nodes, double coordinates[][3] ); 933 934 //! Calculates hex stretch metric 935 /** Sqrt(3) * minimum edge length / maximum diagonal length. 936 Reference --- FIMESH code */ 937 C_FUNC_DEF double v_hex_stretch( int num_nodes, double coordinates[][3] ); 938 939 //! Calculates hex diagonal metric 940 /** Minimum diagonal length / maximum diagonal length. 941 Reference --- Unknown */ 942 C_FUNC_DEF double v_hex_diagonal( int num_nodes, double coordinates[][3] ); 943 944 //! Calculates hex dimension metric 945 /** Pronto-specific characteristic length for stable time step calculation. 946 Char_length = Volume / 2 grad Volume. 947 Reference --- L.M. Taylor, and D.P. Flanagan, Pronto3D - A Three Dimensional Transient 948 Solid Dynamics Program, SAND87-1912, Sandia National Laboratories, 1989. */ 949 C_FUNC_DEF double v_hex_dimension( int num_nodes, double coordinates[][3] ); 950 951 //! Calculates hex oddy metric 952 C_FUNC_DEF double v_hex_oddy( int num_nodes, double coordinates[][3] ); 953 954 //! Calculates hex condition metric 955 /** Average Frobenius condition number of the Jacobian matrix at 8 corners. */ 956 C_FUNC_DEF double v_hex_med_aspect_frobenius( int num_nodes, double coordinates[][3] ); 957 958 //! Calculates hex condition metric 959 /** Maximum Frobenius condition number of the Jacobian matrix at 8 corners. 960 Reference --- P. Knupp, Achieving Finite Element Mesh Quality via 961 Optimization of the Jacobian Matrix Norm and Associated Quantities, 962 Intl. J. Numer. Meth. Engng. 2000, 48:1165-1185. */ 963 C_FUNC_DEF double v_hex_max_aspect_frobenius( int num_nodes, double coordinates[][3] ); 964 //! Calculates hex condition metric. This is a synonym for \ref v_hex_max_aspect_frobenius. 965 C_FUNC_DEF double v_hex_condition( int num_nodes, double coordinates[][3] ); 966 967 //! Calculates hex jacobian metric 968 /** Minimum pointwise volume of local map at 8 corners & center of hex. 969 Reference --- P. Knupp, Achieving Finite Element Mesh Quality via 970 Optimization of the Jacobian Matrix Norm and Associated Quantities, 971 Intl. J. Numer. Meth. Engng. 2000, 48:1165-1185. */ 972 C_FUNC_DEF double v_hex_jacobian( int num_nodes, double coordinates[][3] ); 973 974 //! Calculates hex scaled jacobian metric 975 /** Minimum Jacobian divided by the lengths of the 3 edge vectors. 976 Reference --- P. Knupp, Achieving Finite Element Mesh Quality via 977 Optimization of the Jacobian Matrix Norm and Associated Quantities, 978 Intl. J. Numer. Meth. Engng. 2000, 48:1165-1185. */ 979 C_FUNC_DEF double v_hex_scaled_jacobian( int num_nodes, double coordinates[][3] ); 980 981 //! Calculates hex shear metric 982 /** 3/Mean Ratio of Jacobian Skew matrix. 983 Reference --- P. Knupp, Algebraic Mesh Quality Metrics for 984 Unstructured Initial Meshes, submitted for publication. */ 985 C_FUNC_DEF double v_hex_shear( int num_nodes, double coordinates[][3] ); 986 987 //! Calculates hex shape metric. 988 /** 3/Mean Ratio of weighted Jacobian matrix. 989 Reference --- P. Knupp, Algebraic Mesh Quality Metrics for 990 Unstructured Initial Meshes, submitted for publication. */ 991 C_FUNC_DEF double v_hex_shape( int num_nodes, double coordinates[][3] ); 992 993 //! Calculates hex relative size metric. 994 /** 3/Mean Ratio of weighted Jacobian matrix. 995 Reference --- P. Knupp, Algebraic Mesh Quality Metrics for 996 Unstructured Initial Meshes, submitted for publication. */ 997 C_FUNC_DEF double v_hex_relative_size_squared( int num_nodes, double coordinates[][3] ); 998 999 //! Calculates hex shape-size metric. 1000 /** Product of Shape and Relative Size. 1001 Reference --- P. Knupp, Algebraic Mesh Quality Metrics for 1002 Unstructured Initial Meshes, submitted for publication. */ 1003 C_FUNC_DEF double v_hex_shape_and_size( int num_nodes, double coordinates[][3] ); 1004 1005 //! Calculates hex shear-size metric 1006 /** Product of Shear and Relative Size. 1007 Reference --- P. Knupp, Algebraic Mesh Quality Metrics for 1008 Unstructured Initial Meshes, submitted for publication. */ 1009 C_FUNC_DEF double v_hex_shear_and_size( int num_nodes, double coordinates[][3] ); 1010 1011 //! Calculates hex distortion metric 1012 /** {min(|J|)/actual volume}*parent volume, parent volume = 8 for hex. 1013 Reference --- SDRC/IDEAS Simulation: Finite Element Modeling--User's Guide */ 1014 C_FUNC_DEF double v_hex_distortion( int num_nodes, double coordinates[][3] ); 1015 1016/* individual quality functions for tet elements */ 1017 1018 //! Sets average size (volume) of tet, needed for v_tet_relative_size(...) 1019 C_FUNC_DEF void v_set_tet_size( double size ); 1020 1021 //! Calculates tet edge ratio metric. 1022 /** Hmax / Hmin where Hmax and Hmin are respectively the maximum and the 1023 minimum edge lengths */ 1024 C_FUNC_DEF double v_tet_edge_ratio( int num_nodes, double coordinates[][3] ); 1025 1026 //! Calculates tet radius ratio metric. 1027 /** CR / (3.0 * IR) where CR = circumsphere radius, IR = inscribed sphere radius. 1028 Reference --- V. N. Parthasarathy et al, A comparison of tetrahedron 1029 quality measures, Finite Elem. Anal. Des., Vol 15(1993), 255-261. */ 1030 C_FUNC_DEF double v_tet_radius_ratio( int num_nodes, double coordinates[][3] ); 1031 1032 //! Calculates the radius ratio metric of a positively oriented tet. 1033 /** CR / (3.0 * IR) where CR = circumsphere radius, IR = inscribed sphere radius 1034 if the element is positively-oriented. 1035 Reference --- V. N. Parthasarathy et al, A comparison of tetrahedron 1036 quality measures, Finite Elem. Anal. Des., Vol 15(1993), 255-261. */ 1037 C_FUNC_DEF double v_tet_aspect_beta( int num_nodes, double coordinates[][3] ); 1038 1039 //! Calculates tet aspect ratio metric. 1040 /** Hmax / (2 sqrt(6) r) where Hmax and r respectively denote the greatest edge 1041 length and the inradius of the tetrahedron 1042 Reference --- P. Frey and P.-L. George, Meshing, Hermes (2000). */ 1043 C_FUNC_DEF double v_tet_aspect_ratio( int num_nodes, double coordinates[][3] ); 1044 1045 //! Calculates tet aspect gamma metric. 1046 /** Srms**3 / (8.479670*V) where Srms = sqrt(Sum(Si**2)/6), Si = edge length. 1047 Reference --- V. N. Parthasarathy et al, A comparison of tetrahedron 1048 quality measures, Finite Elem. Anal. Des., Vol 15(1993), 255-261. */ 1049 C_FUNC_DEF double v_tet_aspect_gamma( int num_nodes, double coordinates[][3] ); 1050 1051 //! Calculates tet aspect frobenius metric. 1052 /** Frobenius condition number when the reference element is regular 1053 Reference --- P. Knupp, Achieving Finite Element Mesh Quality via 1054 Optimization of the Jacobian Matrix Norm and Associated Quantities, 1055 Intl. J. Numer. Meth. Engng. 2000, 48:1165-1185. */ 1056 C_FUNC_DEF double v_tet_aspect_frobenius( int num_nodes, double coordinates[][3] ); 1057 1058 //! Calculates tet minimum dihedral angle. 1059 /** Minimum (nonoriented) dihedral angle of a tetrahedron, expressed in degrees. */ 1060 C_FUNC_DEF double v_tet_minimum_angle( int num_nodes, double coordinates[][3] ); 1061 1062 //! Calculates tet collapse ratio metric. 1063 /** Collapse ratio */ 1064 C_FUNC_DEF double v_tet_collapse_ratio( int num_nodes, double coordinates[][3] ); 1065 1066 //! Calculates tet volume. 1067 /** (1/6) * Jacobian at corner node. 1068 Reference --- V. N. Parthasarathy et al, A comparison of tetrahedron 1069 quality measures, Finite Elem. Anal. Des., Vol 15(1993), 255-261. */ 1070 C_FUNC_DEF double v_tet_volume( int num_nodes, double coordinates[][3] ); 1071 1072 //! Calculates tet condition metric. 1073 /** Condition number of the Jacobian matrix at any corner. 1074 Reference --- P. Knupp, Achieving Finite Element Mesh Quality via 1075 Optimization of the Jacobian Matrix Norm and Associated Quantities, 1076 Intl. J. Numer. Meth. Engng. 2000, 48:1165-1185. */ 1077 C_FUNC_DEF double v_tet_condition( int num_nodes, double coordinates[][3] ); 1078 1079 //! Calculates tet jacobian. 1080 /** Minimum pointwise volume at any corner. 1081 Reference --- P. Knupp, Achieving Finite Element Mesh Quality via 1082 Optimization of the Jacobian Matrix Norm and Associated Quantities, 1083 Intl. J. Numer. Meth. Engng. 2000, 48:1165-1185. */ 1084 C_FUNC_DEF double v_tet_jacobian( int num_nodes, double coordinates[][3] ); 1085 1086 //! Calculates tet scaled jacobian. 1087 /** Minimum Jacobian divided by the lengths of 3 edge vectors 1088 Reference --- P. Knupp, Achieving Finite Element Mesh Quality via 1089 Optimization of the Jacobian Matrix Norm and Associated Quantities, 1090 Intl. J. Numer. Meth. Engng. 2000, 48:1165-1185. */ 1091 C_FUNC_DEF double v_tet_scaled_jacobian( int num_nodes, double coordinates[][3] ); 1092 1093 //! Calculates tet shape metric. 1094 /** 3/Mean Ratio of weighted Jacobian matrix. 1095 Reference --- P. Knupp, Algebraic Mesh Quality Metrics for 1096 Unstructured Initial Meshes, submitted for publication. */ 1097 C_FUNC_DEF double v_tet_shape( int num_nodes, double coordinates[][3] ); 1098 1099 //! Calculates tet relative size metric. 1100 /** Min( J, 1/J ), where J is determinant of weighted Jacobian matrix. 1101 Reference --- P. Knupp, Algebraic Mesh Quality Metrics for 1102 Unstructured Initial Meshes, submitted for publication. */ 1103 C_FUNC_DEF double v_tet_relative_size_squared( int num_nodes, double coordinates[][3] ); 1104 1105 //! Calculates tet shape-size metric. 1106 /** Product of Shape and Relative Size. 1107 Reference --- P. Knupp, Algebraic Mesh Quality Metrics for 1108 Unstructured Initial Meshes, submitted for publication. */ 1109 C_FUNC_DEF double v_tet_shape_and_size( int num_nodes, double coordinates[][3] ); 1110 1111 //! Calculates tet distortion metric. 1112 /** {min(|J|)/actual volume}*parent volume, parent volume = 1/6 for tet. 1113 Reference --- SDRC/IDEAS Simulation: Finite Element Modeling--User's Guide */ 1114 C_FUNC_DEF double v_tet_distortion( int num_nodes, double coordinates[][3] ); 1115 1116/* individual quality functions for pyramid elements */ 1117 1118 //! Calculates pyramid volume. 1119 C_FUNC_DEF double v_pyramid_volume( int num_nodes, double coordinates[][3] ); 1120 1121 1122/* individual quality functions for wedge elements */ 1123 1124 //! Calculates wedge volume. 1125 C_FUNC_DEF double v_wedge_volume( int num_nodes, double coordinates[][3] ); 1126 1127 //! Calculates wedge edge ratio metric. 1128 /** Hmax / Hmin where Hmax and Hmin are respectively the maximum and the 1129 minimum edge lengths */ 1130 C_FUNC_DEF double v_wedge_edge_ratio( int num_nodes, double coordinates[][3] ); 1131 1132 //! Calculates wedge max aspect forbenius. 1133 /** max(F_0123, F_1204, F_2015, F_3540, F_4351, F_5432) 1134 Reference --- Adaptation of hex max aspect frobenius */ 1135 C_FUNC_DEF double v_wedge_max_aspect_frobenius( int num_nodes, double coordinates[][3] ); 1136 1137 //! Calculates wedge mean aspect forbenius. 1138 /** 1/6 * (F_0123 + F_1204 + F+2015 + F_3540 + F_4351 + F_5432) 1139 Reference --- Adaptation of hex mean aspect frobenius */ 1140 C_FUNC_DEF double v_wedge_mean_aspect_frobenius( int num_nodes, double coordinates[][3] ); 1141 1142 //! Calculates wedge jacobian metric. 1143 /** min{((L_2 X L_0) * L_3)_k} 1144 Reference --- Adaptation of Tet jacobian metric. */ 1145 C_FUNC_DEF double v_wedge_jacobian( int num_nodes, double coordinates[][3]); 1146 1147 //! Calculates wedge distortion metric. 1148 /** {min(|J|)/actual volume}*parent volume. 1149 Reference --- Adaptation of Hex distortion metric. */ 1150 C_FUNC_DEF double v_wedge_distortion( int num_nodes, double coordinates[][3] ); 1151 1152 //! Calculates the wedge stretch 1153 /** Minimum of the stretch of each quadrilateral face. 1154 Reference -- See quadrilateral stretch */ 1155 C_FUNC_DEF double v_wedge_max_stretch( int num_nodes, double coordinates[][3] ); 1156 1157 //! Calculates wedge scaled jacobian metric. 1158 /** Reference --- Adaptation of Hex and Tet scaled jacobian metric. */ 1159 C_FUNC_DEF double v_wedge_scaled_jacobian( int num_nodes, double coordinates[][3] ); 1160 1161 //! Calculates the wedge shape metric. 1162 /** 3 divided by the minimum mean ratio of the Jacobian matrix at each 1163 element corner. 1164 Reference -- Adaptaation of Hex shape metric. */ 1165 C_FUNC_DEF double v_wedge_shape( int num_nodes, double coordinates[][3] ); 1166 1167 //! Calculates wedge max aspect forbenius. 1168 /** max(F_0123, F_1204, F_2015, F_3540, F_4351, F_5432) 1169 Reference --- Adaptation of hex max aspect frobenius */ 1170 C_FUNC_DEF double v_wedge_condition( int num_nodes, double coordinates[][3] ); 1171 1172 1173/* individual quality functions for knife elements */ 1174 1175 //! Calculates knife volume. 1176 C_FUNC_DEF double v_knife_volume( int num_nodes, double coordinates[][3] ); 1177 1178 1179/* individual quality functions for edge elements */ 1180 1181 //! Calculates edge length. 1182 C_FUNC_DEF double v_edge_length( int num_nodes, double coordinates[][3] ); 1183 1184 1185/* individual quality functions for quad elements */ 1186 //! Sets average size (area) of quad, needed for v_quad_relative_size(...) 1187 C_FUNC_DEF void v_set_quad_size( double size ); 1188 1189 //! Calculates quad edge ratio 1190 /** edge ratio 1191 Reference --- P. P. Pebay, Planar Quadrangle Quality 1192 Measures, Eng. Comp., 2004, 20(2):157-173 */ 1193 C_FUNC_DEF double v_quad_edge_ratio( int num_nodes, double coordinates[][3] ); 1194 1195 //! Calculates quad maximum of edge ratio. 1196 /** Maximum edge length ratio at quad center. 1197 Reference --- J. Robinson, CRE Method of element testing and the 1198 Jacobian shape parameters, Eng. Comput., Vol 4, 1987. */ 1199 C_FUNC_DEF double v_quad_max_edge_ratio( int num_nodes, double coordinates[][3] ); 1200 1201 //! Calculates quad aspect ratio 1202 /** aspect ratio 1203 Reference --- P. P. Pebay, Planar Quadrangle Quality 1204 Measures, Eng. Comp., 2004, 20(2):157-173 */ 1205 C_FUNC_DEF double v_quad_aspect_ratio( int num_nodes, double coordinates[][3] ); 1206 1207 //! Calculates quad radius ratio 1208 /** radius ratio 1209 Reference --- P. P. Pebay, Planar Quadrangle Quality 1210 Measures, Eng. Comp., 2004, 20(2):157-173 */ 1211 C_FUNC_DEF double v_quad_radius_ratio( int num_nodes, double coordinates[][3] ); 1212 1213 //! Calculates quad average Frobenius aspect 1214 /** average Frobenius aspect 1215 Reference --- P. P. Pebay, Planar Quadrangle Quality 1216 Measures, Eng. Comp., 2004, 20(2):157-173 */ 1217 C_FUNC_DEF double v_quad_med_aspect_frobenius( int num_nodes, double coordinates[][3] ); 1218 1219 //! Calculates quad maximum Frobenius aspect 1220 /** average Frobenius aspect 1221 Reference --- P. P. Pebay, Planar Quadrangle Quality 1222 Measures, Eng. Comp., 2004, 20(2):157-173 */ 1223 C_FUNC_DEF double v_quad_max_aspect_frobenius( int num_nodes, double coordinates[][3] ); 1224 1225 //! Calculates quad skew metric. 1226 /** Maximum |cos A| where A is the angle between edges at quad center. 1227 Reference --- J. Robinson, CRE Method of element testing and the 1228 Jacobian shape parameters, Eng. Comput., Vol 4, 1987. */ 1229 C_FUNC_DEF double v_quad_skew( int num_nodes, double coordinates[][3] ); 1230 1231 //! Calculates quad taper metric. 1232 /** Maximum ratio of lengths derived from opposite edges. 1233 Reference --- J. Robinson, CRE Method of element testing and the 1234 Jacobian shape parameters, Eng. Comput., Vol 4, 1987. */ 1235 C_FUNC_DEF double v_quad_taper( int num_nodes, double coordinates[][3] ); 1236 1237 //! Calculates quad warpage metric. 1238 /** Cosine of Minimum Dihedral Angle formed by Planes Intersecting in Diagonals. 1239 Reference --- J. Robinson, CRE Method of element testing and the 1240 Jacobian shape parameters, Eng. Comput., Vol 4, 1987. */ 1241 C_FUNC_DEF double v_quad_warpage( int num_nodes, double coordinates[][3] ); 1242 1243 //! Calculates quad area. 1244 /** Jacobian at quad center. 1245 Reference --- J. Robinson, CRE Method of element testing and the 1246 Jacobian shape parameters, Eng. Comput., Vol 4, 1987. */ 1247 C_FUNC_DEF double v_quad_area( int num_nodes, double coordinates[][3] ); 1248 1249 //! Calculates quad strech metric. 1250 /** Sqrt(2) * minimum edge length / maximum diagonal length. 1251 Reference --- FIMESH code. */ 1252 C_FUNC_DEF double v_quad_stretch( int num_nodes, double coordinates[][3] ); 1253 1254 //! Calculates quad's smallest angle. 1255 /** Smallest included quad angle (degrees). 1256 Reference --- Unknown. */ 1257 C_FUNC_DEF double v_quad_minimum_angle( int num_nodes, double coordinates[][3] ); 1258 1259 //! Calculates quad's largest angle. 1260 /** Largest included quad angle (degrees). 1261 Reference --- Unknown. */ 1262 C_FUNC_DEF double v_quad_maximum_angle( int num_nodes, double coordinates[][3] ); 1263 1264 //! Calculates quad oddy metric. 1265 C_FUNC_DEF double v_quad_oddy( int num_nodes, double coordinates[][3] ); 1266 1267 //! Calculates quad condition number metric. 1268 /** Maximum condition number of the Jacobian matrix at 4 corners. 1269 Reference --- P. Knupp, Achieving Finite Element Mesh Quality via 1270 Optimization of the Jacobian Matrix Norm and Associated Quantities, 1271 Intl. J. Numer. Meth. Engng. 2000, 48:1165-1185. */ 1272 C_FUNC_DEF double v_quad_condition( int num_nodes, double coordinates[][3] ); 1273 1274 //! Calculates quad jacobian. 1275 /** Minimum pointwise volume of local map at 4 corners & center of quad. 1276 Reference --- P. Knupp, Achieving Finite Element Mesh Quality via 1277 Optimization of the Jacobian Matrix Norm and Associated Quantities, 1278 Intl. J. Numer. Meth. Engng. 2000, 48:1165-1185. */ 1279 C_FUNC_DEF double v_quad_jacobian( int num_nodes, double coordinates[][3] ); 1280 1281 //! Calculates quad scaled jacobian. 1282 /** Minimum Jacobian divided by the lengths of the 2 edge vectors. 1283 Reference --- P. Knupp, Achieving Finite Element Mesh Quality via 1284 Optimization of the Jacobian Matrix Norm and Associated Quantities, 1285 Intl. J. Numer. Meth. Engng. 2000, 48:1165-1185. */ 1286 C_FUNC_DEF double v_quad_scaled_jacobian( int num_nodes, double coordinates[][3] ); 1287 1288 //! Calculates quad shear metric. 1289 /** 2/Condition number of Jacobian Skew matrix. 1290 Reference --- P. Knupp, Algebraic Mesh Quality Metrics for 1291 Unstructured Initial Meshes, submitted for publication. */ 1292 C_FUNC_DEF double v_quad_shear( int num_nodes, double coordinates[][3] ); 1293 1294 //! Calculates quad shape metric. 1295 /** 2/Condition number of weighted Jacobian matrix. 1296 Reference --- P. Knupp, Algebraic Mesh Quality Metrics for 1297 Unstructured Initial Meshes, submitted for publication. */ 1298 C_FUNC_DEF double v_quad_shape( int num_nodes, double coordinates[][3] ); 1299 1300 //! Calculates quad relative size metric. 1301 /** Min( J, 1/J ), where J is determinant of weighted Jacobian matrix. 1302 Reference --- P. Knupp, Algebraic Mesh Quality Metrics for 1303 Unstructured Initial Meshes, submitted for publication. */ 1304 C_FUNC_DEF double v_quad_relative_size_squared( int num_nodes, double coordinates[][3] ); 1305 1306 //! Calculates quad shape-size metric. 1307 /** Product of Shape and Relative Size. 1308 Reference --- P. Knupp, Algebraic Mesh Quality Metrics for 1309 Unstructured Initial Meshes, submitted for publication. */ 1310 C_FUNC_DEF double v_quad_shape_and_size( int num_nodes, double coordinates[][3] ); 1311 1312 //! Calculates quad shear-size metric. 1313 /** Product of Shear and Relative Size. 1314 Reference --- P. Knupp, Algebraic Mesh Quality Metrics for 1315 Unstructured Initial Meshes, submitted for publication. */ 1316 C_FUNC_DEF double v_quad_shear_and_size( int num_nodes, double coordinates[][3] ); 1317 1318 //! Calculates quad distortion metric. 1319 /** {min(|J|)/actual area}*parent area, parent area = 4 for quad. 1320 Reference --- SDRC/IDEAS Simulation: Finite Element Modeling--User's Guide */ 1321 C_FUNC_DEF double v_quad_distortion( int num_nodes, double coordinates[][3] ); 1322 1323 1324/* individual quality functions for triangle elements */ 1325 1326 //! Sets average size (area) of tri, needed for v_tri_relative_size(...) 1327 C_FUNC_DEF void v_set_tri_size( double size ); 1328 1329 //! Sets function pointer to calculate triangle normal wrt surface 1330 C_FUNC_DEF void v_set_tri_normal_func( ComputeNormal func ); 1331 1332 //! Calculates triangle metric. 1333 /** edge ratio 1334 Reference --- P. P. Pebay & T. J. Baker, Analysis of Triangle Quality 1335 Measures, AMS Math. Comp., 2003, 72(244):1817-1839 */ 1336 C_FUNC_DEF double v_tri_edge_ratio( int num_nodes, double coordinates[][3] ); 1337 1338 //! Calculates triangle metric. 1339 /** aspect ratio 1340 Reference --- P. P. Pebay & T. J. Baker, Analysis of Triangle Quality 1341 Measures, AMS Math. Comp., 2003, 72(244):1817-1839 */ 1342 C_FUNC_DEF double v_tri_aspect_ratio( int num_nodes, double coordinates[][3] ); 1343 1344 //! Calculates triangle metric. 1345 /** radius ratio 1346 Reference --- P. P. Pebay & T. J. Baker, Analysis of Triangle Quality 1347 Measures, AMS Math. Comp., 2003, 72(244):1817-1839 */ 1348 C_FUNC_DEF double v_tri_radius_ratio( int num_nodes, double coordinates[][3] ); 1349 1350 //! Calculates triangle metric. 1351 /** Frobenius aspect */ 1352 C_FUNC_DEF double v_tri_aspect_frobenius( int num_nodes, double coordinates[][3] ); 1353 1354 //! Calculates triangle metric. 1355 /** Maximum included angle in triangle */ 1356 C_FUNC_DEF double v_tri_area( int num_nodes, double coordinates[][3] ); 1357 1358 //! Calculates triangle metric. 1359 /** Minimum included angle in triangle */ 1360 C_FUNC_DEF double v_tri_minimum_angle( int num_nodes, double coordinates[][3] ); 1361 1362 //! Calculates triangle metric. 1363 /** Maximum included angle in triangle */ 1364 C_FUNC_DEF double v_tri_maximum_angle( int num_nodes, double coordinates[][3] ); 1365 1366 //! Calculates triangle metric. 1367 /** Condition number of the Jacobian matrix. 1368 Reference --- P. Knupp, Achieving Finite Element Mesh Quality via 1369 Optimization of the Jacobian Matrix Norm and Associated Quantities, 1370 Intl. J. Numer. Meth. Engng. 2000, 48:1165-1185. */ 1371 C_FUNC_DEF double v_tri_condition( int num_nodes, double coordinates[][3] ); 1372 1373 //! Calculates triangle metric. 1374 /** Minimum Jacobian divided by the lengths of 2 edge vectors. 1375 Reference --- P. Knupp, Achieving Finite Element Mesh Quality via 1376 Optimization of the Jacobian Matrix Norm and Associated Quantities, 1377 Intl. J. Numer. Meth. Engng. 2000, 48:1165-1185. */ 1378 C_FUNC_DEF double v_tri_scaled_jacobian( int num_nodes, double coordinates[][3] ); 1379 1380 //! Calculates triangle metric. 1381 /** */ 1382 C_FUNC_DEF double v_tri_shear( int num_nodes, double coordinates[][3] ); 1383 1384 //! Calculates triangle metric. 1385 /** Min( J, 1/J ), where J is determinant of weighted Jacobian matrix. 1386 Reference --- P. Knupp, Algebraic Mesh Quality Metrics for 1387 Unstructured Initial Meshes, submitted for publication. */ 1388 C_FUNC_DEF double v_tri_relative_size_squared( int num_nodes, double coordinates[][3] ); 1389 1390 //! Calculates triangle metric. 1391 /** 2/Condition number of weighted Jacobian matrix. 1392 Reference --- P. Knupp, Algebraic Mesh Quality Metrics for 1393 Unstructured Initial Meshes, submitted for publication. */ 1394 C_FUNC_DEF double v_tri_shape( int num_nodes, double coordinates[][3] ); 1395 1396 //! Calculates triangle metric. 1397 /** Product of Shape and Relative Size. 1398 Reference --- P. Knupp, Algebraic Mesh Quality Metrics for 1399 Unstructured Initial Meshes, submitted for publication. */ 1400 C_FUNC_DEF double v_tri_shape_and_size( int num_nodes, double coordinates[][3] ); 1401 1402 //! Calculates triangle metric. 1403 /** {min(|J|)/actual area}*parent area, parent area = 1/2 for triangular element. 1404 Reference --- SDRC/IDEAS Simulation: Finite Element Modeling--User's Guide */ 1405 C_FUNC_DEF double v_tri_distortion( int num_nodes, double coordinates[][3] ); 1406 1407 1408 1409#endif /* __verdict_h */ 1410 1411 1412 1413