Start of boundary handling

applications/utilities/mesh/generation/CV3DMesher/CV3D.C

@@ -162,12 +162,7 @@ void Foam::CV3D::insertGrid()
         {
             for (int k=0; k<nk; k++)
             {
-                point p1
-		(
-		    x0 + i*delta.x(),
-		    y0 + j*delta.y(),
-		    z0 + k*delta.z()
-		);
+                point p1 (x0 + i*delta.x(), y0 + j*delta.y(), z0 + k*delta.z());
 
                 point p2 = p1 + 0.5*delta;
 
@@ -220,18 +215,19 @@ void Foam::CV3D::insertSurfacePointPairs()
         insertSurfaceNearestPointPairs();
     }
 
-//     if (controls_.writeNearestTriangulation)
-//     {
-//         writeFaces("near_allFaces.obj", false);
-//         writeFaces("near_faces.obj", true);
-//         writeTriangles("near_triangles.obj", true);
-//     }
+    if (controls_.writeNearestTriangulation)
+    {   
+        // writeFaces("near_allFaces.obj", false);
+        // writeFaces("near_faces.obj", true);
+        writeTriangles("near_triangles.obj", true);
+    }
 
     if (controls_.insertSurfaceNearPointPairs)
     {
         insertSurfaceNearPointPairs();
     }
 
+    startOfBoundaryConformPointPairs_ = number_of_vertices();
 }
 
 void Foam::CV3D::boundaryConform()
@@ -240,6 +236,18 @@ void Foam::CV3D::boundaryConform()
 
 void Foam::CV3D::removeSurfacePointPairs()
 {
+    for
+    (
+        Triangulation::Finite_vertices_iterator vit = finite_vertices_begin();
+        vit != finite_vertices_end();
+        ++vit
+    )
+    {
+        if (vit->index() >= startOfSurfacePointPairs_)
+        {
+            remove(vit);
+        }
+    }
 }
 
 void Foam::CV3D::write() const

applications/utilities/mesh/generation/CV3DMesher/CV3D.H

@@ -36,7 +36,7 @@ SourceFiles
     insertFeaturePoints.C
     insertSurfaceNearestPointPairs.C
     insertSurfaceNearPointPairs.C
-
+    insertBoundaryConformPointPairs.C
 
 \*---------------------------------------------------------------------------*/
 
@@ -150,6 +150,12 @@ private:
         //- Keep track of the start of the surface point-pairs
         label startOfSurfacePointPairs_;
 
+        //- Keep track of the boundary conform point-pairs
+        //  stored after the insertion of the surface point-pairs in case
+        //  the boundary conform function is called more than once without
+        //  removing and insertin the surface point-pairs
+        label startOfBoundaryConformPointPairs_;
+
     // Private Member Functions
 
         //- Disallow default bitwise copy construct
@@ -168,6 +174,27 @@ private:
         //- Create the initial mesh from the bounding-box
         void insertBoundingBox();
 
+        //- Insert point groups at the feature points.
+        void insertFeaturePoints();
+
+        //- Insert point-pairs at the given set of points using the surface
+        //  normals corresponding to the given set of surface triangles
+        //  and write the inserted point locations to the given file.
+        void insertPointPairs
+        (
+            const DynamicList<point>& nearSurfacePoints,
+            const DynamicList<point>& surfacePoints,
+            const DynamicList<label>& surfaceTris,
+            const fileName fName
+        );
+
+        //- Check to see if dual cell specified by given vertex iterator
+        //  intersects the boundary and hence reqires a point-pair.
+        bool dualCellSurfaceIntersection
+        (
+            const Triangulation::Finite_vertices_iterator& vit
+        ) const;
+
         //- Insert point-pairs at the nearest points on the surface to the 
         //  control vertex of dual-cells which intersect the boundary in order
         //  to provide a boundary-layer mesh.
@@ -179,6 +206,8 @@ private:
         //  insertSurfaceNearestPointPairs.
         void insertSurfaceNearPointPairs();
 
+        void markNearBoundaryPoints();
+
         //- Dual calculation
         void calcDualMesh
         (
@@ -240,13 +269,10 @@ public:
             //  file.  Points outside the geometry are ignored.
             void insertPoints(const fileName& pointFileName);
 
-            //- Create the initial mesh as a regular grid of points.
+            //- Create the initial mesh as a BCC lattice of points.
             //  Points outside the geometry are ignored.
             void insertGrid();
 
-            //- Insert point groups at the feature points.
-            void insertFeaturePoints();
-
             //- Insert all surface point-pairs from
             //  insertSurfaceNearestPointPairs and
             //  findIntersectionForOutsideCentroid

applications/utilities/mesh/generation/CV3DMesher/CV3DIO.C

@@ -193,23 +193,25 @@ void Foam::CV3D::writeMesh(const Time& runTime)
 
     Info << nl << "Temporary hack to produce boundary" << endl;
 
-    for
-    (
-        Triangulation::Finite_vertices_iterator vit = finite_vertices_begin();
-        vit != finite_vertices_end();
-        ++vit
-    )
-    {
-        if
-        (
-            mag(topoint(vit->point()).x()) > 1.5
-         || mag(topoint(vit->point()).y()) > 1.5
-         || mag(topoint(vit->point()).z()) > 1.5   
-        )
-        {
-	    vit->type() = Vb::FAR_POINT;
-        }
-    }
+    // for
+    // (
+    //     Triangulation::Finite_vertices_iterator vit = finite_vertices_begin();
+    //     vit != finite_vertices_end();
+    //     ++vit
+    // )
+    // {
+    //     if
+    //     (
+    //         mag(topoint(vit->point()).x()) > 1.5
+    //      || mag(topoint(vit->point()).y()) > 1.5
+    //      || mag(topoint(vit->point()).z()) > 1.5   
+    //     )
+    //     {
+    // 	    vit->type() = Vb::FAR_POINT;
+    //     }
+    // }
+
+    markNearBoundaryPoints();
 
     pointField points(0);
     faceList faces(0);

applications/utilities/mesh/generation/CV3DMesher/insertBoundaryConformPointPairs.C

@@ -0,0 +1,57 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 1991-2007 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*----------------------------------------------------------------------------*/
+
+#include "CV3D.H"
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+void Foam::CV3D::markNearBoundaryPoints()
+{
+    for
+    (
+        Triangulation::Finite_vertices_iterator vit = finite_vertices_begin();
+        vit != finite_vertices_end();
+        vit++
+    )
+    {
+        if (vit->internalPoint())
+        {
+            point vert(topoint(vit->point()));
+
+            pointIndexHit pHit = 
+                qSurf_.tree().findNearest(vert, 4*controls_.minCellSize2);
+
+            if (pHit.hit())
+            {
+                // vit->setNearBoundary();
+	        vit->type() = Vb::FAR_POINT;
+            }
+        }
+    }
+}
+
+
+// ************************************************************************* //

applications/utilities/mesh/generation/CV3DMesher/insertSurfaceNearestPointPairs.C

@@ -28,9 +28,104 @@ License
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
+bool Foam::CV3D::dualCellSurfaceIntersection
+(
+    const Triangulation::Finite_vertices_iterator& vit
+) const
+{
+  // I think that this needs to be done with facets, but think it through.
+    return false;
+}
+
+
 void Foam::CV3D::insertSurfaceNearestPointPairs()
 {
     Info<< "insertSurfaceNearestPointPairs: " << nl << endl;
+
+    label nSurfacePointsEst = number_of_vertices();
+
+    DynamicList<point> nearSurfacePoints(nSurfacePointsEst);
+    DynamicList<point> surfacePoints(nSurfacePointsEst);
+    DynamicList<label> surfaceTris(nSurfacePointsEst);
+
+    // Local references to surface mesh addressing
+    const pointField& localPoints = qSurf_.localPoints();
+    const labelListList& edgeFaces = qSurf_.edgeFaces();
+    const vectorField& faceNormals = qSurf_.faceNormals();
+    const labelListList& faceEdges = qSurf_.faceEdges();
+
+    for
+    (
+        Triangulation::Finite_vertices_iterator vit = finite_vertices_begin();
+        vit != finite_vertices_end();
+        vit++
+    )
+    {
+        if (vit->internalPoint())
+        {
+	    point vert(topoint(vit->point()));
+
+            pointIndexHit pHit = qSurf_.tree().findNearest
+            (
+                vert,
+                4*controls_.minCellSize2
+            );
+
+            if (pHit.hit())
+            {
+                vit->setNearBoundary();
+
+                // Reference to the nearest triangle
+                const labelledTri& f = qSurf_[pHit.index()];
+
+                // Find where point is on triangle.
+                // Note tolerance needed is relative one
+                // (used in comparing normalized [0..1] triangle coordinates).
+                label nearType, nearLabel;
+                triPointRef
+                (
+                    localPoints[f[0]],
+                    localPoints[f[1]],
+                    localPoints[f[2]]
+                ).classify(pHit.hitPoint(), 1e-6, nearType, nearLabel);
+
+                // If point is on a edge check if it is an internal feature
+
+                bool internalFeatureEdge = false;
+
+                if (nearType == triPointRef::EDGE)
+                {
+                    label edgeI = faceEdges[pHit.index()][nearLabel];
+                    const labelList& eFaces = edgeFaces[edgeI];
+
+                    if
+                    (
+                        eFaces.size() == 2
+                     && (faceNormals[eFaces[0]] & faceNormals[eFaces[1]])
+                       < -0.2
+                    )
+                    {
+                        internalFeatureEdge = true;
+                    }
+                }
+
+               if (!internalFeatureEdge && dualCellSurfaceIntersection(vit))
+                {
+                    nearSurfacePoints.append(vert);
+                    surfacePoints.append(pHit.hitPoint());
+                    surfaceTris.append(pHit.index());
+                }
+	    }
+	}
+    }
+
+    // insertPointPairs
+    // (
+    //     nearSurfacePoints,
+    //     surfacePoints,
+    //     surfaceTris,
+    //     "surfaceNearestIntersections.obj"
+    // );
 }