diff --git a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C index 5157f3c58001cce05f8f33db47eb551df8a7470a..208dba1083342ccfd7639a7e3004658275c13da8 100644 --- a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C +++ b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C @@ -283,7 +283,7 @@ void Foam::conformalVoronoiMesh::insertPoints // ++pit // ) // { - // insertVb(*pit); + // insertPoint(topoint(*pit), Vb::ptInternalPoint); // } label nInserted(number_of_vertices() - preInsertionSize); @@ -318,6 +318,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 @@ -363,7 +367,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(); } @@ -750,6 +753,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()) @@ -770,27 +774,15 @@ void Foam::conformalVoronoiMesh::insertFeaturePoints() << " feature vertices" << endl; } - Info<< "SORT OUT FEATURE POINT DISTRIBUTION AND STORAGE" << endl; - - // featureVertices_.setSize(number_of_vertices()); - - // label featPtI = 0; - - // for - // ( - // Delaunay::Finite_vertices_iterator vit = finite_vertices_begin(); - // vit != finite_vertices_end(); - // vit++ - // ) - // { - // featureVertices_[featPtI] = Vb(vit->point()); + featureVertices_.clear(); - // featureVertices_[featPtI].index() = vit->index(); - - // featureVertices_[featPtI].type() = vit->type(); - - // featPtI++; - // } + forAll(pts, pI) + { + featureVertices_.push_back + ( + Vb(toPoint(pts[pI]), indices[pI], types[pI]) + ); + } constructFeaturePointLocations(); } @@ -1066,14 +1058,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; } @@ -1195,7 +1245,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++; } @@ -1259,7 +1310,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()); @@ -1308,8 +1360,8 @@ bool Foam::conformalVoronoiMesh::distributeBackground() // Remove the entire tessellation reset(); - Info<< "NEED TO MAP FEATURE POINTS" << endl; - // reinsertFeaturePoints(); + // Reinsert feature points, distributing them as necessary. + reinsertFeaturePoints(true); startOfInternalPoints_ = number_of_vertices(); @@ -1326,6 +1378,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]; @@ -1369,6 +1423,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 @@ -1427,6 +1502,8 @@ void Foam::conformalVoronoiMesh::updateSizesAndAlignments ) { storeSizesAndAlignments(storePts); + + timeCheck("Updated sizes and alignments"); } } @@ -1794,6 +1871,8 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh insertBoundingPoints(); + insertFeaturePoints(); + insertInitialPoints(); // Improve the guess that the backgroundMeshDecomposition makes with the @@ -1801,15 +1880,6 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh // better balance the surface conformation load. distributeBackground(); - insertFeaturePoints(); - - Info<< "NEED TO MAP FEATURE POINTS AFTER DISTRIBUTE" << endl; - Info<< "NEED TO ENSURE THAT FEATURE POINTS ARE INSERTED ON THE " - << "CORRECT PROCESSOR" << endl; - - Info<< "NEED TO CHANGE storeSizesAndAlignments AFTER DISTRIBUTE" << endl; - // storeSizesAndAlignments(initPts); - buildSurfaceConformation(rmCoarse); // The introduction of the surface conformation may have distorted the @@ -1822,6 +1892,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); @@ -2174,13 +2259,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(); @@ -2207,14 +2285,20 @@ void Foam::conformalVoronoiMesh::move() updateSizesAndAlignments(pointsToInsert); + buildParallelInterface("move_" + runTime_.timeName()); + + // 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"); } diff --git a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H index ffc87006d6d711c3d433051c1b5fc86b5057d454..8a83f44a416d63feb11cb1c503e7219d2e89dba7 100644 --- a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H +++ b/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(); diff --git a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C index 1b88564fa67810d6d0689becce6535daa58411e7..e129296adca9e9a2dccf590d344263c51db781bd 100644 --- a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C +++ b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C @@ -41,7 +41,11 @@ void Foam::conformalVoronoiMesh::conformToSurface() { // Rebuild, insert and store new surface conformation buildSurfaceConformation(reconfMode); + + storeSurfaceConformation(); } + + // reportSurfaceConformationQuality(); } @@ -411,10 +415,6 @@ void Foam::conformalVoronoiMesh::buildSurfaceConformation << "), stopping iterations" << endl; } } - - // reportSurfaceConformationQuality(); - - storeSurfaceConformation(); } @@ -620,6 +620,8 @@ void Foam::conformalVoronoiMesh::buildParallelInterface receivedVertices, outputName ); + + timeCheck("After buildParallelInterface"); } @@ -1346,6 +1348,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 +1628,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) @@ -1785,6 +1816,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 +1826,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 " diff --git a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H index f620c7cd0841f3948ac8f66157c7a2d0ad8a5246..140c1d7070d3050ec7e5cefb961a7736feb5fb4d 100644 --- a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H +++ b/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 diff --git a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex.H b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex.H index f423d17cf6fa458c9535c2dc18966b0fcad7e3a2..492cf10889867210d3cb2c6cc3553d3414e597af 100644 --- a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex.H +++ b/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 - inline friend bool pointPair - ( - const indexedVertex& v0, - const indexedVertex& v1 - ) - { - return v0.index_ == v1.type_ || v1.index_ == v0.type_; - } + // //- Do the two given vertices consitute a boundary point-pair + // inline friend bool pointPair + // ( + // const indexedVertex& v0, + // const indexedVertex& v1 + // ) + // { + // return v0.index_ == v1.type_ || v1.index_ == v0.type_; + // } - //- Do the three given vertices consitute a boundary triangle - inline friend bool boundaryTriangle - ( - const indexedVertex& v0, - const indexedVertex& v1, - const indexedVertex& v2 - ) - { - return (v0.pairPoint() && pointPair(v1, v2)) - || (v1.pairPoint() && pointPair(v2, v0)) - || (v2.pairPoint() && pointPair(v0, v1)); - } + // //- Do the three given vertices consitute a boundary triangle + // inline friend bool boundaryTriangle + // ( + // const indexedVertex& v0, + // const indexedVertex& v1, + // const indexedVertex& v2 + // ) + // { + // return (v0.pairPoint() && pointPair(v1, v2)) + // || (v1.pairPoint() && pointPair(v2, v0)) + // || (v2.pairPoint() && pointPair(v0, v1)); + // } - //- Do the three given vertices consitute an outside triangle - inline friend bool outsideTriangle - ( - const indexedVertex& v0, - const indexedVertex& v1, - const indexedVertex& v2 - ) - { - return (v0.farPoint() || v0.ppSlave()) - || (v1.farPoint() || v1.ppSlave()) - || (v2.farPoint() || v2.ppSlave()); - } + // //- Do the three given vertices consitute an outside triangle + // inline friend bool outsideTriangle + // ( + // const indexedVertex& v0, + // const indexedVertex& v1, + // const indexedVertex& v2 + // ) + // { + // return (v0.farPoint() || v0.ppSlave()) + // || (v1.farPoint() || v1.ppSlave()) + // || (v2.farPoint() || v2.ppSlave()); + // } // inline void operator=(const Delaunay::Finite_vertices_iterator vit)