Relaxation method experimenting.

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

@@ -202,28 +202,27 @@ void Foam::CV3D::setVertexAlignmentDirections()
 
                     alignmentDirections[2] = nt;
 
-                    Info<< "internal " << vit->internalPoint()
+                    // Info<< "internal " << vit->internalPoint()
-                        << nl << alignmentDirections
+                    //     << nl << alignmentDirections
-                        << nl << "v " << vert + alignmentDirections[0]
+                    //     << nl << "v " << vert + alignmentDirections[0]
-                        << nl << "v " << vert + alignmentDirections[1]
+                    //     << nl << "v " << vert + alignmentDirections[1]
-                        << nl << "v " << vert + alignmentDirections[2]
+                    //     << nl << "v " << vert + alignmentDirections[2]
-                        << nl << "v " << vert
+                    //     << nl << "v " << vert
-                        << nl << "v " << pHit.hitPoint()
+                    //     << nl << "v " << pHit.hitPoint()
-                        << nl << "v " << closestSpokeHitPoint
+                    //     << nl << "v " << closestSpokeHitPoint
-                        << nl << "f 4 1"
+                    //     << nl << "f 4 1"
-                        << nl << "f 4 2"
+                    //     << nl << "f 4 2"
-                        << nl << "f 4 3"
+                    //     << nl << "f 4 3"
-                        << nl << endl;
+                    //     << nl << endl;
                 }
                 else
                 {
-                    // Using on the primary alignment and the closest of the
+                    // Using only the primary alignment
-                    // globalAlignmentDirections
 
-                    alignmentDirections.setSize(0);
+                    alignmentDirections.setSize(1);
+
+                    alignmentDirections[0] = np;
                 }
-
-
             }
             else
             {
@@ -257,6 +256,30 @@ Foam::scalar Foam::CV3D::alignmentDistanceWeight(scalar dist) const
     return w;
 }
 
+
+Foam::scalar Foam::CV3D::faceAreaWeight(scalar faceArea) const
+{
+    scalar fl2 = 0.2;
+
+    scalar fu2 = 1.0;
+
+    scalar m2 = controls_.minCellSize2;
+
+    if (faceArea < fl2*m2)
+    {
+        return 0;
+    }
+    else if (faceArea < fu2*m2)
+    {
+        return faceArea/((fu2 - fl2)*m2) - 1/((fu2/fl2) - 1);
+    }
+    else
+    {
+        return 1;
+    }
+}
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::CV3D::CV3D
@@ -563,19 +586,31 @@ void Foam::CV3D::relaxPoints(const scalar relaxation)
 
             point dualFaceCentre(dualFace.centre(dualVertices));
 
-            scalar weight = 1.0;
+            vector rAB = dVA - dVB;
+
+            scalar rABMag = mag(rAB);
+
+            scalar faceArea = dualFace.mag(dualVertices);
+
+            scalar directStiffness = 2.0*relaxation;
+
+            scalar transverseStiffness = 0.0001*relaxation;
+
+            scalar r0 = 0.9*controls_.minCellSize;
+
+            vector dA = -directStiffness*(1 - r0/rABMag)
+            *faceAreaWeight(faceArea)*rAB;
+
+            vector dT = transverseStiffness*faceAreaWeight(faceArea)
+            *(dualFaceCentre - 0.5*(dVA - dVB));
 
             if (vA->internalPoint())
             {
-               //displacementAccumulator[vA->index()] = vA->index()*vector::one;
+                displacementAccumulator[vA->index()] += dA + dT;
-                displacementAccumulator[vA->index()] +=
-                relaxation*weight*(dualFaceCentre - dVA);
             }
             if (vB->internalPoint())
             {
-               //displacementAccumulator[vB->index()] = vB->index()*vector::one;
+                displacementAccumulator[vB->index()] += -dA + dT;
-                displacementAccumulator[vB->index()] +=
-                relaxation*weight*(dualFaceCentre - dVB);
             }
         }
     }

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

@@ -113,9 +113,6 @@ 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;
@@ -253,6 +250,8 @@ private:
             const vector& n
         );
 
+        inline void insertVb(const Vb& v);
+
         inline void movePoint(const Vertex_handle& vh, const point& p);
 
         //- Create the initial mesh from the bounding-box
@@ -338,6 +337,8 @@ private:
 
         scalar alignmentDistanceWeight(scalar dist) const;
 
+        scalar faceAreaWeight(scalar faceArea) const;
+
 
 public:
 

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

@@ -92,6 +92,26 @@ inline void Foam::CV3D::insertPointPair
 }
 
 
+inline void Foam::CV3D::insertVb(const Vb& v)
+{
+    const Point& Pt(v.point());
+
+    uint nVert = number_of_vertices();
+
+    Vertex_handle vh = insert(Pt);
+
+    if (nVert == number_of_vertices())
+    {
+        FatalErrorIn("Foam::CV3D::insert(const Vb& v")
+        << "Failed to reinsert Vb at " << topoint(Pt)
+            << endl;
+    }
+
+    vh->index() = v.index();
+    vh->type() = v.type();
+}
+
+
 void Foam::CV3D::movePoint(const Vertex_handle& vh, const point& p)
 {
     label vi = vh->index();

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

@@ -45,7 +45,6 @@ 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
     (
@@ -63,17 +62,7 @@ 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);
-}
 
 
 // ************************************************************************* //

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

@@ -578,8 +578,6 @@ void Foam::CV3D::insertFeaturePoints()
         vit++
     )
     {
-        // featureConstrainingVertices_[featPtI] = vit;
-
         featureConstrainingVertices_[featPtI] = Vb(vit->point());
 
         featureConstrainingVertices_[featPtI].index() = vit->index();
@@ -602,24 +600,9 @@ void Foam::CV3D::reinsertFeaturePoints()
 {
     if (featureConstrainingVertices_.size())
     {
-
         forAll(featureConstrainingVertices_, f)
         {
-            const Point& fPt(featureConstrainingVertices_[f].point());
+            insertVb(featureConstrainingVertices_[f]);
-
-            uint nVert = number_of_vertices();
-
-            Vertex_handle vh = insert(fPt);
-
-            if (nVert == number_of_vertices())
-            {
-                FatalErrorIn("Foam::CV3D::reinsertFeaturePoints")
-                    << "Failed to reinsert feature point " << topoint(fPt)
-                    << endl;
-            }
-
-            vh->index() = featureConstrainingVertices_[f].index();
-            vh->type() = featureConstrainingVertices_[f].type();
         }
     }
     else