Merge branch 'cvm' of /home/noisy3/OpenFOAM/OpenFOAM-dev into cvm

src/lagrangian/intermediate/submodels/SubModelBase.C

@@ -90,6 +90,21 @@ const Foam::dictionary& Foam::SubModelBase<CloudType>::coeffDict() const
 }
 
 
+template<class CloudType>
+bool Foam::SubModelBase<CloudType>::defaultCoeffs(const bool printMsg) const
+{
+    bool def = coeffDict_.lookupOrDefault<bool>("defaultCoeffs", false);
+    if (printMsg && def)
+    {
+        Info<< incrIndent;
+        Info<< indent << "Employing default coefficients" << endl;
+        Info<< decrIndent;
+    }
+
+    return def;
+}
+
+
 template<class CloudType>
 CloudType& Foam::SubModelBase<CloudType>::owner()
 {
@@ -105,7 +120,7 @@ bool Foam::SubModelBase<CloudType>::active() const
 
 
 template<class CloudType>
-void  Foam::SubModelBase<CloudType>::cacheFields(const bool)
+void Foam::SubModelBase<CloudType>::cacheFields(const bool)
 {
     // do nothing
 }

src/lagrangian/intermediate/submodels/SubModelBase.H

@@ -63,11 +63,10 @@ protected:
         //- Reference to the cloud dictionary
         const dictionary& dict_;
 
-        //- Reference to the coefficients dictionary
+        //- Coefficients dictionary
         const dictionary& coeffDict_;
 
 
-
 public:
 
     // Constructors
@@ -108,6 +107,9 @@ public:
             //- Return const access to the coefficients dictionary
             const dictionary& coeffDict() const;
 
+            //- Returns true if defaultCoeffs is true and outputs on printMsg
+            bool defaultCoeffs(const bool printMsg) const;
+
             //- Return the model 'active' status - default active = true
             virtual bool active() const;
 

src/lagrangian/spray/submodels/BreakupModel/ETAB/ETAB.C

@@ -35,24 +35,24 @@ Foam::ETAB<CloudType>::ETAB
 )
 :
     BreakupModel<CloudType>(dict, owner, typeName),
-    Cmu_(this->coeffDict().template lookupOrDefault<scalar>("Cmu", 10.0)),
+    Cmu_(10.0),
-    Comega_(this->coeffDict().template lookupOrDefault<scalar>("Comega", 8.0)),
+    Comega_(8.0),
-    k1_(this->coeffDict().template lookupOrDefault<scalar>("k1", 0.2)),
+    k1_(0.2),
-    k2_(this->coeffDict().template lookupOrDefault<scalar>("k2", 0.2)),
+    k2_(0.2),
-    WeCrit_
+    WeCrit_(12.0),
-(
+    WeTransition_(100.0),
-    this->coeffDict().template lookupOrDefault<scalar>("WeCrit", 12.0)
-),
-    WeTransition_
-    (
-        this->coeffDict().template lookupOrDefault<scalar>
-        (
-            "WeTransition",
-            100.0
-        )
-    ),
     AWe_(0.0)
 {
+    if (!this->defaultCoeffs(true))
+    {
+        this->coeffDict().lookup("Cmu") >> Cmu_;
+        this->coeffDict().lookup("Comega") >> Comega_;
+        this->coeffDict().lookup("k1") >> k1_;
+        this->coeffDict().lookup("k2") >> k2_;
+        this->coeffDict().lookup("WeCrit") >> WeCrit_;
+        this->coeffDict().lookup("WeTransition") >> WeTransition_;
+    }
+
     scalar k21 = k2_/k1_;
     AWe_ = (k21*sqrt(WeTransition_) - 1.0)/pow4(WeTransition_);
 

src/lagrangian/spray/submodels/BreakupModel/PilchErdman/PilchErdman.C

@@ -35,9 +35,15 @@ Foam::PilchErdman<CloudType>::PilchErdman
 )
 :
     BreakupModel<CloudType>(dict, owner, typeName),
-    B1_(this->coeffDict().template lookupOrDefault<scalar>("B1", 0.375)),
+    B1_(0.375),
-    B2_(this->coeffDict().template lookupOrDefault<scalar>("B2", 0.236))
+    B2_(0.236)
-{}
+{
+    if (!this->defaultCoeffs(true))
+    {
+        this->coeffDict().lookup("B1") >> B1_;
+        this->coeffDict().lookup("B2") >> B2_;
+    }
+}
 
 
 template <class CloudType>

src/lagrangian/spray/submodels/BreakupModel/ReitzDiwakar/ReitzDiwakar.C

@@ -35,11 +35,19 @@ Foam::ReitzDiwakar<CloudType>::ReitzDiwakar
 )
 :
     BreakupModel<CloudType>(dict, owner, typeName),
-    Cbag_(this->coeffDict().template lookupOrDefault<scalar>("Cbag", 6.0)),
+    Cbag_(6.0),
-    Cb_(this->coeffDict().template lookupOrDefault<scalar>("Cb", 0.785)),
+    Cb_(0.785),
-    Cstrip_(this->coeffDict().template lookupOrDefault<scalar>("Cstrip", 0.5)),
+    Cstrip_(0.5),
-    Cs_(this->coeffDict().template lookupOrDefault<scalar>("Cs", 10.0))
+    Cs_(10.0)
-{}
+{
+    if (!this->defaultCoeffs(true))
+    {
+        this->coeffDict().lookup("Cbag") >> Cbag_;
+        this->coeffDict().lookup("Cb") >> Cb_;
+        this->coeffDict().lookup("Cstrip") >> Cstrip_;
+        this->coeffDict().lookup("Cs") >> Cs_;
+    }
+}
 
 
 template <class CloudType>

src/lagrangian/spray/submodels/BreakupModel/ReitzKHRT/ReitzKHRT.C

@@ -35,20 +35,23 @@ Foam::ReitzKHRT<CloudType>::ReitzKHRT
 )
 :
     BreakupModel<CloudType>(dict, owner, typeName),
-    b0_(this->coeffDict().template lookupOrDefault<scalar>("B0", 0.61)),
+    b0_(0.61),
-    b1_(this->coeffDict().template lookupOrDefault<scalar>("B1", 40.0)),
+    b1_(40.0),
-    cTau_(this->coeffDict().template lookupOrDefault<scalar>("Ctau", 1.0)),
+    cTau_(1.0),
-    cRT_(this->coeffDict().template lookupOrDefault<scalar>("CRT", 0.1)),
+    cRT_(0.1),
-    msLimit_
+    msLimit_(0.03),
-    (
+    weberLimit_(6.0)
-        this->coeffDict().template lookupOrDefault<scalar>
+{
-        (
+    if (!this->defaultCoeffs(true))
-            "msLimit",
+    {
-            0.03
+        this->coeffDict().lookup("B0") >> b0_;
-        )
+        this->coeffDict().lookup("B1") >> b1_;
-    ),
+        this->coeffDict().lookup("Ctau") >> cTau_;
-    weberLimit_(readScalar(this->coeffDict().lookup("WeberLimit")))
+        this->coeffDict().lookup("CRT") >> cRT_;
-{}
+        this->coeffDict().lookup("msLimit") >> msLimit_;
+        this->coeffDict().lookup("WeberLimit") >> weberLimit_;
+    }
+}
 
 
 template <class CloudType>

src/lagrangian/spray/submodels/BreakupModel/SHF/SHF.C

@@ -142,7 +142,7 @@ bool Foam::SHF<CloudType>::update
 {
     bool addChild = false;
 
-    scalar d03 = pow(d, 3);
+    scalar d03 = pow3(d);
     scalar rhopi6 = rho*constant::mathematical::pi/6.0;
     scalar mass0 = nParticle*rhopi6*d03;
     scalar mass = mass0;

src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C

@@ -256,9 +256,21 @@ void Foam::conformalVoronoiMesh::insertPoints
 
         label sizeChange = preDistributionSize - label(points.size());
 
-        if (mag(sizeChange) > 0)
+        // if (mag(sizeChange) > 0)
+        // {
+        //     Pout<< "    distribution points size change " << sizeChange
+        //         << endl;
+        // }
+
+        label totalMagSizeChange = returnReduce
+        (
+            mag(sizeChange), sumOp<label>()
+        );
+
+        if (totalMagSizeChange > 0)
         {
-            Pout<< "    distribution points size change " << sizeChange
+            Info<< "    distribution points size change total "
+                << totalMagSizeChange/2
                 << endl;
         }
 
@@ -272,18 +284,19 @@ void Foam::conformalVoronoiMesh::insertPoints
 
     label nVert = number_of_vertices();
 
-    Info<< "TEMPORARILY USING INDIVIDUAL INSERTION TO DETECT FAILURE" << endl;
     // using the range insert (faster than inserting points one by one)
-    // insert(points.begin(), points.end());
+    insert(points.begin(), points.end());
-    for
+
-    (
+    // Info<< "USING INDIVIDUAL INSERTION TO DETECT FAILURE" << endl;
-        std::list<Point>::iterator pit=points.begin();
+    // for
-        pit != points.end();
+    // (
-        ++pit
+    //     std::list<Point>::iterator pit=points.begin();
-    )
+    //     pit != points.end();
-    {
+    //     ++pit
-        insertVb(*pit);
+    // )
-    }
+    // {
+    //     insertPoint(topoint(*pit), Vb::ptInternalPoint);
+    // }
 
     label nInserted(number_of_vertices() - preInsertionSize);
 
@@ -317,6 +330,10 @@ void Foam::conformalVoronoiMesh::insertPoints
     bool distribute
 )
 {
+    // The pts, indices and types lists must be intact and up-to-date at the
+    // end of this function as they may also be used by other functions
+    // subsequently.
+
     if (Pstream::parRun() && distribute)
     {
         // The link between vertices that form the boundary via pairs cannot be
@@ -362,7 +379,6 @@ void Foam::conformalVoronoiMesh::insertPoints
         if (type > Vb::ptFarPoint)
         {
             // This is a member of a point pair, don't use the type directly
-
             type += number_of_vertices();
         }
 
@@ -749,6 +765,7 @@ void Foam::conformalVoronoiMesh::insertFeaturePoints()
 
     createMixedFeaturePoints(pts, indices, types);
 
+    // Insert the created points, distributing to the appropriate processor
     insertPoints(pts, indices, types, true);
 
     if(cvMeshControls().objOutput())
@@ -769,27 +786,15 @@ void Foam::conformalVoronoiMesh::insertFeaturePoints()
             << " feature vertices" << endl;
     }
 
-    Info<< "SORT OUT FEATURE POINT DISTRIBUTION AND STORAGE" << endl;
+    featureVertices_.clear();
 
-    // featureVertices_.setSize(number_of_vertices());
+    forAll(pts, pI)
-
+    {
-    // label featPtI = 0;
+        featureVertices_.push_back
-
+        (
-    // for
+            Vb(toPoint(pts[pI]), indices[pI], types[pI])
-    // (
+        );
-    //     Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
+    }
-    //     vit != finite_vertices_end();
-    //     vit++
-    // )
-    // {
-    //     featureVertices_[featPtI] = Vb(vit->point());
-
-    //     featureVertices_[featPtI].index() = vit->index();
-
-    //     featureVertices_[featPtI].type() = vit->type();
-
-    //     featPtI++;
-    // }
 
     constructFeaturePointLocations();
 }
@@ -1065,14 +1070,72 @@ void Foam::conformalVoronoiMesh::constructFeaturePointLocations()
 }
 
 
-void Foam::conformalVoronoiMesh::reinsertFeaturePoints()
+void Foam::conformalVoronoiMesh::reinsertFeaturePoints(bool distribute)
 {
     Info<< nl << "Reinserting stored feature points" << endl;
 
-    forAll(featureVertices_, f)
+    label preReinsertionSize(number_of_vertices());
+
+    if (distribute)
     {
-        insertVb(featureVertices_[f]);
+        DynamicList<Foam::point> pointsToInsert;
+        DynamicList<label> indices;
+        DynamicList<label> types;
+
+        for
+        (
+            std::list<Vb>::iterator vit=featureVertices_.begin();
+            vit != featureVertices_.end();
+            ++vit
+        )
+        {
+            pointsToInsert.append(topoint(vit->point()));
+            indices.append(vit->index());
+            types.append(vit->type());
+        }
+
+        // Insert distributed points
+        insertPoints(pointsToInsert, indices, types, true);
+
+        // Save points in new distribution
+        featureVertices_.clear();
+
+        forAll(pointsToInsert, pI)
+        {
+            featureVertices_.push_back
+            (
+                Vb(toPoint(pointsToInsert[pI]), indices[pI], types[pI])
+            );
+        }
     }
+    else
+    {
+        for
+        (
+            std::list<Vb>::iterator vit=featureVertices_.begin();
+            vit != featureVertices_.end();
+            ++vit
+        )
+        {
+            // Assuming that all of the reinsertions are pair points, and that
+            // the index and type are relative, i.e. index 0 and type relative
+            // to it.
+            insertPoint
+            (
+                vit->point(),
+                vit->index() + number_of_vertices(),
+                vit->type() + number_of_vertices()
+            );
+        }
+    }
+
+    Info<< "    Reinserted "
+        << returnReduce
+        (
+            label(number_of_vertices()) - preReinsertionSize,
+            sumOp<label>()
+        )
+        << " vertices" << endl;
 }
 
 
@@ -1167,12 +1230,6 @@ void Foam::conformalVoronoiMesh::insertInitialPoints()
         writePoints("initialPoints.obj", true);
     }
 
-    Info<< "NEED TO CHANGE storeSizesAndAlignments AFTER DISTRIBUTE" << endl;
-    Info<< "NEED TO MAP FEATURE POINTS AFTER DISTRIBUTE" << endl;
-    Info<< "NEED TO ENSURE THAT FEATURE POINTS ARE INSERTED ON THE "
-        << "CORRECT PRCOESSOR" << endl;
-
-    storeSizesAndAlignments(initPts);
 }
 
 
@@ -1200,7 +1257,8 @@ bool Foam::conformalVoronoiMesh::distributeBackground()
             vit++
         )
         {
-            if (vit->real())
+            // Only store real vertices that are not feature vertices
+            if (vit->real() && vit->index() >= startOfInternalPoints_)
             {
                 nRealVertices++;
             }
@@ -1264,7 +1322,8 @@ bool Foam::conformalVoronoiMesh::distributeBackground()
             vit++
         )
         {
-            if (vit->real())
+            // Only store real vertices that are not feature vertices
+            if (vit->real() && vit->index() >= startOfInternalPoints_)
             {
                 Foam::point v = topoint(vit->point());
 
@@ -1313,8 +1372,8 @@ bool Foam::conformalVoronoiMesh::distributeBackground()
         // Remove the entire tessellation
         reset();
 
-        Info<< "NEED TO MAP FEATURE POINTS" << endl;
+        // Reinsert feature points, distributing them as necessary.
-        // reinsertFeaturePoints();
+        reinsertFeaturePoints(true);
 
         startOfInternalPoints_ = number_of_vertices();
 
@@ -1331,6 +1390,8 @@ bool Foam::conformalVoronoiMesh::distributeBackground()
             forAll(cellVertices[cI], cVPI)
             {
                 pointsToInsert.append(cellVertices[cI][cVPI]);
+
+                // All insertions relative to index of zero
                 indices.append(0);
 
                 label type = cellVertexTypes[cI][cVPI];
@@ -1374,6 +1435,27 @@ bool Foam::conformalVoronoiMesh::distributeBackground()
 }
 
 
+void Foam::conformalVoronoiMesh::storeSizesAndAlignments()
+{
+    std::list<Point> storePts;
+
+    for
+    (
+        Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
+        vit != finite_vertices_end();
+        vit++
+    )
+    {
+        if (vit->internalPoint())
+        {
+            storePts.push_back(vit->point());
+        }
+    }
+
+    storeSizesAndAlignments(storePts);
+}
+
+
 void Foam::conformalVoronoiMesh::storeSizesAndAlignments
 (
     const std::list<Point>& storePts
@@ -1432,6 +1514,8 @@ void Foam::conformalVoronoiMesh::updateSizesAndAlignments
     )
     {
         storeSizesAndAlignments(storePts);
+
+        timeCheck("Updated sizes and alignments");
     }
 }
 
@@ -1799,6 +1883,8 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh
 
     insertBoundingPoints();
 
+    insertFeaturePoints();
+
     insertInitialPoints();
 
     // Improve the guess that the backgroundMeshDecomposition makes with the
@@ -1806,8 +1892,6 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh
     // better balance the surface conformation load.
     distributeBackground();
 
-    insertFeaturePoints();
-
     buildSurfaceConformation(rmCoarse);
 
     // The introduction of the surface conformation may have distorted the
@@ -1820,6 +1904,21 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh
         buildParallelInterface("rebuild");
     }
 
+    // Do not store the surface conformation until after it has been
+    // (potentially) redistributed.
+    storeSurfaceConformation();
+
+    // Use storeSizesAndAlignments with no feed points because all background
+    // points may have been distributed.  It is a requirement that none of the
+    // preceding functions requires look up of sizes or alignments from the
+    // Delaunay vertices, i.e. setVertexSizeAndAlignment cannot be called
+    // before this point.
+    storeSizesAndAlignments();
+
+    // Report any Delaunay vertices that do not think that they are in the
+    // domain the processor they are on.
+    reportProcessorOccupancy();
+
     if(cvMeshControls().objOutput())
     {
         writePoints("allInitialPoints.obj", false);
@@ -1897,8 +1996,6 @@ void Foam::conformalVoronoiMesh::move()
 
     std::list<Point> pointsToInsert;
 
-    label pointsAdded = 0;
-
     for
     (
         Delaunay::Finite_edges_iterator eit = finite_edges_begin();
@@ -1906,18 +2003,18 @@ void Foam::conformalVoronoiMesh::move()
         ++eit
     )
     {
+        Cell_handle c = eit->first;
+        Vertex_handle vA = c->vertex(eit->second);
+        Vertex_handle vB = c->vertex(eit->third);
+
         if
         (
-            eit->first->vertex(eit->second)->internalOrBoundaryPoint()
+            vA->internalOrBoundaryPoint()
-         && eit->first->vertex(eit->third)->internalOrBoundaryPoint()
+         && vB->internalOrBoundaryPoint()
         )
         {
             face dualFace = buildDualFace(eit);
 
-            Cell_handle c = eit->first;
-            Vertex_handle vA = c->vertex(eit->second);
-            Vertex_handle vB = c->vertex(eit->third);
-
             Foam::point dVA = topoint(vA->point());
             Foam::point dVB = topoint(vB->point());
 
@@ -2056,7 +2153,6 @@ void Foam::conformalVoronoiMesh::move()
                                 toPoint(0.5*(dVA + dVB))
                             );
 
-                            pointsAdded++;
                         }
                     }
                     else if
@@ -2133,8 +2229,8 @@ void Foam::conformalVoronoiMesh::move()
         }
     }
 
-    vector totalDisp = sum(displacementAccumulator);
+    vector totalDisp = gSum(displacementAccumulator);
-    scalar totalDist = sum(mag(displacementAccumulator));
+    scalar totalDist = gSum(mag(displacementAccumulator));
 
     // Relax the calculated displacement
     displacementAccumulator *= relaxation;
@@ -2172,13 +2268,6 @@ void Foam::conformalVoronoiMesh::move()
         }
     }
 
-    // Write the mesh before clearing it.  Beware that writeMesh destroys the
-    // indexing of the tessellation.  Do not filter the dual faces.
-    if (runTime_.outputTime())
-    {
-        writeMesh(runTime_.timeName(), false);
-    }
-
     // Remove the entire tessellation
     reset();
 
@@ -2192,11 +2281,6 @@ void Foam::conformalVoronoiMesh::move()
 
     insertPoints(pointsToInsert);
 
-    label pointsRemoved =
-        displacementAccumulator.size()
-      - number_of_vertices()
-      + pointsAdded;
-
     timeCheck("Internal points inserted");
 
     conformToSurface();
@@ -2205,14 +2289,16 @@ void Foam::conformalVoronoiMesh::move()
 
     updateSizesAndAlignments(pointsToInsert);
 
+    // Write the intermediate mesh, do not filter the dual faces.
+    if (runTime_.outputTime())
+    {
+        writeMesh(runTime_.timeName(), false);
+    }
+
     Info<< nl
         << "Total displacement = " << totalDisp << nl
         << "Total distance = " << totalDist << nl
-        << "Points added = " << pointsAdded << nl
-        << "Points removed = " << pointsRemoved
         << endl;
-
-    timeCheck("Updated sizes and alignments (if required), end of move");
 }
 
 

src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H

@@ -138,7 +138,7 @@ private:
 
         //- Store the feature constraining points to be reinserted after a
         //  triangulation clear
-        List<Vb> featureVertices_;
+        std::list<Vb> featureVertices_;
 
         //- Storing the locations of all of the features to be conformed to.
         //  Single pointField required by the featurePointTree.
@@ -240,17 +240,32 @@ private:
         //- Return the required alignment directions at the given location
         tensor requiredAlignment(const Foam::point& pt) const;
 
-        //- Insert point and return its index
+        //- Insert Foam::point and return its auto-generated index
         inline label insertPoint
         (
-            const Foam::point& pt,
+            const Foam::point& p,
             const label type
         );
 
-        //- Insert point and return its index
+        //- Insert Point and return its auto-generated index
+        inline label insertPoint
+        (
+            const Point& P,
+            const label type
+        );
+
+        //- Insert Foam::point with specified index and type
         inline void insertPoint
         (
-            const Foam::point& pt,
+            const Foam::point& p,
+            const label index,
+            const label type
+        );
+
+        //- Insert Point with specified index and type
+        inline void insertPoint
+        (
+            const Point& P,
             const label index,
             const label type
         );
@@ -295,10 +310,6 @@ private:
             DynamicList<label>& types
         );
 
-        //- Insert a Vb (a typedef of CGAL::indexedVertex<K>) with the
-        //  possibility of an offset for the index and the type
-        inline void insertVb(const Vb& v, label offset = 0);
-
         //- Insert pairs of points on the surface with the given normals, at the
         //  specified spacing
         void insertSurfacePointPairs
@@ -420,7 +431,7 @@ private:
         void constructFeaturePointLocations();
 
         //- Reinsert stored feature point defining points
-        void reinsertFeaturePoints();
+        void reinsertFeaturePoints(bool distribute = false);
 
         //- Demand driven construction of octree for feature points
         const indexedOctree<treeDataPoint>& featurePointTree() const;
@@ -445,6 +456,12 @@ private:
         //  referred vertices, so the parallel interface may need rebuilt.
         bool distributeBackground();
 
+        //- Store data for sizeAndAlignmentLocations_, storedSizes_ and
+        //  storedAlignments_ and initialise the sizeAndAlignmentTreePtr_,
+        //  determining the appropriate sizeAndAlignmentLocations_
+        //  automatically
+        void storeSizesAndAlignments();
+
         //- Store data for sizeAndAlignmentLocations_, storedSizes_ and
         //  storedAlignments_ and initialise the sizeAndAlignmentTreePtr_
         void storeSizesAndAlignments(const std::list<Point>& storePts);
@@ -587,6 +604,9 @@ private:
             label& hitSurfaceLargest
         ) const;
 
+        //- Write out vertex-processor occupancy information for debugging
+        void reportProcessorOccupancy();
+
         //- Write out debugging information about the surface conformation
         //  quality
         void reportSurfaceConformationQuality();

src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C

@@ -1457,65 +1457,65 @@ Foam::labelHashSet Foam::conformalVoronoiMesh::checkPolyMeshQuality
 
     List<polyPatch*> patches(patchStarts.size());
 
-    // label nValidPatches = 0;
+    label nValidPatches = 0;
-
-    // forAll(patches, p)
-    // {
-    //     if (patchTypes[p] == processorPolyPatch::typeName)
-    //     {
-    //         // Do not create empty processor patches
-
-    //         if (patchSizes[p] > 0)
-    //         {
-    //             patches[nValidPatches] = new processorPolyPatch
-    //             (
-    //                 patchNames[p],
-    //                 patchSizes[p],
-    //                 patchStarts[p],
-    //                 nValidPatches,
-    //                 pMesh.boundaryMesh(),
-    //                 Pstream::myProcNo(),
-    //                 procNeighbours[p]
-    //             );
-
-    //             nValidPatches++;
-    //         }
-    //     }
-    //     else
-    //     {
-    //         patches[nValidPatches] = polyPatch::New
-    //         (
-    //             patchTypes[p],
-    //             patchNames[p],
-    //             patchSizes[p],
-    //             patchStarts[p],
-    //             nValidPatches,
-    //             pMesh.boundaryMesh()
-    //         ).ptr();
-
-    //         nValidPatches++;
-    //     }
-    // }
-
-    // patches.setSize(nValidPatches);
-
-    // pMesh.addPatches(patches);
-
-    Info<< "ADDPATCHES NOT IN PARALLEL" << endl;
 
     forAll(patches, p)
     {
-        patches[p] = new polyPatch
+        if (patchTypes[p] == processorPolyPatch::typeName)
-        (
+        {
-            patchNames[p],
+            // Do not create empty processor patches
-            patchSizes[p],
+
-            patchStarts[p],
+            if (patchSizes[p] > 0)
-            p,
+            {
-            pMesh.boundaryMesh()
+                patches[nValidPatches] = new processorPolyPatch
-        );
+                (
+                    patchNames[p],
+                    patchSizes[p],
+                    patchStarts[p],
+                    nValidPatches,
+                    pMesh.boundaryMesh(),
+                    Pstream::myProcNo(),
+                    procNeighbours[p]
+                );
+
+                nValidPatches++;
+            }
+        }
+        else
+        {
+            patches[nValidPatches] = polyPatch::New
+            (
+                patchTypes[p],
+                patchNames[p],
+                patchSizes[p],
+                patchStarts[p],
+                nValidPatches,
+                pMesh.boundaryMesh()
+            ).ptr();
+
+            nValidPatches++;
+        }
     }
 
-    pMesh.addPatches(patches, false);
+    patches.setSize(nValidPatches);
+
+    pMesh.addPatches(patches);
+
+    // Info<< "ADDPATCHES NOT IN PARALLEL" << endl;
+
+    // forAll(patches, p)
+    // {
+    //     patches[p] = new polyPatch
+    //     (
+    //         patchNames[p],
+    //         patchSizes[p],
+    //         patchStarts[p],
+    //         p,
+    //         pMesh.boundaryMesh()
+    //     );
+    // }
+
+    // pMesh.addPatches(patches, false);
 
     // pMesh.overrideCellCentres(cellCentres);
 

src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C

@@ -36,12 +36,18 @@ void Foam::conformalVoronoiMesh::conformToSurface()
     {
         // Reinsert stored surface conformation
         reinsertSurfaceConformation();
+
+        buildParallelInterface("move_" + runTime_.timeName());
     }
     else
     {
         // Rebuild, insert and store new surface conformation
         buildSurfaceConformation(reconfMode);
+
+        storeSurfaceConformation();
     }
+
+    // reportSurfaceConformationQuality();
 }
 
 
@@ -411,10 +417,6 @@ void Foam::conformalVoronoiMesh::buildSurfaceConformation
                 << "), stopping iterations" << endl;
         }
     }
-
-    // reportSurfaceConformationQuality();
-
-    storeSurfaceConformation();
 }
 
 
@@ -620,6 +622,8 @@ void Foam::conformalVoronoiMesh::buildParallelInterface
         receivedVertices,
         outputName
     );
+
+    timeCheck("After buildParallelInterface");
 }
 
 
@@ -1346,6 +1350,27 @@ void Foam::conformalVoronoiMesh::dualCellLargestSurfaceIncursion
 }
 
 
+void Foam::conformalVoronoiMesh::reportProcessorOccupancy()
+{
+    for
+    (
+        Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
+        vit != finite_vertices_end();
+        vit++
+    )
+    {
+        if (vit->real())
+        {
+            if (!positionOnThisProc(topoint(vit->point())))
+            {
+                Pout<< topoint(vit->point()) << " is not on this processor "
+                    << endl;
+            }
+        }
+    }
+}
+
+
 void Foam::conformalVoronoiMesh::reportSurfaceConformationQuality()
 {
     Info<< nl << "Check surface conformation quality" << endl;
@@ -1605,6 +1630,14 @@ void Foam::conformalVoronoiMesh::buildEdgeLocationTree
 
 void Foam::conformalVoronoiMesh::buildSizeAndAlignmentTree() const
 {
+    if (sizeAndAlignmentLocations_.empty())
+    {
+        FatalErrorIn("buildSizeAndAlignmentTree()")
+            << "sizeAndAlignmentLocations empty, must be populated before "
+            << "sizeAndAlignmentTree can be built."
+            << exit(FatalError);
+    }
+
     treeBoundBox overallBb
     (
         geometryToConformTo_.globalBounds().extend(rndGen_, 1e-4)
@@ -1755,7 +1788,14 @@ void Foam::conformalVoronoiMesh::storeSurfaceConformation()
         vit++
     )
     {
-        if (!vit->referred() && vit->pairPoint())
+        // Store points that are not referred, part of a pair, but not feature
+        // points
+        if
+        (
+            !vit->referred()
+         && vit->pairPoint()
+         && vit->index() >= startOfInternalPoints_
+        )
         {
             surfaceConformationVertices_.push_back
             (
@@ -1785,6 +1825,9 @@ void Foam::conformalVoronoiMesh::reinsertSurfaceConformation()
 
     label preReinsertionSize(number_of_vertices());
 
+    // It is assumed that the stored surface conformation is on the correct
+    // processor and does not need distributed
+
     for
     (
         std::list<Vb>::iterator vit=surfaceConformationVertices_.begin();
@@ -1792,7 +1835,14 @@ void Foam::conformalVoronoiMesh::reinsertSurfaceConformation()
         ++vit
     )
     {
-        insertVb(*vit, number_of_vertices());
+        // Assuming that all of the reinsertions are pair points, and that the
+        // index and type are relative, i.e. index 0 and type relative to it.
+        insertPoint
+        (
+            vit->point(),
+            vit->index() + number_of_vertices(),
+            vit->type() + number_of_vertices()
+        );
     }
 
     Info<< "    Reinserted "

src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H

@@ -215,14 +215,24 @@ inline Foam::label Foam::conformalVoronoiMesh::insertPoint
     const Foam::point& p,
     const label type
 )
+{
+    return insertPoint(toPoint(p), type);
+}
+
+
+inline Foam::label Foam::conformalVoronoiMesh::insertPoint
+(
+    const Point& P,
+    const label type
+)
 {
     uint nVert = number_of_vertices();
 
-    Vertex_handle vh = insert(toPoint(p));
+    Vertex_handle vh = insert(P);
 
     if (nVert == number_of_vertices())
     {
-        Pout << "Failed to insert point " << p << endl;
+        Pout << "Failed to insert point " << topoint(P) << endl;
     }
     else
     {
@@ -240,14 +250,25 @@ inline void Foam::conformalVoronoiMesh::insertPoint
     const label index,
     const label type
 )
+{
+    insertPoint(toPoint(p), index, type);
+}
+
+
+inline void Foam::conformalVoronoiMesh::insertPoint
+(
+    const Point& P,
+    const label index,
+    const label type
+)
 {
     uint nVert = number_of_vertices();
 
-    Vertex_handle vh = insert(toPoint(p));
+    Vertex_handle vh = insert(P);
 
     if (nVert == number_of_vertices())
     {
-        Pout << "Failed to insert point " << p << endl;
+        Pout << "Failed to insert point " << topoint(P) << endl;
     }
     else
     {
@@ -311,37 +332,6 @@ inline void Foam::conformalVoronoiMesh::createPointPair
 }
 
 
-inline void Foam::conformalVoronoiMesh::insertVb(const Vb& v, label offset)
-{
-    const Point& Pt(v.point());
-
-    uint nVert = number_of_vertices();
-
-    Vertex_handle vh = insert(Pt);
-
-    if (nVert == number_of_vertices())
-    {
-        FatalErrorIn("Foam::conformalVoronoiMesh::insertVb(const Vb& v")
-            << "Failed to reinsert Vb at " << topoint(Pt)
-            << nl
-            << exit(FatalError);
-    }
-    else
-    {
-        vh->index() = v.index() + offset;
-
-        if (v.pairPoint())
-        {
-            vh->type() = v.type() + offset;
-        }
-        else
-        {
-            vh->type() = v.type();
-        }
-    }
-}
-
-
 inline bool Foam::conformalVoronoiMesh::isBoundaryDualFace
 (
     const Delaunay::Finite_edges_iterator& eit

src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshIO.C

@@ -26,7 +26,6 @@ License
 #include "conformalVoronoiMesh.H"
 #include "IOstreams.H"
 #include "OFstream.H"
-#include "zeroGradientPointPatchField.H"
 #include "pointMesh.H"
 #include "pointFields.H"
 
@@ -312,7 +311,7 @@ void Foam::conformalVoronoiMesh::writeMesh
         writeObjMesh(points, faces, word(meshName + ".obj"));
     }
 
-    polyMesh mesh
+    fvMesh mesh
     (
         IOobject
         (
@@ -372,9 +371,10 @@ void Foam::conformalVoronoiMesh::writeMesh
 
     patches.setSize(nValidPatches);
 
-    Info<< "ADDPATCHES NOT IN PARALLEL" << endl;
+    mesh.addFvPatches(patches);
-    // mesh.addPatches(patches);
+
-    mesh.addPatches(patches, false);
+    // Info<< "ADDPATCHES NOT IN PARALLEL" << endl;
+    // mesh.addPatches(patches, false);
 
     if (!mesh.write())
     {
@@ -404,11 +404,9 @@ void Foam::conformalVoronoiMesh::writeMesh
 
     // cellCs.write();
 
-    Info<< "DISABLED WRITING OF CELL SIZE AND PROTRUSION SET" << endl;
+    writeCellSizes(mesh);
 
-    // writeCellSizes(mesh);
+    findRemainingProtrusionSet(mesh);
-
-    // findRemainingProtrusionSet(mesh);
 }
 
 
@@ -466,7 +464,7 @@ void Foam::conformalVoronoiMesh::writeCellSizes
             ),
             mesh,
             dimensionedScalar("cellSize", dimLength, 0),
-            zeroGradientPointPatchField<scalar>::typeName
+            zeroGradientFvPatchScalarField::typeName
         );
 
         scalarField& cellSize = targetCellSize.internalField();
@@ -492,7 +490,7 @@ void Foam::conformalVoronoiMesh::writeCellSizes
         //     ),
         //     mesh,
         //     dimensionedScalar("cellVolume", dimLength, 0),
-        //     zeroGradientPointPatchField<scalar>::typeName
+        //     zeroGradientFvPatchScalarField::typeName
         // );
 
         // targetCellVolume.internalField() = pow3(cellSize);
@@ -511,7 +509,7 @@ void Foam::conformalVoronoiMesh::writeCellSizes
         //     ),
         //     mesh,
         //     dimensionedScalar("cellVolume", dimVolume, 0),
-        //     zeroGradientPointPatchField<scalar>::typeName
+        //     zeroGradientFvPatchScalarField::typeName
         // );
 
         // actualCellVolume.internalField() = mesh.cellVolumes();
@@ -530,7 +528,7 @@ void Foam::conformalVoronoiMesh::writeCellSizes
         //     ),
         //     mesh,
         //     dimensionedScalar("cellSize", dimLength, 0),
-        //     zeroGradientPointPatchField<scalar>::typeName
+        //     zeroGradientFvPatchScalarField::typeName
         // );
 
         // equivalentCellSize.internalField() = pow

src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex.H

@@ -342,6 +342,7 @@ public:
         return internalPoint(type) || ppMaster(index, type);
     }
 
+
     //- Is point near the boundary or part of the boundary definition
     inline bool nearOrOnBoundary() const
     {
@@ -349,43 +350,43 @@ public:
     }
 
 
-    //- Do the two given vertices consitute a boundary point-pair
+    // //- Do the two given vertices consitute a boundary point-pair
-    inline friend bool pointPair
+    // inline friend bool pointPair
-    (
+    // (
-        const indexedVertex& v0,
+    //     const indexedVertex& v0,
-        const indexedVertex& v1
+    //     const indexedVertex& v1
-    )
+    // )
-    {
+    // {
-        return v0.index_ == v1.type_ || v1.index_ == v0.type_;
+    //     return v0.index_ == v1.type_ || v1.index_ == v0.type_;
-    }
+    // }
 
 
-    //- Do the three given vertices consitute a boundary triangle
+    // //- Do the three given vertices consitute a boundary triangle
-    inline friend bool boundaryTriangle
+    // inline friend bool boundaryTriangle
-    (
+    // (
-        const indexedVertex& v0,
+    //     const indexedVertex& v0,
-        const indexedVertex& v1,
+    //     const indexedVertex& v1,
-        const indexedVertex& v2
+    //     const indexedVertex& v2
-    )
+    // )
-    {
+    // {
-        return (v0.pairPoint() && pointPair(v1, v2))
+    //     return (v0.pairPoint() && pointPair(v1, v2))
-            || (v1.pairPoint() && pointPair(v2, v0))
+    //         || (v1.pairPoint() && pointPair(v2, v0))
-            || (v2.pairPoint() && pointPair(v0, v1));
+    //         || (v2.pairPoint() && pointPair(v0, v1));
-    }
+    // }
 
 
-    //- Do the three given vertices consitute an outside triangle
+    // //- Do the three given vertices consitute an outside triangle
-    inline friend bool outsideTriangle
+    // inline friend bool outsideTriangle
-    (
+    // (
-        const indexedVertex& v0,
+    //     const indexedVertex& v0,
-        const indexedVertex& v1,
+    //     const indexedVertex& v1,
-        const indexedVertex& v2
+    //     const indexedVertex& v2
-    )
+    // )
-    {
+    // {
-        return (v0.farPoint() || v0.ppSlave())
+    //     return (v0.farPoint() || v0.ppSlave())
-            || (v1.farPoint() || v1.ppSlave())
+    //         || (v1.farPoint() || v1.ppSlave())
-            || (v2.farPoint() || v2.ppSlave());
+    //         || (v2.farPoint() || v2.ppSlave());
-    }
+    // }
 
 
     // inline void operator=(const Delaunay::Finite_vertices_iterator vit)

tutorials/mesh/cvMesh/flange/Allrun

@@ -8,7 +8,7 @@ cd ${0%/*} || exit 1    # run from this directory
 # mesh from that.
 cp system/controlDict-generatePoints system/controlDict
 cp system/cvMeshDict-generatePoints system/cvMeshDict
-runApplication surfaceFeatureExtract constant/triSurface/flange.obj flange -includedAngle 125
+runApplication surfaceFeatureExtract constant/triSurface/flange.obj flange -includedAngle 150
 runApplication cvMesh
 
 # Use pre-generated aligned points (constant/internalDelaunayVertices)