src/dynamicMesh/dynamicTopoFvMesh/tetMetrics/tetMetrics.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | foam-extend: Open Source CFD
   4    \\    /   O peration     | Version:     3.2
   5     \\  /    A nd           | Web:         http://www.foam-extend.org
   6      \\/     M anipulation  | For copyright notice see file Copyright
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of foam-extend.
  10
  11     foam-extend is free software: you can redistribute it and/or modify it
  12     under the terms of the GNU General Public License as published by the
  13     Free Software Foundation, either version 3 of the License, or (at your
  14     option) any later version.
  15
  16     foam-extend is distributed in the hope that it will be useful, but
  17     WITHOUT ANY WARRANTY; without even the implied warranty of
  18     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19     General Public License for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
  23
  24 Class
  25     tetMetrics
  26
  27 Description
  28     Implementation of tetrahedral mesh-quality metrics.
  29
  30     Most measures are referenced from:
  31     Jonathan R. Shewchuk. What is a good linear finite element?
  32     Interpolation, Conditioning, Anisotropy and Quality Measures.
  33     Eleventh International Meshing Roundtable (Ithaca, NY), pp 115-126, 2002.
  34
  35     NOTE:
  36     For metrics involving volume the function assumes points (p0-p1-p2)
  37     are in counter-clockwise fashion when viewed from p3.
  38
  39 Author
  40     Sandeep Menon
  41     University of Massachusetts Amherst
  42     All rights reserved
  43
  44 \*----------------------------------------------------------------------------*/
  45
  46 #include "tetMetrics.H"
  47 #include "addToMemberFunctionSelectionTable.H"
  48
  49 namespace Foam
  50 {
  51
  52 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  53
  54 defineTypeNameAndDebug(Knupp,0);
  55 defineTypeNameAndDebug(Dihedral,0);
  56 defineTypeNameAndDebug(cubicMeanRatio,0);
  57 defineTypeNameAndDebug(Frobenius,0);
  58 defineTypeNameAndDebug(PGH,0);
  59 defineTypeNameAndDebug(CSG,0);
  60
  61
  62 addToMemberFunctionSelectionTable(tetMetric, Knupp, metric, Point);
  63 addToMemberFunctionSelectionTable(tetMetric, Dihedral, metric, Point);
  64 addToMemberFunctionSelectionTable(tetMetric, cubicMeanRatio, metric, Point);
  65 addToMemberFunctionSelectionTable(tetMetric, Frobenius, metric, Point);
  66 addToMemberFunctionSelectionTable(tetMetric, PGH, metric, Point);
  67 addToMemberFunctionSelectionTable(tetMetric, CSG, metric, Point);
  68
  69
  70 // Enumeration for tets
  71 label Dihedral::tetEnum[6][4] =
  72 {
  73     {0,1,2,3},
  74     {0,2,3,1},
  75     {0,3,1,2},
  76     {1,2,0,3},
  77     {1,3,0,2},
  78     {2,3,0,1}
  79 };
  80
  81 // * * * * * * * * * * * * * Static Members Functions * * * * * * * * * * *  //
  82
  83 // Tetrahedral mesh-quality metric suggested by Knupp [2003].
  84 scalar Knupp::metric
  85 (
  86     const point& p0,
  87     const point& p1,
  88     const point& p2,
  89     const point& p3
  90 )
  91 {
  92     // Obtain signed tet volume
  93     scalar V = (1.0/6.0)*(((p1 - p0) ^ (p2 - p0)) & (p3 - p0));
  94
  95     // Obtain the magSqr edge-lengths
  96     scalar Le = ((p1-p0) & (p1-p0))
  97               + ((p2-p0) & (p2-p0))
  98               + ((p3-p0) & (p3-p0))
  99               + ((p2-p1) & (p2-p1))
 100               + ((p3-p1) & (p3-p1))
 101               + ((p3-p2) & (p3-p2));
 102
 103     // Return signed quality
 104     return sign(V)*((24.96100588*::cbrt(V*V))/Le);
 105 }
 106
 107
 108 // Minimum dihedral angle among six edges of the tetrahedron. Normalized
 109 // by 70.529 degrees (equilateral tet) and signed by volume.
 110 scalar Dihedral::metric
 111 (
 112     const point& p0,
 113     const point& p1,
 114     const point& p2,
 115     const point& p3
 116 )
 117 {
 118     scalar minAngle = 0.0;
 119     FixedList<scalar,6> cosAngles(1.0);
 120     FixedList<vector,4> pts(vector::zero);
 121
 122     // Assign point-positions
 123     pts[0] = p0;
 124     pts[1] = p1;
 125     pts[2] = p2;
 126     pts[3] = p3;
 127
 128     // Permute over all six edges
 129     for (label i = 0; i < 6; i++)
 130     {
 131         // Normalize the axis
 132         vector v0 = (pts[tetEnum[i][1]] - pts[tetEnum[i][0]])
 133                     /mag(pts[tetEnum[i][1]] - pts[tetEnum[i][0]]);
 134
 135         // Obtain plane-vectors
 136         vector v1 = (pts[tetEnum[i][2]] - pts[tetEnum[i][0]]);
 137         vector v2 = (pts[tetEnum[i][3]] - pts[tetEnum[i][0]]);
 138
 139         v1 -= (v1 & v0)*v0;
 140         v2 -= (v2 & v0)*v0;
 141
 142         cosAngles[i] = acos((v1/mag(v1)) & (v2/mag(v2)));
 143
 144         // Compute minimum angle on-the-fly
 145         if (i > 0)
 146         {
 147             minAngle = minAngle < cosAngles[i] ? minAngle : cosAngles[i];
 148         }
 149         else
 150         {
 151             minAngle = cosAngles[0];
 152         }
 153     }
 154
 155     // Compute signed volume and multiply by the normalized angle
 156     return sign(((p1 - p0) ^ (p2 - p0)) & (p3 - p0))*(minAngle/1.2309632);
 157 }
 158
 159
 160 // Cubic Mean Ratio Tetrahedral mesh metric
 161 // Liu,A. and Joe, B., “On the shape of tetrahedra from bisection”
 162 // Mathematics of Computation, Vol. 63, 1994, pp. 141–154.
 163 scalar cubicMeanRatio::metric
 164 (
 165     const point& p0,
 166     const point& p1,
 167     const point& p2,
 168     const point& p3
 169 )
 170 {
 171     // Obtain signed tet volume
 172     scalar V = (1.0/6.0)*(((p1 - p0) ^ (p2 - p0)) & (p3 - p0));
 173
 174     // Obtain the magSqr edge-lengths
 175     scalar Le = ((p1-p0) & (p1-p0))
 176               + ((p2-p0) & (p2-p0))
 177               + ((p3-p0) & (p3-p0))
 178               + ((p2-p1) & (p2-p1))
 179               + ((p3-p1) & (p3-p1))
 180               + ((p3-p2) & (p3-p2));
 181
 182     // Return signed quality
 183     return  sign(V)*((15552.0*V*V)/(Le*Le*Le));
 184 }
 185
 186
 187 // Tetrahedral mesh-metric based on the Frobenius Condition Number
 188 // Patrick M. Knupp. Matrix Norms & the Condition Number: A General Framework
 189 // to Improve Mesh Quality via Node-Movement. Eighth International Meshing
 190 // Roundtable (Lake Tahoe, California), pages 13–22, October 1999.
 191 scalar Frobenius::metric
 192 (
 193     const point& p0,
 194     const point& p1,
 195     const point& p2,
 196     const point& p3
 197 )
 198 {
 199     // Obtain signed tet volume
 200     scalar V = (1.0/6.0)*(((p1 - p0) ^ (p2 - p0)) & (p3 - p0));
 201
 202     // Obtain the magSqr edge-lengths
 203     scalar Le = ((p1-p0) & (p1-p0))
 204               + ((p2-p0) & (p2-p0))
 205               + ((p3-p0) & (p3-p0))
 206               + ((p2-p1) & (p2-p1))
 207               + ((p3-p1) & (p3-p1))
 208               + ((p3-p2) & (p3-p2));
 209
 210     // Compute magSqr of face-areas
 211     scalar A = magSqr(0.5*((p1-p0) ^ (p2-p0)))
 212              + magSqr(0.5*((p1-p0) ^ (p3-p0)))
 213              + magSqr(0.5*((p2-p0) ^ (p3-p0)))
 214              + magSqr(0.5*((p3-p1) ^ (p2-p1)));
 215
 216     // Return signed quality
 217     return 3.67423461*(V/sqrt((Le/6.0)*(A/4.0)));
 218 }
 219
 220
 221 // Tetrahedral mesh-metric suggested by:
 222 // V. N. Parthasarathy, C. M. Graichen, and A. F. Hathaway.
 223 // Fast Evaluation & Improvement of Tetrahedral 3-D Grid Quality. [1991]
 224 scalar PGH::metric
 225 (
 226     const point& p0,
 227     const point& p1,
 228     const point& p2,
 229     const point& p3
 230 )
 231 {
 232     // Obtain signed tet volume
 233     scalar V = (1.0/6.0)*(((p1 - p0) ^ (p2 - p0)) & (p3 - p0));
 234
 235     // Obtain the magSqr edge-lengths
 236     scalar Le = ((p1-p0) & (p1-p0))
 237               + ((p2-p0) & (p2-p0))
 238               + ((p3-p0) & (p3-p0))
 239               + ((p2-p1) & (p2-p1))
 240               + ((p3-p1) & (p3-p1))
 241               + ((p3-p2) & (p3-p2));
 242
 243     // Return signed quality
 244     return 8.48528137*(V/pow(Le/4.0,1.5));
 245 }
 246
 247
 248 // Metric suggested by:
 249 // Hugues L. de Cougny, Mark S. Shephard, and Marcel K. Georges.
 250 // Explicit Node Point Smoothing Within Octree. Technical Report 10-1990,
 251 // Scientific Computation Research Center, Rensselaer Polytechnic Institute,
 252 // Troy, New York. [1990]
 253 scalar CSG::metric
 254 (
 255     const point& p0,
 256     const point& p1,
 257     const point& p2,
 258     const point& p3
 259 )
 260 {
 261     // Obtain signed tet volume
 262     scalar V = (1.0/6.0)*(((p1 - p0) ^ (p2 - p0)) & (p3 - p0));
 263
 264     // Compute magSqr of face-areas
 265     scalar A = magSqr(0.5*((p1-p0) ^ (p2-p0)))
 266              + magSqr(0.5*((p1-p0) ^ (p3-p0)))
 267              + magSqr(0.5*((p2-p0) ^ (p3-p0)))
 268              + magSqr(0.5*((p3-p1) ^ (p2-p1)));
 269
 270     // Return signed quality
 271     return 6.83852117*(V/pow(A,0.75));
 272 }
 273
 274
 275 } // End namespace Foam
 276
 277 // ************************************************************************* //