diff --git a/applications/utilities/mesh/generation/CV3DMesher/CV3D.C b/applications/utilities/mesh/generation/CV3DMesher/CV3D.C index aa0c7e3e09..d4c449fbb0 100644 --- a/applications/utilities/mesh/generation/CV3DMesher/CV3D.C +++ b/applications/utilities/mesh/generation/CV3DMesher/CV3D.C @@ -65,6 +65,198 @@ void Foam::CV3D::reinsertPoints(const pointField& points) } +void Foam::CV3D::setVertexAlignmentDirections() +{ + for + ( + Triangulation::Finite_vertices_iterator vit = finite_vertices_begin(); + vit != finite_vertices_end(); + vit++ + ) + { + if (vit->internalOrBoundaryPoint()) + { + List alignmentDirections(vit->alignmentDirections()); + + point vert(topoint(vit->point())); + + pointIndexHit pHit = qSurf_.tree().findNearest + ( + vert, + controls_.nearWallAlignedDist2 + ); + + if (pHit.hit()) + { + // Primary alignment + vector np = qSurf_.faceNormals()[pHit.index()]; + + // Generate equally spaced 'spokes' in a circle normal to the + // direction from the vertex to the closest point on the surface + // and look for a secondary intersection. + + vector d = pHit.hitPoint() - vert; + + tensor R = rotationTensor(vector(0,0,1), np); + + label s = 36; + + scalar closestSpokeHitDistance = GREAT; + + point closestSpokeHitPoint = point(GREAT,GREAT,GREAT); + + label closestSpokeHitDistanceIndex = -1; + + for(label i = 0; i < s; i++) + { + vector spoke + ( + Foam::cos(i*mathematicalConstant::twoPi/s), + Foam::sin(i*mathematicalConstant::twoPi/s), + 0 + ); + + spoke *= controls_.nearWallAlignedDist; + + spoke = R & spoke; + + pointIndexHit spokeHit; + + // internal spoke + + spokeHit = qSurf_.tree().findLine + ( + vert, + vert + spoke + ); + + if (spokeHit.hit()) + { + scalar spokeHitDistance = mag + ( + spokeHit.hitPoint() - vert + ); + + if (spokeHitDistance < closestSpokeHitDistance) + { + closestSpokeHitPoint = spokeHit.hitPoint(); + + closestSpokeHitDistanceIndex = spokeHit.index(); + } + } + + //external spoke + + point mirrorVert = vert + 2*d; + + spokeHit = qSurf_.tree().findLine + ( + mirrorVert, + mirrorVert + spoke + ); + + if (spokeHit.hit()) + { + scalar spokeHitDistance = mag + ( + spokeHit.hitPoint() - mirrorVert + ); + + if (spokeHitDistance < closestSpokeHitDistance) + { + closestSpokeHitPoint = spokeHit.hitPoint(); + + closestSpokeHitDistanceIndex = spokeHit.index(); + } + } + } + + if (closestSpokeHitDistanceIndex > -1) + { + // Auxiliary alignment generated by spoke intersection + // normal. + + vector na = + qSurf_.faceNormals()[closestSpokeHitDistanceIndex]; + + // Secondary alignment + vector ns = np ^ na; + + if (mag(ns) < SMALL) + { + FatalErrorIn("Foam::CV3D::setVertexAlignmentDirections") + << "Parallel normals detected in spoke search." + << nl << exit(FatalError); + } + + ns /= mag(ns); + + // Tertiary alignment + vector nt = ns ^ np; + + alignmentDirections.setSize(3); + + alignmentDirections[0] = np; + + alignmentDirections[1] = ns; + + alignmentDirections[2] = nt; + + Info<< "internal " << vit->internalPoint() + << nl << alignmentDirections + << nl << "v " << vert + alignmentDirections[0] + << nl << "v " << vert + alignmentDirections[1] + << nl << "v " << vert + alignmentDirections[2] + << nl << "v " << vert + << nl << "v " << pHit.hitPoint() + << nl << "v " << closestSpokeHitPoint + << nl << "f 4 1" + << nl << "f 4 2" + << nl << "f 4 3" + << nl << endl; + } + else + { + // Using on the primary alignment and the closest of the + // globalAlignmentDirections + + alignmentDirections.setSize(0); + } + + + } + else + { + alignmentDirections.setSize(0); + } + } + } +} + + +Foam::scalar Foam::CV3D::alignmentDistanceWeight(scalar dist) const +{ + scalar w; + + scalar x = dist/controls_.nearWallAlignedDist; + + if (x < 0.5) + { + w = 1; + } + + else if (x < 1) + { + w = 2*(1 - x); + } + else + { + w = 0; + } + + return w; +} + // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // Foam::CV3D::CV3D @@ -208,36 +400,6 @@ void Foam::CV3D::insertGrid() Random rndGen(1321); scalar pert = controls_.randomPerturbation*cmptMin(delta); - // for (int i=0; iindex() = nVert++; - // } - // } - // } - // } - - // Info<< nVert - startOfInternalPoints_ << " vertices inserted" << endl; - std::vector initialPoints; for (int i=0; iindex() = nVert++; } } } @@ -341,12 +504,10 @@ void Foam::CV3D::relaxPoints(const scalar relaxation) dualVerti++; } - // else - // { - // cit->cellIndex() = -1; - // } } + setVertexAlignmentDirections(); + dualVertices.setSize(dualVerti); // loop around the Delaunay edges to construct the dual faces. diff --git a/applications/utilities/mesh/generation/CV3DMesher/CV3D.H b/applications/utilities/mesh/generation/CV3DMesher/CV3D.H index e711b17b25..d1d4defd2e 100644 --- a/applications/utilities/mesh/generation/CV3DMesher/CV3D.H +++ b/applications/utilities/mesh/generation/CV3DMesher/CV3D.H @@ -113,6 +113,9 @@ public: //- Square of nearWallAlignedDist scalar nearWallAlignedDist2; + //- Alignment direction vectors + vectorField globalAlignmentDirections; + //- Chose if the cell orientation should relax during the iterations // or remain fixed to the x-y directions Switch relaxOrientation; @@ -331,6 +334,10 @@ private: labelList& patchStarts ); + void setVertexAlignmentDirections(); + + scalar alignmentDistanceWeight(scalar dist) const; + public: diff --git a/applications/utilities/mesh/generation/CV3DMesher/CV3DIO.C b/applications/utilities/mesh/generation/CV3DMesher/CV3DIO.C index 7199b4a424..450d9d143b 100644 --- a/applications/utilities/mesh/generation/CV3DMesher/CV3DIO.C +++ b/applications/utilities/mesh/generation/CV3DMesher/CV3DIO.C @@ -166,6 +166,7 @@ void Foam::CV3D::writeMesh(bool writeToTimestep) patchStarts ); + writeDual(points, faces, "dualMesh.obj"); IOobject io ( @@ -198,11 +199,19 @@ void Foam::CV3D::writeMesh(bool writeToTimestep) polyMesh pMesh ( io, - points, - faces, - owner, - neighbour + xferMove(points), + xferMove(faces), + xferMove(owner), + xferMove(neighbour) ); + // polyMesh pMesh + // ( + // io, + // points, + // faces, + // owner, + // neighbour + // ); List patches(patchStarts.size()); @@ -227,7 +236,7 @@ void Foam::CV3D::writeMesh(bool writeToTimestep) << exit(FatalError); } - writeDual(points, faces, "dualMesh.obj"); + } // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * // diff --git a/applications/utilities/mesh/generation/CV3DMesher/CV3DMesher.C b/applications/utilities/mesh/generation/CV3DMesher/CV3DMesher.C index 0139f6cd17..60a60bb545 100644 --- a/applications/utilities/mesh/generation/CV3DMesher/CV3DMesher.C +++ b/applications/utilities/mesh/generation/CV3DMesher/CV3DMesher.C @@ -73,7 +73,7 @@ int main(int argc, char *argv[]) << "Read surface with " << surf.size() << " triangles from file " << args.args()[1] << nl << endl; - // Read and triangulation + // Read and triangulation // ~~~~~~~~~~~~~~~~~~~~~~ CV3D mesh(runTime, surf); diff --git a/applications/utilities/mesh/generation/CV3DMesher/controls.C b/applications/utilities/mesh/generation/CV3DMesher/controls.C index 6675af0dfe..6bf7d89e18 100644 --- a/applications/utilities/mesh/generation/CV3DMesher/controls.C +++ b/applications/utilities/mesh/generation/CV3DMesher/controls.C @@ -45,6 +45,7 @@ Foam::CV3D::controls::controls(const IOdictionary& controlDict) readScalar(controlDict.lookup("nearWallAlignedDist"))*minCellSize ), nearWallAlignedDist2(Foam::sqr(nearWallAlignedDist)), + globalAlignmentDirections(controlDict.lookup("globalAlignmentDirections")), relaxOrientation(controlDict.lookup("relaxOrientation")), insertSurfaceNearestPointPairs ( @@ -62,7 +63,17 @@ Foam::CV3D::controls::controls(const IOdictionary& controlDict) writeFinalTriangulation(controlDict.lookup("writeFinalTriangulation")), randomiseInitialGrid(controlDict.lookup("randomiseInitialGrid")), randomPerturbation(readScalar(controlDict.lookup("randomPerturbation"))) -{} +{ + if(globalAlignmentDirections.size() < 3) + { + FatalErrorIn("Foam::CV3D::controls::controls") + << "Too many alignment directions specified - 3 at maximum." + << exit(FatalError); + + } + + globalAlignmentDirections /= mag(globalAlignmentDirections); +} // ************************************************************************* // diff --git a/applications/utilities/mesh/generation/CV3DMesher/indexedVertex.H b/applications/utilities/mesh/generation/CV3DMesher/indexedVertex.H index 382f19d4b6..1eecba4877 100644 --- a/applications/utilities/mesh/generation/CV3DMesher/indexedVertex.H +++ b/applications/utilities/mesh/generation/CV3DMesher/indexedVertex.H @@ -35,6 +35,8 @@ Description #define indexedVertex_H #include +#include "List.H" +#include "vector.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -63,6 +65,8 @@ class indexedVertex // Lowest numbered is inside one (master). int type_; + Foam::List alignmentDirections_; + public: @@ -89,28 +93,32 @@ public: : Vb(), index_(INTERNAL_POINT), - type_(INTERNAL_POINT) + type_(INTERNAL_POINT), + alignmentDirections_(0) {} indexedVertex(const Point& p) : Vb(p), index_(INTERNAL_POINT), - type_(INTERNAL_POINT) + type_(INTERNAL_POINT), + alignmentDirections_(0) {} indexedVertex(const Point& p, Cell_handle f) : Vb(f, p), index_(INTERNAL_POINT), - type_(INTERNAL_POINT) + type_(INTERNAL_POINT), + alignmentDirections_(0) {} indexedVertex(Cell_handle f) : Vb(f), index_(INTERNAL_POINT), - type_(INTERNAL_POINT) + type_(INTERNAL_POINT), + alignmentDirections_(0) {} @@ -135,6 +143,16 @@ public: return type_; } + Foam::List& alignmentDirections() + { + return alignmentDirections_; + } + + const Foam::List& alignmentDirections() const + { + return alignmentDirections_; + } + inline bool uninitialised() const { return type_ == INTERNAL_POINT && index_ == INTERNAL_POINT;