diff --git a/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict b/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict
index 791de6ec441faf0a359af527300e592a18576b3d..7d1ad76db3007982d4977face4a0b095f9697d57 100644
--- a/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict
+++ b/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict
@@ -478,7 +478,7 @@ meshQualityControls
// resulting volume after snapping (by approximation) is larger than
// minVolCollapseRatio times old volume (i.e. not collapsed to flat cell).
// If <0 : delete always.
- //minVolCollapseRatio 0.5;
+ //minVolCollapseRatio 0.1;
// Advanced
diff --git a/src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshTools.C b/src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshTools.C
index 3af85508be1304af48032919b768bf6d95a80803..4945eb0b098d73b50721c055a6c9001fa10e10ed 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshTools.C
+++ b/src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshTools.C
@@ -97,7 +97,6 @@ Foam::tmp<Foam::scalarField> Foam::polyMeshTools::faceSkewness
{
const labelList& own = mesh.faceOwner();
const labelList& nei = mesh.faceNeighbour();
- const faceList& fcs = mesh.faces();
const polyBoundaryMesh& pbm = mesh.boundaryMesh();
tmp<scalarField> tskew(new scalarField(mesh.nFaces()));
@@ -154,31 +153,16 @@ Foam::tmp<Foam::scalarField> Foam::polyMeshTools::faceSkewness
{
label faceI = pp.start() + i;
- vector Cpf = fCtrs[faceI] - cellCtrs[own[faceI]];
-
- vector normal = fAreas[faceI];
- normal /= mag(normal) + VSMALL;
- vector d = normal*(normal & Cpf);
-
-
- // Skewness vector
- vector sv =
- Cpf
- - ((fAreas[faceI] & Cpf)/((fAreas[faceI] & d) + VSMALL))*d;
- vector svHat = sv/(mag(sv) + VSMALL);
-
- // Normalisation distance calculated as the approximate distance
- // from the face centre to the edge of the face in the direction
- // of the skewness
- scalar fd = 0.4*mag(d) + VSMALL;
- const face& f = fcs[faceI];
- forAll(f, pi)
- {
- fd = max(fd, mag(svHat & (p[f[pi]] - fCtrs[faceI])));
- }
+ skew[faceI] = primitiveMeshTools::boundaryFaceSkewness
+ (
+ mesh,
+ p,
+ fCtrs,
+ fAreas,
- // Normalised skewness
- skew[faceI] = mag(sv)/fd;
+ faceI,
+ cellCtrs[own[faceI]]
+ );
}
}
}
diff --git a/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCheck/primitiveMeshTools.C b/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCheck/primitiveMeshTools.C
index e94c8edbf77a1709d63377a11c6fb98e5765e067..241f4eef8ed91e131a37c6e5ae9288ac1087e8a1 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCheck/primitiveMeshTools.C
+++ b/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCheck/primitiveMeshTools.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2012-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -64,6 +64,44 @@ Foam::scalar Foam::primitiveMeshTools::faceSkewness
return mag(sv)/fd;
}
+Foam::scalar Foam::primitiveMeshTools::boundaryFaceSkewness
+(
+ const primitiveMesh& mesh,
+ const pointField& p,
+ const vectorField& fCtrs,
+ const vectorField& fAreas,
+
+ const label faceI,
+ const point& ownCc
+)
+{
+ vector Cpf = fCtrs[faceI] - ownCc;
+
+ vector normal = fAreas[faceI];
+ normal /= mag(normal) + VSMALL;
+ vector d = normal*(normal & Cpf);
+
+
+ // Skewness vector
+ vector sv =
+ Cpf
+ - ((fAreas[faceI] & Cpf)/((fAreas[faceI] & d) + VSMALL))*d;
+ vector svHat = sv/(mag(sv) + VSMALL);
+
+ // Normalisation distance calculated as the approximate distance
+ // from the face centre to the edge of the face in the direction
+ // of the skewness
+ scalar fd = 0.4*mag(d) + VSMALL;
+ const face& f = mesh.faces()[faceI];
+ forAll(f, pi)
+ {
+ fd = max(fd, mag(svHat & (p[f[pi]] - fCtrs[faceI])));
+ }
+
+ // Normalised skewness
+ return mag(sv)/fd;
+}
+
Foam::scalar Foam::primitiveMeshTools::faceOrthogonality
(
@@ -119,7 +157,6 @@ Foam::tmp<Foam::scalarField> Foam::primitiveMeshTools::faceSkewness
{
const labelList& own = mesh.faceOwner();
const labelList& nei = mesh.faceNeighbour();
- const faceList& fcs = mesh.faces();
tmp<scalarField> tskew(new scalarField(mesh.nFaces()));
scalarField& skew = tskew();
@@ -145,31 +182,15 @@ Foam::tmp<Foam::scalarField> Foam::primitiveMeshTools::faceSkewness
for (label faceI = mesh.nInternalFaces(); faceI < mesh.nFaces(); faceI++)
{
- vector Cpf = fCtrs[faceI] - cellCtrs[own[faceI]];
-
- vector normal = fAreas[faceI];
- normal /= mag(normal) + VSMALL;
- vector d = normal*(normal & Cpf);
-
-
- // Skewness vector
- vector sv =
- Cpf
- - ((fAreas[faceI] & Cpf)/((fAreas[faceI] & d) + VSMALL))*d;
- vector svHat = sv/(mag(sv) + VSMALL);
-
- // Normalisation distance calculated as the approximate distance
- // from the face centre to the edge of the face in the direction
- // of the skewness
- scalar fd = 0.4*mag(d) + VSMALL;
- const face& f = fcs[faceI];
- forAll(f, pi)
- {
- fd = max(fd, mag(svHat & (p[f[pi]] - fCtrs[faceI])));
- }
-
- // Normalised skewness
- skew[faceI] = mag(sv)/fd;
+ skew[faceI] = boundaryFaceSkewness
+ (
+ mesh,
+ p,
+ fCtrs,
+ fAreas,
+ faceI,
+ cellCtrs[own[faceI]]
+ );
}
return tskew;
diff --git a/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCheck/primitiveMeshTools.H b/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCheck/primitiveMeshTools.H
index 5e3b2c07fd4b03aaacab2f82fca09c83187e74c9..170d089e117a9be1831d5b0c00733919b0102a2d 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCheck/primitiveMeshTools.H
+++ b/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCheck/primitiveMeshTools.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2012-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -48,32 +48,6 @@ namespace Foam
class primitiveMeshTools
{
-
-protected:
-
- // Protected Member Functions
-
- //- Skewness of single face
- static scalar faceSkewness
- (
- const primitiveMesh& mesh,
- const pointField& p,
- const vectorField& fCtrs,
- const vectorField& fAreas,
-
- const label faceI,
- const point& ownCc,
- const point& neiCc
- );
-
- //- Orthogonality of single face
- static scalar faceOrthogonality
- (
- const point& ownCc,
- const point& neiCc,
- const vector& s
- );
-
public:
//- Generate non-orthogonality field (internal faces only)
@@ -154,6 +128,41 @@ public:
const PackedBoolList& internalOrCoupledFace
);
+
+ // Helpers: single face check
+
+ //- Skewness of single face
+ static scalar faceSkewness
+ (
+ const primitiveMesh& mesh,
+ const pointField& p,
+ const vectorField& fCtrs,
+ const vectorField& fAreas,
+
+ const label faceI,
+ const point& ownCc,
+ const point& neiCc
+ );
+
+ //- Skewness of single boundary face
+ static scalar boundaryFaceSkewness
+ (
+ const primitiveMesh& mesh,
+ const pointField& p,
+ const vectorField& fCtrs,
+ const vectorField& fAreas,
+
+ const label faceI,
+ const point& ownCc
+ );
+
+ //- Orthogonality of single face
+ static scalar faceOrthogonality
+ (
+ const point& ownCc,
+ const point& neiCc,
+ const vector& s
+ );
};
diff --git a/src/dynamicMesh/motionSmoother/motionSmootherCheck.C b/src/dynamicMesh/motionSmoother/motionSmootherCheck.C
index 77ce89ca71e2a7d238d973b5af6c13a766a59b37..3efb62a995ae9f66cb46e096950f17a0869c9957 100644
--- a/src/dynamicMesh/motionSmoother/motionSmootherCheck.C
+++ b/src/dynamicMesh/motionSmoother/motionSmootherCheck.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -239,6 +239,7 @@ bool Foam::motionSmoother::checkMesh
maxIntSkew,
maxBounSkew,
mesh,
+ mesh.points(),
mesh.cellCentres(),
mesh.faceCentres(),
mesh.faceAreas(),
@@ -481,14 +482,14 @@ bool Foam::motionSmoother::checkMesh
(
readScalar(dict.lookup("minArea", true))
);
- const scalar maxIntSkew
- (
- readScalar(dict.lookup("maxInternalSkewness", true))
- );
- const scalar maxBounSkew
- (
- readScalar(dict.lookup("maxBoundarySkewness", true))
- );
+ //const scalar maxIntSkew
+ //(
+ // readScalar(dict.lookup("maxInternalSkewness", true))
+ //);
+ //const scalar maxBounSkew
+ //(
+ // readScalar(dict.lookup("maxBoundarySkewness", true))
+ //);
const scalar minWeight
(
readScalar(dict.lookup("minFaceWeight", true))
@@ -615,27 +616,30 @@ bool Foam::motionSmoother::checkMesh
nWrongFaces = nNewWrongFaces;
}
- if (maxIntSkew > 0 || maxBounSkew > 0)
- {
- meshGeom.checkFaceSkewness
- (
- report,
- maxIntSkew,
- maxBounSkew,
- checkFaces,
- baffles,
- &wrongFaces
- );
-
- label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
-
- Info<< " faces with skewness > "
- << setw(3) << maxIntSkew
- << " (internal) or " << setw(3) << maxBounSkew
- << " (boundary) : " << nNewWrongFaces-nWrongFaces << endl;
- nWrongFaces = nNewWrongFaces;
- }
+ //- Note: cannot check the skewness without the points and don't want
+ // to store them on polyMeshGeometry.
+ //if (maxIntSkew > 0 || maxBounSkew > 0)
+ //{
+ // meshGeom.checkFaceSkewness
+ // (
+ // report,
+ // maxIntSkew,
+ // maxBounSkew,
+ // checkFaces,
+ // baffles,
+ // &wrongFaces
+ // );
+ //
+ // label nNewWrongFaces = returnReduce(wrongFaces.size(),sumOp<label>());
+ //
+ // Info<< " faces with skewness > "
+ // << setw(3) << maxIntSkew
+ // << " (internal) or " << setw(3) << maxBounSkew
+ // << " (boundary) : " << nNewWrongFaces-nWrongFaces << endl;
+ //
+ // nWrongFaces = nNewWrongFaces;
+ //}
if (minWeight >= 0 && minWeight < 1)
{
diff --git a/src/dynamicMesh/motionSmoother/polyMeshGeometry/polyMeshGeometry.C b/src/dynamicMesh/motionSmoother/polyMeshGeometry/polyMeshGeometry.C
index d3b69d392f391a8d554d3557d677b006a3d8902e..57725dc57c95aaf24cabc25133b7a70b9b110afb 100644
--- a/src/dynamicMesh/motionSmoother/polyMeshGeometry/polyMeshGeometry.C
+++ b/src/dynamicMesh/motionSmoother/polyMeshGeometry/polyMeshGeometry.C
@@ -29,6 +29,7 @@ License
#include "tetrahedron.H"
#include "syncTools.H"
#include "unitConversion.H"
+#include "primitiveMeshTools.H"
namespace Foam
{
@@ -286,29 +287,6 @@ Foam::scalar Foam::polyMeshGeometry::checkNonOrtho
}
-Foam::scalar Foam::polyMeshGeometry::calcSkewness
-(
- const point& ownCc,
- const point& neiCc,
- const point& fc
-)
-{
- scalar dOwn = mag(fc - ownCc);
- scalar dNei = mag(fc - neiCc);
-
- point faceIntersection =
- ownCc*dNei/(dOwn+dNei+VSMALL)
- + neiCc*dOwn/(dOwn+dNei+VSMALL);
-
- return
- mag(fc - faceIntersection)
- / (
- mag(neiCc-ownCc)
- + VSMALL
- );
-}
-
-
// Create the neighbour pyramid - it will have positive volume
bool Foam::polyMeshGeometry::checkFaceTet
(
@@ -1008,6 +986,7 @@ bool Foam::polyMeshGeometry::checkFaceSkewness
const scalar internalSkew,
const scalar boundarySkew,
const polyMesh& mesh,
+ const pointField& points,
const vectorField& cellCentres,
const vectorField& faceCentres,
const vectorField& faceAreas,
@@ -1024,14 +1003,8 @@ bool Foam::polyMeshGeometry::checkFaceSkewness
const polyBoundaryMesh& patches = mesh.boundaryMesh();
// Calculate coupled cell centre
- pointField neiCc(mesh.nFaces()-mesh.nInternalFaces());
-
- for (label faceI = mesh.nInternalFaces(); faceI < mesh.nFaces(); faceI++)
- {
- neiCc[faceI-mesh.nInternalFaces()] = cellCentres[own[faceI]];
- }
- syncTools::swapBoundaryFacePositions(mesh, neiCc);
-
+ pointField neiCc;
+ syncTools::swapBoundaryCellPositions(mesh, cellCentres, neiCc);
scalar maxSkew = 0;
@@ -1043,11 +1016,16 @@ bool Foam::polyMeshGeometry::checkFaceSkewness
if (mesh.isInternalFace(faceI))
{
- scalar skewness = calcSkewness
+ scalar skewness = primitiveMeshTools::faceSkewness
(
+ mesh,
+ points,
+ faceCentres,
+ faceAreas,
+
+ faceI,
cellCentres[own[faceI]],
- cellCentres[nei[faceI]],
- faceCentres[faceI]
+ cellCentres[nei[faceI]]
);
// Check if the skewness vector is greater than the PN vector.
@@ -1073,11 +1051,16 @@ bool Foam::polyMeshGeometry::checkFaceSkewness
}
else if (patches[patches.whichPatch(faceI)].coupled())
{
- scalar skewness = calcSkewness
+ scalar skewness = primitiveMeshTools::faceSkewness
(
+ mesh,
+ points,
+ faceCentres,
+ faceAreas,
+
+ faceI,
cellCentres[own[faceI]],
- neiCc[faceI-mesh.nInternalFaces()],
- faceCentres[faceI]
+ neiCc[faceI-mesh.nInternalFaces()]
);
// Check if the skewness vector is greater than the PN vector.
@@ -1103,21 +1086,17 @@ bool Foam::polyMeshGeometry::checkFaceSkewness
}
else
{
- // Boundary faces: consider them to have only skewness error.
- // (i.e. treat as if mirror cell on other side)
-
- vector faceNormal = faceAreas[faceI];
- faceNormal /= mag(faceNormal) + ROOTVSMALL;
-
- vector dOwn = faceCentres[faceI] - cellCentres[own[faceI]];
-
- vector dWall = faceNormal*(faceNormal & dOwn);
+ scalar skewness = primitiveMeshTools::boundaryFaceSkewness
+ (
+ mesh,
+ points,
+ faceCentres,
+ faceAreas,
- point faceIntersection = cellCentres[own[faceI]] + dWall;
+ faceI,
+ cellCentres[own[faceI]]
+ );
- scalar skewness =
- mag(faceCentres[faceI] - faceIntersection)
- /(2*mag(dWall) + ROOTVSMALL);
// Check if the skewness vector is greater than the PN vector.
// This does not cause trouble but is a good indication of a poor
@@ -1148,12 +1127,18 @@ bool Foam::polyMeshGeometry::checkFaceSkewness
label face1 = baffles[i].second();
const point& ownCc = cellCentres[own[face0]];
+ const point& neiCc = cellCentres[own[face1]];
- scalar skewness = calcSkewness
+ scalar skewness = primitiveMeshTools::faceSkewness
(
+ mesh,
+ points,
+ faceCentres,
+ faceAreas,
+
+ face0,
ownCc,
- cellCentres[own[face1]],
- faceCentres[face0]
+ neiCc
);
// Check if the skewness vector is greater than the PN vector.
@@ -2361,30 +2346,30 @@ bool Foam::polyMeshGeometry::checkFaceTets
}
-bool Foam::polyMeshGeometry::checkFaceSkewness
-(
- const bool report,
- const scalar internalSkew,
- const scalar boundarySkew,
- const labelList& checkFaces,
- const List<labelPair>& baffles,
- labelHashSet* setPtr
-) const
-{
- return checkFaceSkewness
- (
- report,
- internalSkew,
- boundarySkew,
- mesh_,
- cellCentres_,
- faceCentres_,
- faceAreas_,
- checkFaces,
- baffles,
- setPtr
- );
-}
+//bool Foam::polyMeshGeometry::checkFaceSkewness
+//(
+// const bool report,
+// const scalar internalSkew,
+// const scalar boundarySkew,
+// const labelList& checkFaces,
+// const List<labelPair>& baffles,
+// labelHashSet* setPtr
+//) const
+//{
+// return checkFaceSkewness
+// (
+// report,
+// internalSkew,
+// boundarySkew,
+// mesh_,
+// cellCentres_,
+// faceCentres_,
+// faceAreas_,
+// checkFaces,
+// baffles,
+// setPtr
+// );
+//}
bool Foam::polyMeshGeometry::checkFaceWeights
diff --git a/src/dynamicMesh/motionSmoother/polyMeshGeometry/polyMeshGeometry.H b/src/dynamicMesh/motionSmoother/polyMeshGeometry/polyMeshGeometry.H
index 2ef8604cac8714a45711323d82f2cbc4f475d574..405564a175273efd8d273725ef121558bc816f60 100644
--- a/src/dynamicMesh/motionSmoother/polyMeshGeometry/polyMeshGeometry.H
+++ b/src/dynamicMesh/motionSmoother/polyMeshGeometry/polyMeshGeometry.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -229,6 +229,7 @@ public:
const scalar internalSkew,
const scalar boundarySkew,
const polyMesh& mesh,
+ const pointField& points,
const vectorField& cellCentres,
const vectorField& faceCentres,
const vectorField& faceAreas,
@@ -372,16 +373,6 @@ public:
labelHashSet* setPtr
) const;
- bool checkFaceSkewness
- (
- const bool report,
- const scalar internalSkew,
- const scalar boundarySkew,
- const labelList& checkFaces,
- const List<labelPair>& baffles,
- labelHashSet* setPtr
- ) const;
-
bool checkFaceWeights
(
const bool report,
diff --git a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriverShrink.C b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriverShrink.C
index 3a7c4fb2496a98fcb0befd6c099dc9961c51b6ae..bed0e3c52fcd03d6e64fe3d93052b46984ca894f 100644
--- a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriverShrink.C
+++ b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriverShrink.C
@@ -336,7 +336,7 @@ void Foam::autoLayerDriver::smoothNormals
) const
{
// Get smoothly varying internal normals field.
- Info<< "shrinkMeshDistance : Smoothing normals ..." << endl;
+ Info<< "shrinkMeshDistance : Smoothing normals in interior ..." << endl;
const fvMesh& mesh = meshRefiner_.mesh();
const edgeList& edges = mesh.edges();
@@ -1161,6 +1161,11 @@ void Foam::autoLayerDriver::medialAxisSmoothingInfo
scalar mDist = medialDist[pointI];
if (wDist2 < sqr(SMALL) && mDist < SMALL)
+ //- Note: maybe less strict:
+ //(
+ // wDist2 < sqr(meshRefiner_.mergeDistance())
+ // && mDist < meshRefiner_.mergeDistance()
+ //)
{
medialRatio[pointI] = 0.0;
}