diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.C index d78011c9db4ce3fbb157feffd44999abe19cdace..7eb8216b00f1e4033d7f7dc005226394289e4fe4 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.C @@ -122,7 +122,7 @@ Foam::DelaunayMesh<Triangulation>::DelaunayMesh if (indices.headerOk()) { vh->index() = indices[ptI]; - vertexCount()++; + vertexCount_++; } else { diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.H index 1fbf19c1d6957cdc69fabcf678fd13071faac0f8..63d227fccac877066f5230354f56f23dfef8b736 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.H +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.H @@ -201,7 +201,6 @@ public: inline void resetCellCount(); inline label vertexCount() const; - inline label& vertexCount(); inline void resetVertexCount(); diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMeshI.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMeshI.H index bf2277e356da608452c1b3c47cf8db1d90bc2a08..bad25b2b835832d9b301f735631df7d4ef1804a3 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMeshI.H +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMeshI.H @@ -124,12 +124,6 @@ Foam::label Foam::DelaunayMesh<Triangulation>::vertexCount() const return vertexCount_; } -template<class Triangulation> -Foam::label& Foam::DelaunayMesh<Triangulation>::vertexCount() -{ - return vertexCount_; -} - template<class Triangulation> void Foam::DelaunayMesh<Triangulation>::resetVertexCount() diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.C index 37125bb76b4f2cbd8f38f5fce425e6f898d95dc8..a033e5076d54e8180830dda7b9cf55f67c484983 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.C @@ -40,11 +40,15 @@ Foam::cellAspectRatioControl::cellAspectRatioControl aspectRatioDict_.lookupOrDefault<vector> ( "aspectRatioDirection", - vector(0, 0, 0) + vector::zero ) ) { - Info<< nl << "Cell Aspect Ratio Control" << nl + // Normalise the direction + aspectRatioDirection_ /= mag(aspectRatioDirection_) + SMALL; + + Info<< nl + << "Cell Aspect Ratio Control" << nl << " Ratio : " << aspectRatio_ << nl << " Direction : " << aspectRatioDirection_ << endl; @@ -66,22 +70,20 @@ void Foam::cellAspectRatioControl::updateCellSizeAndFaceArea scalar& targetCellSize ) const { - const scalar cosAngle = mag - ( - vectorTools::cosPhi(alignmentDir, aspectRatioDirection_) - ); + const scalar cosAngle = + mag(vectorTools::cosPhi(alignmentDir, aspectRatioDirection_)); // Change target face area based on aspect ratio - targetFaceArea - += targetFaceArea - *(aspectRatio_ - 1.0) - *(1.0 - cosAngle); + targetFaceArea += + targetFaceArea + *(aspectRatio_ - 1.0) + *(1.0 - cosAngle); // Change target cell size based on aspect ratio - targetCellSize - += targetCellSize - *(aspectRatio_ - 1.0) - *cosAngle; + targetCellSize += + targetCellSize + *(aspectRatio_ - 1.0) + *cosAngle; alignmentDir *= 0.5*targetCellSize; } @@ -95,16 +97,15 @@ void Foam::cellAspectRatioControl::updateDeltaVector vector& delta ) const { - const scalar cosAngle = mag - ( - vectorTools::cosPhi(alignmentDir, aspectRatioDirection_) - ); - - delta += 0.5 - *delta - *cosAngle - *(targetCellSize/rABMag) - *(aspectRatio_ - 1.0); + const scalar cosAngle = + mag(vectorTools::cosPhi(alignmentDir, aspectRatioDirection_)); + + delta += + 0.5 + *delta + *cosAngle + *(targetCellSize/rABMag) + *(aspectRatio_ - 1.0); } diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.H index cdc7510f23d4045987bae3d8b64084fcadf1ab26..b3964432b1af9b6384580826b708dffdc51f080e 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.H +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.H @@ -56,7 +56,7 @@ class cellAspectRatioControl const scalar aspectRatio_; - const vector aspectRatioDirection_; + vector aspectRatioDirection_; // Private Member Functions @@ -73,10 +73,7 @@ public: // Constructors //- Construct from dictionary - cellAspectRatioControl - ( - const dictionary& motionDict - ); + cellAspectRatioControl(const dictionary& motionDict); //- Destructor diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C index 62d0c02afe63e1249e0300fd95f54a09e8a18223..3ba5b8c9e733793b5a9b7669bdf84b62a20e3efb 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C @@ -570,159 +570,8 @@ void Foam::conformalVoronoiMesh::buildCellSizeAndAlignmentMesh() } -void Foam::conformalVoronoiMesh::storeSizesAndAlignments() -{ - DynamicList<Point> storePts(number_of_vertices()); - - for - ( - Delaunay::Finite_vertices_iterator vit = finite_vertices_begin(); - vit != finite_vertices_end(); - vit++ - ) - { - if (vit->internalPoint()) - { - storePts.append(vit->point()); - } - } - - storePts.shrink(); - - storeSizesAndAlignments(storePts); -} - - -void Foam::conformalVoronoiMesh::storeSizesAndAlignments -( - const List<Point>& storePts -) -{ -// timeCheck("Start of storeSizesAndAlignments"); -// -// Info << nl << "Store size and alignment" << endl; -// -// sizeAndAlignmentLocations_.setSize(storePts.size()); -// -// storedSizes_.setSize(sizeAndAlignmentLocations_.size()); -// -// storedAlignments_.setSize(sizeAndAlignmentLocations_.size()); -// -// label i = 0; -// -// //checkCellSizing(); -// -// for -// ( -// List<Point>::const_iterator pit = storePts.begin(); -// pit != storePts.end(); -// ++pit -// ) -// { -// pointFromPoint pt = topoint(*pit); -// -//// storedAlignments_[i] = requiredAlignment(pt); -//// -//// storedSizes_[i] = cellShapeControls().cellSize(pt); -// -// cellShapeControls().cellSizeAndAlignment -// ( -// pt, -// storedSizes_[i], -// storedAlignments_[i] -// ); -// -// i++; -// } -// -// timeCheck("Sizes and alignments calculated, build tree"); -// -// buildSizeAndAlignmentTree(); -// -// timeCheck("Size and alignment tree built"); -} - - -void Foam::conformalVoronoiMesh::updateSizesAndAlignments -( - const List<Point>& storePts -) -{ - // This function is only used in serial, the background redistribution - // triggers this when unbalance is detected in parallel. - - if - ( - !Pstream::parRun() - && runTime_.run() - && runTime_.timeIndex() - ) - { - storeSizesAndAlignments(storePts); - - timeCheck("Updated sizes and alignments"); - } -} - - -const Foam::indexedOctree<Foam::treeDataPoint>& -Foam::conformalVoronoiMesh::sizeAndAlignmentTree() const -{ - if (sizeAndAlignmentTreePtr_.empty()) - { - buildSizeAndAlignmentTree(); - } - - return sizeAndAlignmentTreePtr_(); -} - - void Foam::conformalVoronoiMesh::setVertexSizeAndAlignment() { -// Info<< nl << "Looking up target cell alignment and size" << endl; -// -// const indexedOctree<treeDataPoint>& tree = sizeAndAlignmentTree(); -// -// for -// ( -// Delaunay::Finite_vertices_iterator vit = finite_vertices_begin(); -// vit != finite_vertices_end(); -// vit++ -// ) -// { -// if -// ( -// vit->internalOrBoundaryPoint() -// || vit->referredInternalOrBoundaryPoint() -// ) -// { -// pointFromPoint pt = topoint(vit->point()); -// -// pointIndexHit info = tree.findNearest(pt, sqr(GREAT)); -// -// if (info.hit()) -// { -// vit->alignment() = storedAlignments_[info.index()]; -// -// vit->targetCellSize() = storedSizes_[info.index()]; -// } -// else -// { -// WarningIn -// ( -// "void " -// "Foam::conformalVoronoiMesh::setVertexSizeAndAlignment()" -// ) -// << "Point " << pt << " did not find a nearest point " -// << " for alignment and size lookup." << endl; -// -// vit->alignment() = cellShapeControls().cellAlignment(pt); -// -// vit->targetCellSize() = cellShapeControls().cellSize(pt); -// } -// } -// } - Info<< nl << "Calculating target cell alignment and size" << endl; for @@ -744,44 +593,6 @@ void Foam::conformalVoronoiMesh::setVertexSizeAndAlignment() ); } } - -// OFstream str(runTime_.path()/"alignments_internal.obj"); -// -// for -// ( -// Finite_vertices_iterator vit = finite_vertices_begin(); -// vit != finite_vertices_end(); -// ++vit -// ) -// { -// if (!vit->farPoint()) -// { -// // Write alignments -// const tensor& alignment = vit->alignment(); -// pointFromPoint pt = topoint(vit->point()); -// -// if -// ( -// alignment.x() == triad::unset[0] -// || alignment.y() == triad::unset[0] -// || alignment.z() == triad::unset[0] -// ) -// { -// Info<< "Bad alignment = " << vit->info(); -// -// vit->alignment() = tensor::I; -// -// Info<< "New alignment = " << vit->info(); -// -// continue; -// } -// -// meshTools::writeOBJ(str, pt, alignment.x() + pt); -// meshTools::writeOBJ(str, pt, alignment.y() + pt); -// meshTools::writeOBJ(str, pt, alignment.z() + pt); -// } -// } - } @@ -1039,10 +850,6 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh featurePointLocations_(), edgeLocationTreePtr_(), surfacePtLocationTreePtr_(), - sizeAndAlignmentLocations_(), - storedSizes_(), - storedAlignments_(), - sizeAndAlignmentTreePtr_(), surfaceConformationVertices_(), initialPointsMethod_ ( @@ -1133,10 +940,6 @@ void Foam::conformalVoronoiMesh::initialiseForMotion() // (potentially) redistributed. storeSurfaceConformation(); - // Use storeSizesAndAlignments with no feed points because all background - // points may have been distributed. - storeSizesAndAlignments(); - // Report any Delaunay vertices that do not think that they are in the // domain the processor they are on. // reportProcessorOccupancy(); @@ -1836,8 +1639,6 @@ void Foam::conformalVoronoiMesh::move() writeMesh(time().timeName()); } - updateSizesAndAlignments(pointsToInsert); - Info<< nl << "Total displacement = " << totalDisp << nl << "Total distance = " << totalDist << nl diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H index ac11d1a8fabbbc8ec14d1b51c951b4521476fd80..777f277d58a1d7c6a89e8a37e8d700470455f4c0 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H @@ -177,19 +177,6 @@ private: mutable DynamicList<Foam::point> existingSurfacePtLocations_; - //- Store locations where the cell size and alignments will be - // pre-calculated and looked up - pointField sizeAndAlignmentLocations_; - - //- Stored cell size at sizeAndAlignmentLocations_ - scalarField storedSizes_; - - //- Stored alignments at sizeAndAlignmentLocations_ - tensorField storedAlignments_; - - //- Search tree for size and alignment lookup points - mutable autoPtr<indexedOctree<treeDataPoint> > sizeAndAlignmentTreePtr_; - //- Store the surface and feature edge conformation locations to be // reinserted List<Vb> surfaceConformationVertices_; @@ -271,40 +258,6 @@ private: //- Return the local maximum surface protrusion distance inline scalar maxSurfaceProtrusion(const Foam::point& pt) const; - //- Insert Point and return its auto-generated index - inline bool insertPoint - ( - const Point& P, - const indexedVertexEnum::vertexType type - ); - - //- Insert Foam::point with specified index and type - inline bool insertPoint - ( - const Foam::point& p, - const indexedVertexEnum::vertexType type - ); - - //- Insert Point with specified index, type and original processor - inline bool insertReferredPoint - ( - const Point& P, - const label index, - const indexedVertexEnum::vertexType type, - const label processor - ); - - inline bool insertReferredPoint(const Vb& P); - - //- Insert Foam::point with specified index, type and original processor - inline bool insertReferredPoint - ( - const Foam::point& p, - const label index, - const indexedVertexEnum::vertexType type, - const label processor - ); - //- Insert Delaunay vertices using the CGAL range insertion method, // optionally check processor occupancy and distribute to other // processors @@ -539,22 +492,6 @@ private: void buildCellSizeAndAlignmentMesh(); - //- 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 List<Point>& storePts); - - //- Restore the sizes and alignments if required - void updateSizesAndAlignments(const List<Point>& storePts); - - //- Demand driven construction of octree for and alignment points - const indexedOctree<treeDataPoint>& sizeAndAlignmentTree() const; - //- Set the size and alignment data for each vertex void setVertexSizeAndAlignment(); @@ -750,9 +687,6 @@ private: const DynamicList<Foam::point>& existingSurfacePtLocations ) const; - //- Build or rebuild the sizeAndAlignmentTree - void buildSizeAndAlignmentTree() const; - //- Process the surface conformation locations to decide which surface // and edge conformation locations to add void addSurfaceAndEdgeHits @@ -1063,10 +997,6 @@ public: // const List<scalar>& radiusSqrs // ) const; - typedef K::Vector_3 CGALVector; - - inline CGALVector toCGALVector(const Foam::vector& v) const; - // Access diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C index 7c07d12e67c435b97b38d121bf3ffe1fd7fe2d1f..5a45b8ab4f7f0ed22a8b821f1d29b108cfd82338 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C @@ -211,11 +211,6 @@ void Foam::conformalVoronoiMesh::checkDuals() List<Point> duals(number_of_finite_cells()); - typedef CGAL::Exact_predicates_exact_constructions_kernel EK2; - typedef CGAL::Regular_triangulation_euclidean_traits_3<EK2> EK; - typedef CGAL::Cartesian_converter<baseK::Kernel, EK2> To_exact; - typedef CGAL::Cartesian_converter<EK2, baseK::Kernel> Back_from_exact; - // PackedBoolList bPoints(number_of_finite_cells()); // indexDualVertices(duals, bPoints); @@ -1787,7 +1782,6 @@ void Foam::conformalVoronoiMesh::indexDualVertices OBJstream snapping1("snapToSurface1.obj"); OBJstream snapping2("snapToSurface2.obj"); OFstream tetToSnapTo("tetsToSnapTo.obj"); - label offset = 0; for ( diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C index 8be496816f11bbd9eb8cf95b1be6b14eda63b6f3..43feda36545c50ef926888cef64020a36a44df96 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C @@ -76,10 +76,6 @@ void Foam::conformalVoronoiMesh::conformToSurface() { sync(decomposition_().procBounds()); } - - // Use storeSizesAndAlignments with no feed points because all - // background points may have been distributed. - storeSizesAndAlignments(); } // Do not store the surface conformation until after it has been @@ -2055,38 +2051,6 @@ void Foam::conformalVoronoiMesh::buildSurfacePtLocationTree } -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) - ); - - overallBb.min() -= Foam::point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL); - overallBb.max() += Foam::point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL); - - sizeAndAlignmentTreePtr_.reset - ( - new indexedOctree<treeDataPoint> - ( - treeDataPoint(sizeAndAlignmentLocations_), - overallBb, // overall search domain - 10, // max levels - 20.0, // maximum ratio of cubes v.s. cells - 100.0 // max. duplicity; n/a since no bounding boxes. - ) - ); -} - - void Foam::conformalVoronoiMesh::addSurfaceAndEdgeHits ( const Foam::point& vit, diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H index a912515fb311b84588d33f190f6e51ca5b510dc5..9fa34b82d389f50cff6d725762c0e69b8ff15ed2 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H @@ -236,94 +236,6 @@ inline Foam::scalar Foam::conformalVoronoiMesh::maxSurfaceProtrusion } -inline bool Foam::conformalVoronoiMesh::insertPoint -( - const Foam::point& p, - const indexedVertexEnum::vertexType type -) -{ - return insertPoint(toPoint<Point>(p), type); -} - - -inline bool Foam::conformalVoronoiMesh::insertPoint -( - const Point& P, - const indexedVertexEnum::vertexType type -) -{ - uint nVert = number_of_vertices(); - - Vertex_handle vh = insert(P); - - bool pointInserted = true; - - if (nVert == number_of_vertices()) - { - Pout<< "Failed to insert point : " << topoint(P) - << " of type " << type << endl; - pointInserted = false; - } - else - { - vh->index() = getNewVertexIndex(); - vh->type() = type; - } - - return pointInserted; -} - - -inline bool Foam::conformalVoronoiMesh::insertReferredPoint(const Vb& P) -{ - return insertReferredPoint(P.point(), P.index(), P.type(), P.procIndex()); -} - - -inline bool Foam::conformalVoronoiMesh::insertReferredPoint -( - const Foam::point& p, - const label index, - const indexedVertexEnum::vertexType type, - const label processor -) -{ - return insertReferredPoint(toPoint<Point>(p), index, type, processor); -} - - -inline bool Foam::conformalVoronoiMesh::insertReferredPoint -( - const Point& P, - const label index, - const indexedVertexEnum::vertexType type, - const label processor -) -{ - uint nVert = number_of_vertices(); - - Vertex_handle vh = insert(P); - - bool pointInserted = true; - - if (nVert == number_of_vertices()) - { - Pout<< "Failed to insert point " << topoint(P) - << " type: " << type << " index: " << index - << " proc: " << processor << endl; - pointInserted = false; - } - else - { - vh->index() = index; - vh->type() = type; - vh->procIndex() = processor; - } - - return pointInserted; -} - - inline void Foam::conformalVoronoiMesh::createPointPair ( const scalar ppDist, @@ -663,13 +575,6 @@ inline bool Foam::conformalVoronoiMesh::isProcBoundaryEdge // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // -inline Foam::conformalVoronoiMesh::CGALVector -Foam::conformalVoronoiMesh::toCGALVector(const Foam::vector& v) const -{ - return CGALVector(v.x(), v.y(), v.z()); -} - - inline const Foam::Time& Foam::conformalVoronoiMesh::time() const { return runTime_; diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCellI.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCellI.H index 7875d67d392e25ebefbe6d127bf05797112ac9a4..a812fe466455221b262f56ce17a8185d7383ee45 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCellI.H +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCellI.H @@ -108,17 +108,15 @@ int CGAL::indexedCell<Gt, Cb>::cellIndex() const #ifdef CGAL_INEXACT + template<class Gt, class Cb> const Foam::point& CGAL::indexedCell<Gt, Cb>::dual() { - // if (Foam::foamyHexMeshChecks::coplanarTet(*this, 1e-20) == 0) - // { - // Do exact calc - // } - return reinterpret_cast<const Foam::point&>(this->circumcenter()); } + #else + template<class Gt, class Cb> const Foam::point CGAL::indexedCell<Gt, Cb>::dual() { @@ -131,6 +129,7 @@ int CGAL::indexedCell<Gt, Cb>::cellIndex() const CGAL::to_double(P.z()) ); } + #endif diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/pointFeatureEdgesTypes.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/pointFeatureEdgesTypes.H index 3dd3761c59eedad38b32423461aebbf81dea2502..4277f0d3b78acf546414814c0865a533cac26214 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/pointFeatureEdgesTypes.H +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/pointFeatureEdgesTypes.H @@ -25,8 +25,7 @@ Class pointFeatureEdgesTypes Description - struct for holding information on the types of feature edges attached to - feature points + Holds information on the types of feature edges attached to feature points. \*---------------------------------------------------------------------------*/ diff --git a/etc/caseDicts/foamyHexMeshDict b/etc/caseDicts/foamyHexMeshDict index c866c832adc595f601186d31873d2f7469f09410..cb8ac1517bbc62cf2c8ae50b99b2ea73b20f8531 100644 --- a/etc/caseDicts/foamyHexMeshDict +++ b/etc/caseDicts/foamyHexMeshDict @@ -5,12 +5,19 @@ | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + object foamyHexMeshDict; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // #inputMode merge; surfaceConformation { - locationInMesh (0 0 0); pointPairDistanceCoeff 0.1; mixedFeaturePointPPDistanceCoeff 5.0; featurePointExclusionDistanceCoeff 0.65; @@ -22,11 +29,11 @@ surfaceConformation featurePointControls { - specialiseFeaturePoints on; - edgeAiming on; - guardFeaturePoints off; - snapFeaturePoints off; - circulateEdges off; + specialiseFeaturePoints on; + edgeAiming on; + guardFeaturePoints off; + snapFeaturePoints off; + circulateEdges off; } conformationControls @@ -115,20 +122,18 @@ polyMeshFiltering } +backgroundMeshDecomposition +{ + minLevels 1; + sampleResolution 4; + spanScale 20; + maxCellWeightCoeff 20; +} + + meshQualityControls { - maxNonOrtho 65; - maxBoundarySkewness 50; - maxInternalSkewness 10; - maxConcave 80; - minVol -1E30; - minTetQuality 1e-30; - minArea -1; - minTwist 0.02; - minDeterminant 0.001; - minFaceWeight 0.02; - minVolRatio 0.01; - minTriangleTwist -1; + #include "meshQualityDict" } diff --git a/src/OpenFOAM/meshes/meshShapes/face/face.C b/src/OpenFOAM/meshes/meshShapes/face/face.C index 21576535e9aad00966b892337190341a071ef953..41650aad6c16d94bb3d16308ca740dae47871a44 100644 --- a/src/OpenFOAM/meshes/meshShapes/face/face.C +++ b/src/OpenFOAM/meshes/meshShapes/face/face.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -829,4 +829,26 @@ Foam::label Foam::face::trianglesQuads } +Foam::label Foam::longestEdge(const face& f, const pointField& pts) +{ + const edgeList& eds = f.edges(); + + label longestEdgeI = -1; + scalar longestEdgeLength = -SMALL; + + forAll(eds, edI) + { + scalar edgeLength = eds[edI].mag(pts); + + if (edgeLength > longestEdgeLength) + { + longestEdgeI = edI; + longestEdgeLength = edgeLength; + } + } + + return longestEdgeI; +} + + // ************************************************************************* // diff --git a/src/OpenFOAM/meshes/meshShapes/face/face.H b/src/OpenFOAM/meshes/meshShapes/face/face.H index c481258a3cc2d1803cee7239354c32950eec1fc5..adedf991bd79e260bab805a6d5e8fb1e9f2b5af7 100644 --- a/src/OpenFOAM/meshes/meshShapes/face/face.H +++ b/src/OpenFOAM/meshes/meshShapes/face/face.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -408,6 +408,12 @@ public: }; +// Global functions + +//- Find the longest edge on a face. Face point labels index into pts. +label longestEdge(const face& f, const pointField& pts); + + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // } // End namespace Foam diff --git a/src/OpenFOAM/meshes/polyMesh/zones/ZoneMesh/ZoneMesh.C b/src/OpenFOAM/meshes/polyMesh/zones/ZoneMesh/ZoneMesh.C index 33e17117c460f755e423d61dd0ffc31e3acfc113..eeec0974d65c630fd9214bea00922d217f827a28 100644 --- a/src/OpenFOAM/meshes/polyMesh/zones/ZoneMesh/ZoneMesh.C +++ b/src/OpenFOAM/meshes/polyMesh/zones/ZoneMesh/ZoneMesh.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -492,7 +492,7 @@ bool Foam::ZoneMesh<ZoneType, MeshType>::checkParallelSync << " of type " << zones[zoneI].type() << " is not correctly synchronised" << " across coupled boundaries." - << " (coupled faces are either not both " + << " (coupled faces are either not both" << " present in set or have same flipmap)" << endl; } } diff --git a/src/OpenFOAM/meshes/polyMesh/zones/pointZone/pointZone.C b/src/OpenFOAM/meshes/polyMesh/zones/pointZone/pointZone.C index ff15dcea3b247615594dd44c0021f78b225da7e0..30d1161bb7c56f33b4b51289cc64b318fe4d05d7 100644 --- a/src/OpenFOAM/meshes/polyMesh/zones/pointZone/pointZone.C +++ b/src/OpenFOAM/meshes/polyMesh/zones/pointZone/pointZone.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -154,8 +154,15 @@ bool Foam::pointZone::checkParallelSync(const bool report) const forAll(maxZone, pointI) { - // Check point in zone on both sides - if (maxZone[pointI] != minZone[pointI]) + // Check point in same (or no) zone on all processors + if + ( + ( + maxZone[pointI] != -1 + || minZone[pointI] != labelMax + ) + && (maxZone[pointI] != minZone[pointI]) + ) { if (report && !error) { @@ -167,7 +174,8 @@ bool Foam::pointZone::checkParallelSync(const bool report) const << (minZone[pointI] == labelMax ? -1 : minZone[pointI]) << " on some processors and in zone " << maxZone[pointI] - << " on some other processors." + << " on some other processors." << nl + << "(suppressing further warnings)" << endl; } error = true; diff --git a/src/OpenFOAM/primitives/Pair/Pair.H b/src/OpenFOAM/primitives/Pair/Pair.H index 12c16fc6c3fbf42ee0f04dc9508b0c30e9bf05dd..6715bffd7e587631a867fbc25dd115e57e530d6a 100644 --- a/src/OpenFOAM/primitives/Pair/Pair.H +++ b/src/OpenFOAM/primitives/Pair/Pair.H @@ -108,12 +108,6 @@ public: return this->operator[](1); } - //- Return reverse pair - inline Pair<Type> reversePair() const - { - return Pair<Type>(second(), first()); - } - //- Return other inline const Type& other(const Type& a) const { @@ -147,11 +141,11 @@ public: // - -1: same pair, but reversed order static inline int compare(const Pair<Type>& a, const Pair<Type>& b) { - if (a[0] == b[0] && a[1] == b[1]) + if (a == b) { return 1; } - else if (a[0] == b[1] && a[1] == b[0]) + else if (a == reverse(b)) { return -1; } @@ -160,33 +154,64 @@ public: return 0; } } +}; - // Friend Operators +template<class Type> +Pair<Type> reverse(const Pair<Type>& p) +{ + return Pair<Type>(p.second(), p.first()); +} - friend bool operator==(const Pair<Type>& a, const Pair<Type>& b) - { - return (a.first() == b.first() && a.second() == b.second()); - } - friend bool operator!=(const Pair<Type>& a, const Pair<Type>& b) - { - return !(a == b); - } +template<class Type> +bool operator==(const Pair<Type>& a, const Pair<Type>& b) +{ + return (a.first() == b.first() && a.second() == b.second()); +} - friend bool operator<(const Pair<Type>& a, const Pair<Type>& b) - { - return - ( - a.first() < b.first() - || - ( - !(b.first() < a.first()) - && a.second() < b.second() - ) - ); - } -}; + +template<class Type> +bool operator!=(const Pair<Type>& a, const Pair<Type>& b) +{ + return !(a == b); +} + + +template<class Type> +bool operator<(const Pair<Type>& a, const Pair<Type>& b) +{ + return + ( + a.first() < b.first() + || + ( + !(b.first() < a.first()) + && a.second() < b.second() + ) + ); +} + + +template<class Type> +bool operator<=(const Pair<Type>& a, const Pair<Type>& b) +{ + return !(b < a); +} + + +template<class Type> +bool operator>(const Pair<Type>& a, const Pair<Type>& b) +{ + return (b < a); +} + + +template<class Type> +bool operator>=(const Pair<Type>& a, const Pair<Type>& b) +{ + return !(a < b); +} // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/src/OpenFOAM/primitives/Tuple2/Tuple2.H b/src/OpenFOAM/primitives/Tuple2/Tuple2.H index b7aae7598ce719d7037840c5f4527358c8dc91e1..26ca000575bba9c812ea540b3cf6b5fd76592ca4 100644 --- a/src/OpenFOAM/primitives/Tuple2/Tuple2.H +++ b/src/OpenFOAM/primitives/Tuple2/Tuple2.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -47,20 +47,6 @@ namespace Foam template<class Type1, class Type2> class Tuple2; -template<class Type1, class Type2> -inline bool operator== -( - const Tuple2<Type1, Type2>&, - const Tuple2<Type1, Type2>& -); - -template<class Type1, class Type2> -inline bool operator!= -( - const Tuple2<Type1, Type2>&, - const Tuple2<Type1, Type2>& -); - template<class Type1, class Type2> inline Istream& operator>>(Istream&, Tuple2<Type1, Type2>&); @@ -129,27 +115,6 @@ public: return s_; } - //- Return reverse pair - inline Tuple2<Type2, Type1> reverseTuple2() const - { - return Tuple2<Type2, Type1>(second(), first()); - } - - - // Friend Operators - - friend bool operator== <Type1, Type2> - ( - const Tuple2<Type1, Type2>& a, - const Tuple2<Type1, Type2>& b - ); - - friend bool operator!= <Type1, Type2> - ( - const Tuple2<Type1, Type2>& a, - const Tuple2<Type1, Type2>& b - ); - // IOstream operators @@ -169,6 +134,14 @@ public: }; +//- Return reverse of a tuple2 +template<class Type1, class Type2> +inline Tuple2<Type2, Type1> reverse(const Tuple2<Type1, Type2>& t) +{ + return Tuple2<Type2, Type1>(t.second(), t.first()); +} + + template<class Type1, class Type2> inline bool operator== ( diff --git a/src/OpenFOAM/primitives/triad/triad.C b/src/OpenFOAM/primitives/triad/triad.C index ac997d8bd27f49ddce3957bceb1a9dcbd8435c3d..a642463f422c7407e6f517aa2516321c0eb29fe9 100644 --- a/src/OpenFOAM/primitives/triad/triad.C +++ b/src/OpenFOAM/primitives/triad/triad.C @@ -282,6 +282,11 @@ void Foam::triad::align(const vector& v) Foam::triad Foam::triad::sortxyz() const { + if (!this->set()) + { + return *this; + } + triad t; if diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.C b/src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.C index 70c2e7c558cafe5179d88c91ab29cf28330a6b99..6ea1e748b052e3205fc31ec17f2b2b2faf24ebde 100644 --- a/src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.C +++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.C @@ -42,32 +42,6 @@ defineTypeNameAndDebug(edgeCollapser, 0); } -Foam::label Foam::edgeCollapser::longestEdge -( - const face& f, - const pointField& pts -) -{ - const edgeList& eds = f.edges(); - - label longestEdgeI = -1; - scalar longestEdgeLength = -SMALL; - - forAll(eds, edI) - { - scalar edgeLength = eds[edI].mag(pts); - - if (edgeLength > longestEdgeLength) - { - longestEdgeI = edI; - longestEdgeLength = edgeLength; - } - } - - return longestEdgeI; -} - - Foam::HashSet<Foam::label> Foam::edgeCollapser::checkBadFaces ( const polyMesh& mesh, diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.H b/src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.H index 2393458cf77be0add69297f9589e5b71d8b9e180..dab8bd298b0becb2bdae6c9c9b9ccd58e4e45d3e 100644 --- a/src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.H +++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.H @@ -262,9 +262,6 @@ public: // Check - //- Find the longest edge in a face - static label longestEdge(const face& f, const pointField& pts); - //- Calls motionSmoother::checkMesh and returns a set of bad faces static HashSet<label> checkBadFaces ( diff --git a/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementBaffles.C b/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementBaffles.C index 6b9d1074c55438384b3e5af61b06033fe576df6c..16528d3ea781c8057f3a68bf7be6a791c8cbd44f 100644 --- a/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementBaffles.C +++ b/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementBaffles.C @@ -450,7 +450,7 @@ Foam::Map<Foam::labelPair> Foam::meshRefinement::getZoneBafflePatches labelPair patches = zPatches; if (fZone.flipMap()[i]) { - patches = patches.reversePair(); + patches = reverse(patches); } if (!bafflePatch.insert(faceI, patches)) diff --git a/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C b/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C index 98b16e9f6b852e005aa4294a00e579d6f8572c69..f0277b2990ffa1cc2ed6dd89db768f70b02ae05a 100644 --- a/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C +++ b/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C @@ -605,8 +605,10 @@ void Foam::refinementSurfaces::setMinLevelFields // searchableBox (size=6) if (geom.regions().size() > 1 && geom.globalSize() > 10) { - // Representative local coordinates - const pointField ctrs(geom.coordinates()); + // Representative local coordinates and bounding sphere + pointField ctrs; + scalarField radiusSqr; + geom.boundingSpheres(ctrs, radiusSqr); labelList minLevelField(ctrs.size(), -1); { @@ -614,12 +616,7 @@ void Foam::refinementSurfaces::setMinLevelFields // distributed surface where local indices differ from global // ones (needed for getRegion call) List<pointIndexHit> info; - geom.findNearest - ( - ctrs, - scalarField(ctrs.size(), sqr(GREAT)), - info - ); + geom.findNearest(ctrs, radiusSqr, info); // Get per element the region labelList region; @@ -628,7 +625,7 @@ void Foam::refinementSurfaces::setMinLevelFields // From the region get the surface-wise refinement level forAll(minLevelField, i) { - if (info[i].hit()) + if (info[i].hit()) //Note: should not be necessary { minLevelField[i] = minLevel(surfI, region[i]); } diff --git a/src/meshTools/searchableSurface/searchableBox.C b/src/meshTools/searchableSurface/searchableBox.C index 7d639b99f8697205b13c4d6159a395b6630c65a9..5970e8257597703f908c667596a6a4ec8eecb8fc 100644 --- a/src/meshTools/searchableSurface/searchableBox.C +++ b/src/meshTools/searchableSurface/searchableBox.C @@ -249,6 +249,41 @@ Foam::tmp<Foam::pointField> Foam::searchableBox::coordinates() const } +void Foam::searchableBox::boundingSpheres +( + pointField& centres, + scalarField& radiusSqr +) const +{ + centres.setSize(size()); + radiusSqr.setSize(size()); + radiusSqr = 0.0; + + const pointField pts(treeBoundBox::points()); + const faceList& fcs = treeBoundBox::faces; + + forAll(fcs, i) + { + const face& f = fcs[i]; + + centres[i] = f.centre(pts); + forAll(f, fp) + { + const point& pt = pts[f[fp]]; + + radiusSqr[i] = Foam::max + ( + radiusSqr[i], + Foam::magSqr(pt-centres[i]) + ); + } + } + + // Add a bit to make sure all points are tested inside + radiusSqr += Foam::sqr(SMALL); +} + + Foam::tmp<Foam::pointField> Foam::searchableBox::points() const { return treeBoundBox::points(); diff --git a/src/meshTools/searchableSurface/searchableBox.H b/src/meshTools/searchableSurface/searchableBox.H index 40e5c5d7110ff3fe5a5c5a79cdd6b72c2c66e73e..d677c1b87de1f59b4c884d1b1fd74c8638095450 100644 --- a/src/meshTools/searchableSurface/searchableBox.H +++ b/src/meshTools/searchableSurface/searchableBox.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -129,6 +129,14 @@ public: // Usually the element centres (should be of length size()). virtual tmp<pointField> coordinates() const; + //- Get bounding spheres (centre and radius squared), one per element. + // Any point on element is guaranteed to be inside. + virtual void boundingSpheres + ( + pointField& centres, + scalarField& radiusSqr + ) const; + //- Get the points that define the surface. virtual tmp<pointField> points() const; diff --git a/src/meshTools/searchableSurface/searchableCylinder.C b/src/meshTools/searchableSurface/searchableCylinder.C index 70d674d521b580acae0845f51890be762d07ce73..7657e883c98acb6c2163fde2ac20d664a94846da 100644 --- a/src/meshTools/searchableSurface/searchableCylinder.C +++ b/src/meshTools/searchableSurface/searchableCylinder.C @@ -47,6 +47,23 @@ Foam::tmp<Foam::pointField> Foam::searchableCylinder::coordinates() const } +void Foam::searchableCylinder::boundingSpheres +( + pointField& centres, + scalarField& radiusSqr +) const +{ + centres.setSize(1); + centres[0] = 0.5*(point1_ + point2_); + + radiusSqr.setSize(1); + radiusSqr[0] = Foam::magSqr(point1_-centres[0]) + Foam::sqr(radius_); + + // Add a bit to make sure all points are tested inside + radiusSqr += Foam::sqr(SMALL); +} + + Foam::tmp<Foam::pointField> Foam::searchableCylinder::points() const { tmp<pointField> tPts(new pointField(2)); diff --git a/src/meshTools/searchableSurface/searchableCylinder.H b/src/meshTools/searchableSurface/searchableCylinder.H index 0bb570e4806b5ee41ae936bd12d05aa6381a39d7..ccf6850dc1be487f7be92f8e02dbb6a9ed1a62c5 100644 --- a/src/meshTools/searchableSurface/searchableCylinder.H +++ b/src/meshTools/searchableSurface/searchableCylinder.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -154,6 +154,14 @@ public: // Usually the element centres (should be of length size()). virtual tmp<pointField> coordinates() const; + //- Get bounding spheres (centre and radius squared), one per element. + // Any point on element is guaranteed to be inside. + virtual void boundingSpheres + ( + pointField& centres, + scalarField& radiusSqr + ) const; + //- Get the points that define the surface. virtual tmp<pointField> points() const; diff --git a/src/meshTools/searchableSurface/searchablePlane.C b/src/meshTools/searchableSurface/searchablePlane.C index a76821d45740e53e4b571668977360b6bf229481..e6c9668b1f134558c972d7301a29afe64d648cab 100644 --- a/src/meshTools/searchableSurface/searchablePlane.C +++ b/src/meshTools/searchableSurface/searchablePlane.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -134,6 +134,20 @@ const Foam::wordList& Foam::searchablePlane::regions() const } +void Foam::searchablePlane::boundingSpheres +( + pointField& centres, + scalarField& radiusSqr +) const +{ + centres.setSize(1); + centres[0] = refPoint(); + + radiusSqr.setSize(1); + radiusSqr[0] = Foam::sqr(GREAT); +} + + void Foam::searchablePlane::findNearest ( const pointField& samples, diff --git a/src/meshTools/searchableSurface/searchablePlane.H b/src/meshTools/searchableSurface/searchablePlane.H index cc64844d0df9b9da1844ecae8e94057901dbde2c..a79f93f90292ae931e00d15fcb42c1c199a69733 100644 --- a/src/meshTools/searchableSurface/searchablePlane.H +++ b/src/meshTools/searchableSurface/searchablePlane.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -130,6 +130,15 @@ public: return tCtrs; } + //- Get bounding spheres (centre and radius squared), one per element. + // Any point on element is guaranteed to be inside. + // Note: radius limited to sqr(GREAT) + virtual void boundingSpheres + ( + pointField& centres, + scalarField& radiusSqr + ) const; + //- Get the points that define the surface. virtual tmp<pointField> points() const { diff --git a/src/meshTools/searchableSurface/searchablePlate.C b/src/meshTools/searchableSurface/searchablePlate.C index 37404db05117444141a9a024425f66ad2816ac78..0ed10bd675c123ca4d77374a0d8546c85646d149 100644 --- a/src/meshTools/searchableSurface/searchablePlate.C +++ b/src/meshTools/searchableSurface/searchablePlate.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -274,6 +274,71 @@ const Foam::wordList& Foam::searchablePlate::regions() const } +Foam::tmp<Foam::pointField> Foam::searchablePlate::coordinates() const +{ + return tmp<pointField>(new pointField(1, origin_ + 0.5*span_)); +} + + +void Foam::searchablePlate::boundingSpheres +( + pointField& centres, + scalarField& radiusSqr +) const +{ + centres.setSize(1); + centres[0] = origin_ + 0.5*span_; + + radiusSqr.setSize(1); + radiusSqr[0] = Foam::magSqr(0.5*span_); + + // Add a bit to make sure all points are tested inside + radiusSqr += Foam::sqr(SMALL); +} + + +Foam::tmp<Foam::pointField> Foam::searchablePlate::points() const +{ + tmp<pointField> tPts(new pointField(4)); + pointField& pts = tPts(); + + pts[0] = origin_; + pts[2] = origin_ + span_; + + if (span_.x() < span_.y() && span_.x() < span_.z()) + { + pts[1] = origin_ + point(0, span_.y(), 0); + pts[3] = origin_ + point(0, 0, span_.z()); + } + else if (span_.y() < span_.z()) + { + pts[1] = origin_ + point(span_.x(), 0, 0); + pts[3] = origin_ + point(0, 0, span_.z()); + } + else + { + pts[1] = origin_ + point(span_.x(), 0, 0); + pts[3] = origin_ + point(0, span_.y(), 0); + } + + return tPts; +} + + +bool Foam::searchablePlate::overlaps(const boundBox& bb) const +{ + return + ( + (origin_.x() + span_.x()) >= bb.min().x() + && origin_.x() <= bb.max().x() + && (origin_.y() + span_.y()) >= bb.min().y() + && origin_.y() <= bb.max().y() + && (origin_.z() + span_.z()) >= bb.min().z() + && origin_.z() <= bb.max().z() + ); +} + + void Foam::searchablePlate::findNearest ( const pointField& samples, diff --git a/src/meshTools/searchableSurface/searchablePlate.H b/src/meshTools/searchableSurface/searchablePlate.H index 86b0785077370effb75ecb83c13091d2506c6fcf..af2d732ff635ecbc0aaa4f1ce3272138568e08aa 100644 --- a/src/meshTools/searchableSurface/searchablePlate.H +++ b/src/meshTools/searchableSurface/searchablePlate.H @@ -145,53 +145,21 @@ public: //- Get representative set of element coordinates // Usually the element centres (should be of length size()). - virtual tmp<pointField> coordinates() const - { - tmp<pointField> tCtrs(new pointField(1, origin_ + 0.5*span_)); - return tCtrs; - } - - //- Get the points that define the surface. - virtual tmp<pointField> points() const - { - tmp<pointField> tPts(new pointField(4)); - pointField& pts = tPts(); + virtual tmp<pointField> coordinates() const; - pts[0] = origin_; - pts[2] = origin_ + span_; - - if (span_.x() < span_.y() && span_.x() < span_.z()) - { - pts[1] = origin_ + point(0, span_.y(), 0); - pts[3] = origin_ + point(0, 0, span_.z()); - } - else if (span_.y() < span_.z()) - { - pts[1] = origin_ + point(span_.x(), 0, 0); - pts[3] = origin_ + point(0, 0, span_.z()); - } - else - { - pts[1] = origin_ + point(span_.x(), 0, 0); - pts[3] = origin_ + point(0, span_.y(), 0); - } + //- Get bounding spheres (centre and radius squared), one per element. + // Any point on element is guaranteed to be inside. + virtual void boundingSpheres + ( + pointField& centres, + scalarField& radiusSqr + ) const; - return tPts; - } + //- Get the points that define the surface. + virtual tmp<pointField> points() const; //- Does any part of the surface overlap the supplied bound box? - virtual bool overlaps(const boundBox& bb) const - { - return - ( - (origin_.x() + span_.x()) >= bb.min().x() - && origin_.x() <= bb.max().x() - && (origin_.y() + span_.y()) >= bb.min().y() - && origin_.y() <= bb.max().y() - && (origin_.z() + span_.z()) >= bb.min().z() - && origin_.z() <= bb.max().z() - ); - } + virtual bool overlaps(const boundBox& bb) const; // Multiple point queries. diff --git a/src/meshTools/searchableSurface/searchableSphere.C b/src/meshTools/searchableSurface/searchableSphere.C index c367c6281babd1d176a5845146cdce897d6b6c93..e4eae4ab27452e40e4dcb1977754ebc5b13d9f1f 100644 --- a/src/meshTools/searchableSurface/searchableSphere.C +++ b/src/meshTools/searchableSurface/searchableSphere.C @@ -185,6 +185,23 @@ const Foam::wordList& Foam::searchableSphere::regions() const +void Foam::searchableSphere::boundingSpheres +( + pointField& centres, + scalarField& radiusSqr +) const +{ + centres.setSize(1); + centres[0] = centre_; + + radiusSqr.setSize(1); + radiusSqr[0] = Foam::sqr(radius_); + + // Add a bit to make sure all points are tested inside + radiusSqr += Foam::sqr(SMALL); +} + + void Foam::searchableSphere::findNearest ( const pointField& samples, diff --git a/src/meshTools/searchableSurface/searchableSphere.H b/src/meshTools/searchableSurface/searchableSphere.H index 5d46f8ad4b5128e3d6748e2a8d3009042a6d8ab0..50900f936185c9a884fca3ef3d338452d2be1aee 100644 --- a/src/meshTools/searchableSurface/searchableSphere.H +++ b/src/meshTools/searchableSurface/searchableSphere.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -60,7 +60,7 @@ private: //- Centre point const point centre_; - //- Radius squared + //- Radius const scalar radius_; //- Names of regions @@ -139,6 +139,14 @@ public: return tCtrs; } + //- Get bounding spheres (centre and radius squared), one per element. + // Any point on element is guaranteed to be inside. + virtual void boundingSpheres + ( + pointField& centres, + scalarField& radiusSqr + ) const; + //- Get the points that define the surface. virtual tmp<pointField> points() const { diff --git a/src/meshTools/searchableSurface/searchableSurface.H b/src/meshTools/searchableSurface/searchableSurface.H index 90224ae9ce54c517b3cdf589467b3533607d930e..6a2adcc9466561b1366d6b642a15c5e6ee56cf28 100644 --- a/src/meshTools/searchableSurface/searchableSurface.H +++ b/src/meshTools/searchableSurface/searchableSurface.H @@ -190,6 +190,14 @@ public: // Usually the element centres (should be of length size()). virtual tmp<pointField> coordinates() const = 0; + //- Get bounding spheres (centre and radius squared), one per element. + // Any point on element is guaranteed to be inside. + virtual void boundingSpheres + ( + pointField& centres, + scalarField& radiusSqr + ) const = 0; + //- Get the points that define the surface. virtual tmp<pointField> points() const = 0; diff --git a/src/meshTools/searchableSurface/searchableSurfaceCollection.C b/src/meshTools/searchableSurface/searchableSurfaceCollection.C index 3a02066afae137a1d251eaf982bd5db3cf3a82da..3aad85f72d7bfd552b681aa6cb1f359bd25e89f7 100644 --- a/src/meshTools/searchableSurface/searchableSurfaceCollection.C +++ b/src/meshTools/searchableSurface/searchableSurfaceCollection.C @@ -354,6 +354,42 @@ Foam::searchableSurfaceCollection::coordinates() const } +void Foam::searchableSurfaceCollection::boundingSpheres +( + pointField& centres, + scalarField& radiusSqr +) const +{ + centres.setSize(size()); + radiusSqr.setSize(centres.size()); + + // Append individual coordinates + label coordI = 0; + + forAll(subGeom_, surfI) + { + scalar maxScale = cmptMax(scale_[surfI]); + + pointField subCentres; + scalarField subRadiusSqr; + subGeom_[surfI].boundingSpheres(subCentres, subRadiusSqr); + + forAll(subCentres, i) + { + centres[coordI++] = transform_[surfI].globalPosition + ( + cmptMultiply + ( + subCentres[i], + scale_[surfI] + ) + ); + radiusSqr[coordI++] = maxScale*subRadiusSqr[i]; + } + } +} + + Foam::tmp<Foam::pointField> Foam::searchableSurfaceCollection::points() const { diff --git a/src/meshTools/searchableSurface/searchableSurfaceCollection.H b/src/meshTools/searchableSurface/searchableSurfaceCollection.H index 53e677eaa656ce4adcf6b5b8c8119bd55f0f441d..dd8a21664b9823f4f27d52dd24ea8866c96e438d 100644 --- a/src/meshTools/searchableSurface/searchableSurfaceCollection.H +++ b/src/meshTools/searchableSurface/searchableSurfaceCollection.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -176,6 +176,14 @@ public: // Usually the element centres (should be of length size()). virtual tmp<pointField> coordinates() const; + //- Get bounding spheres (centre and radius squared), one per element. + // Any point on element is guaranteed to be inside. + virtual void boundingSpheres + ( + pointField& centres, + scalarField& radiusSqr + ) const; + //- Get the points that define the surface. virtual tmp<pointField> points() const; diff --git a/src/meshTools/searchableSurface/searchableSurfaceWithGaps.H b/src/meshTools/searchableSurface/searchableSurfaceWithGaps.H index 1344cda2d9ed5fc4da35dcea0e3ca7a57e1284e2..d07d331c8a56a442c5b1b9c435205aebf98f90a6 100644 --- a/src/meshTools/searchableSurface/searchableSurfaceWithGaps.H +++ b/src/meshTools/searchableSurface/searchableSurfaceWithGaps.H @@ -156,6 +156,17 @@ public: return surface().coordinates(); } + //- Get bounding spheres (centre and radius squared), one per element. + // Any point on element is guaranteed to be inside. + virtual void boundingSpheres + ( + pointField& centres, + scalarField& radiusSqr + ) const + { + surface().boundingSpheres(centres, radiusSqr); + } + //- Get the points that define the surface. virtual tmp<pointField> points() const { diff --git a/src/meshTools/searchableSurface/triSurfaceMesh.C b/src/meshTools/searchableSurface/triSurfaceMesh.C index 125450603b06524ee9b9509224cffe2cc2cf1fa9..7b6b8ad6ad9785eb428197e4a10b2062993e74d5 100644 --- a/src/meshTools/searchableSurface/triSurfaceMesh.C +++ b/src/meshTools/searchableSurface/triSurfaceMesh.C @@ -145,10 +145,12 @@ bool Foam::triSurfaceMesh::addFaceToEdge bool Foam::triSurfaceMesh::isSurfaceClosed() const { + const pointField& pts = triSurface::points(); + // Construct pointFaces. Let's hope surface has compact point // numbering ... labelListList pointFaces; - invertManyToMany(points()().size(), *this, pointFaces); + invertManyToMany(pts.size(), *this, pointFaces); // Loop over all faces surrounding point. Count edges emanating from point. // Every edge should be used by two faces exactly. @@ -241,7 +243,9 @@ Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io, const triSurface& s) minQuality_(-1), surfaceClosed_(-1) { - bounds() = boundBox(points()); + const pointField& pts = triSurface::points(); + + bounds() = boundBox(pts); } @@ -287,7 +291,9 @@ Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io) minQuality_(-1), surfaceClosed_(-1) { - bounds() = boundBox(points()); + const pointField& pts = triSurface::points(); + + bounds() = boundBox(pts); } @@ -347,7 +353,9 @@ Foam::triSurfaceMesh::triSurfaceMesh triSurface::scalePoints(scaleFactor); } - bounds() = boundBox(points()); + const pointField& pts = triSurface::points(); + + bounds() = boundBox(pts); // Have optional minimum quality for normal calculation if (dict.readIfPresent("minQuality", minQuality_) && minQuality_ > 0) @@ -393,6 +401,34 @@ Foam::tmp<Foam::pointField> Foam::triSurfaceMesh::coordinates() const } +void Foam::triSurfaceMesh::boundingSpheres +( + pointField& centres, + scalarField& radiusSqr +) const +{ + centres = coordinates(); + radiusSqr.setSize(size()); + radiusSqr = 0.0; + + const pointField& pts = triSurface::points(); + + forAll(*this, faceI) + { + const labelledTri& f = triSurface::operator[](faceI); + const point& fc = centres[faceI]; + forAll(f, fp) + { + const point& pt = pts[f[fp]]; + radiusSqr[faceI] = max(radiusSqr[faceI], Foam::magSqr(fc-pt)); + } + } + + // Add a bit to make sure all points are tested inside + radiusSqr += Foam::sqr(SMALL); +} + + Foam::tmp<Foam::pointField> Foam::triSurfaceMesh::points() const { return triSurface::points(); @@ -606,6 +642,7 @@ void Foam::triSurfaceMesh::getNormal ) const { const triSurface& s = static_cast<const triSurface&>(*this); + const pointField& pts = s.points(); normal.setSize(info.size()); @@ -621,9 +658,9 @@ void Foam::triSurfaceMesh::getNormal if (info[i].hit()) { label faceI = info[i].index(); - normal[i] = s[faceI].normal(points()); + normal[i] = s[faceI].normal(pts); - scalar qual = s[faceI].tri(points()).quality(); + scalar qual = s[faceI].tri(pts).quality(); if (qual < minQuality_) { @@ -633,11 +670,11 @@ void Foam::triSurfaceMesh::getNormal forAll(fFaces, j) { label nbrI = fFaces[j]; - scalar nbrQual = s[nbrI].tri(points()).quality(); + scalar nbrQual = s[nbrI].tri(pts).quality(); if (nbrQual > qual) { qual = nbrQual; - normal[i] = s[nbrI].normal(points()); + normal[i] = s[nbrI].normal(pts); } } } @@ -662,7 +699,7 @@ void Foam::triSurfaceMesh::getNormal //normal[i] = faceNormals()[faceI]; //- Uncached - normal[i] = s[faceI].normal(points()); + normal[i] = s[faceI].normal(pts); normal[i] /= mag(normal[i]) + VSMALL; } else diff --git a/src/meshTools/searchableSurface/triSurfaceMesh.H b/src/meshTools/searchableSurface/triSurfaceMesh.H index fa7e9e511c769a0d37b7eb6d9b9a885ecaa367cd..e6ecd39594aa03cb67c7e97be15ab93203cdf9a4 100644 --- a/src/meshTools/searchableSurface/triSurfaceMesh.H +++ b/src/meshTools/searchableSurface/triSurfaceMesh.H @@ -187,6 +187,14 @@ public: // Usually the element centres (should be of length size()). virtual tmp<pointField> coordinates() const; + //- Get bounding spheres (centre and radius squared). Any point + // on surface is guaranteed to be inside. + virtual void boundingSpheres + ( + pointField& centres, + scalarField& radiusSqr + ) const; + //- Get the points that define the surface. virtual tmp<pointField> points() const;