diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/controlMeshRefinement/controlMeshRefinement.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/controlMeshRefinement/controlMeshRefinement.C index b100d8eebb..8d954799a1 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/controlMeshRefinement/controlMeshRefinement.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/controlMeshRefinement/controlMeshRefinement.C @@ -295,14 +295,14 @@ void Foam::controlMeshRefinement::initialMeshPopulation shapeController_.minimumCellSize() ); } - else if (maxPriority == controlFunction.maxPriority()) - { - vertices[vI].targetCellSize() = max - ( - min(sizes[vI], size), - shapeController_.minimumCellSize() - ); - } +// else if (maxPriority == controlFunction.maxPriority()) +// { +// vertices[vI].targetCellSize() = max +// ( +// min(sizes[vI], size), +// shapeController_.minimumCellSize() +// ); +// } else { vertices[vI].targetCellSize() = max diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C index e580febdd8..8224492921 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C @@ -1487,8 +1487,13 @@ void Foam::conformalVoronoiMesh::move() // Save displacements to file. if (foamyHexMeshControls().objOutput() && time().outputTime()) { - Pout<< "Writing point displacement vectors to file." << endl; - OFstream str("displacements_" + runTime_.timeName() + ".obj"); + Info<< "Writing point displacement vectors to file." << endl; + OFstream str + ( + time().path() + + "displacements_" + runTime_.timeName() + + ".obj" + ); for ( diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C index db77f62f5f..f5a023592c 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C @@ -550,53 +550,130 @@ void Foam::conformalVoronoiMesh::buildSurfaceConformation() } } +// for +// ( +// Delaunay::Finite_edges_iterator eit = finite_edges_begin(); +// eit != finite_edges_end(); +// ++eit +// ) +// { +// Cell_handle c = eit->first; +// Vertex_handle vA = c->vertex(eit->second); +// Vertex_handle vB = c->vertex(eit->third); +// +// if +// ( +// vA->referred() +// || vB->referred() +// ) +// { +// continue; +// } +// +// if +// ( +// (vA->internalPoint() && vA->referred()) +// || (vB->internalPoint() && vB->referred()) +// ) +// { +// continue; +// } +// +// if +// ( +// (vA->internalPoint() && vB->externalBoundaryPoint()) +// || (vB->internalPoint() && vA->externalBoundaryPoint()) +// || (vA->internalBoundaryPoint() && vB->internalBoundaryPoint()) +// ) +// { +// pointIndexHitAndFeatureDynList surfaceIntersections(0.5*AtoV); +// pointIndexHit surfHit; +// label hitSurface; +// +// geometryToConformTo_.findSurfaceNearest +// ( +// ( +// vA->internalPoint() +// ? topoint(vA->point()) +// : topoint(vB->point()) +// ), +// magSqr(topoint(vA->point()) - topoint(vB->point())), +// surfHit, +// hitSurface +// ); +// +// if (surfHit.hit()) +// { +// surfaceIntersections.append +// ( +// pointIndexHitAndFeature(surfHit, hitSurface) +// ); +// +// addSurfaceAndEdgeHits +// ( +// ( +// vA->internalPoint() +// ? topoint(vA->point()) +// : topoint(vB->point()) +// ), +// surfaceIntersections, +// surfacePtReplaceDistCoeffSqr, +// edgeSearchDistCoeffSqr, +// surfaceHits, +// featureEdgeHits, +// surfaceToTreeShape, +// edgeToTreeShape, +// surfacePtToEdgePtDist, +// false +// ); +// } +// } +// } + for ( - Delaunay::Finite_edges_iterator eit = finite_edges_begin(); - eit != finite_edges_end(); - ++eit + Delaunay::Finite_cells_iterator cit = finite_cells_begin(); + cit != finite_cells_end(); + ++cit ) { - Cell_handle c = eit->first; - Vertex_handle vA = c->vertex(eit->second); - Vertex_handle vB = c->vertex(eit->third); - - if - ( - vA->referred() - || vB->referred() - ) + label nInternal = 0; + for (label vI = 0; vI < 4; vI++) { - continue; + if (cit->vertex(vI)->internalPoint()) + { + nInternal++; + } } - if - ( - (vA->internalPoint() && vA->referred()) - || (vB->internalPoint() && vB->referred()) - ) + if (nInternal == 1 && cit->hasBoundaryPoint()) + //if (cit->boundaryDualVertex() && !cit->hasReferredPoint()) { - continue; - } + const Foam::point& pt = cit->dual(); + + const Foam::point cellCentre = + topoint + ( + CGAL::centroid + ( + CGAL::Tetrahedron_3 + ( + cit->vertex(0)->point(), + cit->vertex(1)->point(), + cit->vertex(2)->point(), + cit->vertex(3)->point() + ) + ) + ); - if - ( - (vA->internalPoint() && vB->externalBoundaryPoint()) - || (vB->internalPoint() && vA->externalBoundaryPoint()) - ) - { pointIndexHitAndFeatureDynList surfaceIntersections(0.5*AtoV); pointIndexHit surfHit; label hitSurface; - geometryToConformTo_.findSurfaceNearest + geometryToConformTo_.findSurfaceNearestIntersection ( - ( - vA->internalPoint() - ? topoint(vA->point()) - : topoint(vB->point()) - ), - magSqr(topoint(vA->point()) - topoint(vB->point())), + cellCentre, + pt, surfHit, hitSurface ); @@ -610,11 +687,7 @@ void Foam::conformalVoronoiMesh::buildSurfaceConformation() addSurfaceAndEdgeHits ( - ( - vA->internalPoint() - ? topoint(vA->point()) - : topoint(vB->point()) - ), + pt, surfaceIntersections, surfacePtReplaceDistCoeffSqr, edgeSearchDistCoeffSqr, @@ -1613,6 +1686,7 @@ void Foam::conformalVoronoiMesh::limitDisplacement // Do not allow infinite recursion if (callCount > 7) { + displacement = vector::zero; return; } @@ -1661,6 +1735,7 @@ void Foam::conformalVoronoiMesh::limitDisplacement if (magSqr(pt - surfHit.hitPoint()) <= searchDistanceSqr) { // Cannot limit displacement, point closer than tolerance + displacement = vector::zero; return; } } @@ -1748,17 +1823,19 @@ bool Foam::conformalVoronoiMesh::nearSurfacePoint const Foam::point& pt = pHit.first().hitPoint(); pointIndexHit closePoint; - const bool closeToSurfacePt = pointIsNearSurfaceLocation(pt, closePoint); if ( closeToSurfacePt - && mag(pt - closePoint.hitPoint()) > pointPairDistance(pt) + && ( + magSqr(pt - closePoint.hitPoint()) + > sqr(pointPairDistance(pt)) + ) ) { - const scalar cosAngle - = angleBetweenSurfacePoints(pt, closePoint.hitPoint()); + const scalar cosAngle = + angleBetweenSurfacePoints(pt, closePoint.hitPoint()); // @todo make this tolerance run-time selectable? if (cosAngle < searchAngleOppositeSurface) @@ -1768,11 +1845,6 @@ bool Foam::conformalVoronoiMesh::nearSurfacePoint const scalar searchDist = targetCellSize(closePoint.hitPoint()); - if (searchDist < SMALL) - { - Pout<< "WARNING: SMALL CELL SIZE" << endl; - } - geometryToConformTo_.findSurfaceNearest ( closePoint.hitPoint(), @@ -1789,15 +1861,15 @@ bool Foam::conformalVoronoiMesh::nearSurfacePoint norm ); - const vector nA = norm[0]; + const vector& nA = norm[0]; pointIndexHit oppositeHit; label oppositeSurfaceHit = -1; geometryToConformTo_.findSurfaceNearestIntersection ( - closePoint.hitPoint() + SMALL*nA, - closePoint.hitPoint() + mag(pt - closePoint.hitPoint())*nA, + closePoint.hitPoint() + 0.5*pointPairDistance(pt)*nA, + closePoint.hitPoint() + 5*targetCellSize(pt)*nA, oppositeHit, oppositeSurfaceHit ); @@ -1808,21 +1880,10 @@ bool Foam::conformalVoronoiMesh::nearSurfacePoint pHit.first() = oppositeHit; pHit.second() = oppositeSurfaceHit; -// appendToSurfacePtTree(pHit.first().hitPoint()); -// surfaceToTreeShape.append -// ( -// existingSurfacePtLocations_.size() - 1 -// ); - return !closeToSurfacePt; } } } - else - { -// appendToSurfacePtTree(pt); -// surfaceToTreeShape.append(existingSurfacePtLocations_.size() - 1); - } return closeToSurfacePt; } diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H index 64534cb839..e9e03a6330 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H @@ -230,13 +230,13 @@ inline void Foam::conformalVoronoiMesh::createPointPair { vector ppDistn = ppDist*n; - const Foam::point internalPt = surfPt - ppDistn; - const Foam::point externalPt = surfPt + ppDistn; +// const Foam::point internalPt = surfPt - ppDistn; +// const Foam::point externalPt = surfPt + ppDistn; - bool internalInside = geometryToConformTo_.inside(internalPt); - bool externalOutside = geometryToConformTo_.outside(externalPt); +// bool internalInside = geometryToConformTo_.inside(internalPt); +// bool externalOutside = geometryToConformTo_.outside(externalPt); - if (internalInside && externalOutside) +// if (internalInside && externalOutside) { pts.append ( @@ -259,15 +259,24 @@ inline void Foam::conformalVoronoiMesh::createPointPair Pstream::myProcNo() ) ); + +// if (ptPair) +// { +// ptPairs_.addPointPair +// ( +// pts[pts.size() - 2].index(), +// pts[pts.size() - 1].index() // external 0 -> slave +// ); +// } } - else - { - Info<< "Warning: point pair not inside/outside" << nl - << " surfPt = " << surfPt << nl - << " internal = " << internalPt << " " << internalInside << nl - << " external = " << externalPt << " " << externalOutside - << endl; - } +// else +// { +// Info<< "Warning: point pair not inside/outside" << nl +// << " surfPt = " << surfPt << nl +// << " internal = " << internalPt << " " << internalInside << nl +// << " external = " << externalPt << " " << externalOutside +// << endl; +// } } diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H index 7a95b036a9..82bbe43b95 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H @@ -183,6 +183,7 @@ public: inline bool hasSeedPoint() const; inline bool hasInternalPoint() const; + inline bool hasBoundaryPoint() const; inline bool hasConstrainedPoint() const; 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 8e035dbc7b..2f5d076138 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCellI.H +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCellI.H @@ -228,6 +228,19 @@ inline bool CGAL::indexedCell::hasInternalPoint() const } +template +inline bool CGAL::indexedCell::hasBoundaryPoint() const +{ + return + ( + this->vertex(0)->boundaryPoint() + || this->vertex(1)->boundaryPoint() + || this->vertex(2)->boundaryPoint() + || this->vertex(3)->boundaryPoint() + ); +} + + template inline bool CGAL::indexedCell::hasConstrainedPoint() const { diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.C index 0360cf0d83..73bd6038d9 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.C @@ -683,14 +683,18 @@ Foam::Field Foam::conformationSurfaces::wellInOutSide vector hitDir = info[0].rawPoint() - samplePts[i]; hitDir /= mag(hitDir) + SMALL; - if + pointIndexHit surfHit; + label hitSurface; + + findSurfaceNearestIntersection ( - findSurfaceAnyIntersection - ( - samplePts[i], - info[0].rawPoint() - 1e-3*mag(hitDir)*(hitDir) - ) - ) + samplePts[i], + info[0].rawPoint(), + surfHit, + hitSurface + ); + + if (surfHit.hit() && hitSurface != surfaces_[s]) { continue; } diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/rayShooting/rayShooting.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/rayShooting/rayShooting.C index ebecad88e1..66dc10737d 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/rayShooting/rayShooting.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/rayShooting/rayShooting.C @@ -37,6 +37,90 @@ namespace Foam defineTypeNameAndDebug(rayShooting, 0); addToRunTimeSelectionTable(initialPointsMethod, rayShooting, dictionary); + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +void rayShooting::splitLine +( + const line& l, + const scalar& pert, + DynamicList& initialPoints +) const +{ + Foam::point midPoint(l.centre()); + const scalar localCellSize(cellShapeControls().cellSize(midPoint)); + const scalar lineLength(l.mag()); + + const scalar minDistFromSurfaceSqr + ( + minimumSurfaceDistanceCoeffSqr_ + *sqr(localCellSize) + ); + + if + ( + magSqr(midPoint - l.start()) > minDistFromSurfaceSqr + && magSqr(midPoint - l.end()) > minDistFromSurfaceSqr + ) + { + // Add extra points if line length is much bigger than local cell size +// if (lineLength > 4.0*localCellSize) +// { +// splitLine +// ( +// line(l.start(), midPoint), +// pert, +// initialPoints +// ); +// +// splitLine +// ( +// line(midPoint, l.end()), +// pert, +// initialPoints +// ); +// } + + if (randomiseInitialGrid_) + { + Foam::point newPt + ( + midPoint.x() + pert*(rndGen().scalar01() - 0.5), + midPoint.y() + pert*(rndGen().scalar01() - 0.5), + midPoint.z() + pert*(rndGen().scalar01() - 0.5) + ); + + if + ( + !geometryToConformTo().findSurfaceAnyIntersection + ( + midPoint, + newPt + ) + ) + { + midPoint = newPt; + } + else + { + WarningIn + ( + "rayShooting::splitLine" + "(" + " const line&," + " const scalar&," + " DynamicList&" + ")" + ) << "Point perturbation crosses a surface. Not inserting." + << endl; + } + } + + initialPoints.append(toPoint(midPoint)); + } +} + + // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // rayShooting::rayShooting @@ -75,7 +159,7 @@ List rayShooting::initialPoints() const const searchableSurfaces& surfaces = geometryToConformTo().geometry(); const labelList& surfacesToConformTo = geometryToConformTo().surfaces(); - const scalar maxRayLength = surfaces.bounds().mag(); + const scalar maxRayLength(surfaces.bounds().mag()); // Initialise points list label initialPointsSize = 0; @@ -94,8 +178,7 @@ List rayShooting::initialPoints() const const pointField& faceCentres = faceCentresTmp(); Info<< " Shoot rays from " << s.name() << nl - << " nRays = " << faceCentres.size() << endl; - + << " nRays = " << faceCentres.size() << endl; forAll(faceCentres, fcI) { @@ -110,9 +193,11 @@ List rayShooting::initialPoints() const continue; } - const scalar pert = + const scalar pert + ( randomPerturbationCoeff_ - *cellShapeControls().cellSize(fC); + *cellShapeControls().cellSize(fC) + ); pointIndexHit surfHitStart; label hitSurfaceStart; @@ -181,27 +266,12 @@ List rayShooting::initialPoints() const } } - Foam::point midPoint(l.centre()); - - const scalar minDistFromSurfaceSqr = - minimumSurfaceDistanceCoeffSqr_ - *sqr(cellShapeControls().cellSize(midPoint)); - - if (randomiseInitialGrid_) - { - midPoint.x() += pert*(rndGen().scalar01() - 0.5); - midPoint.y() += pert*(rndGen().scalar01() - 0.5); - midPoint.z() += pert*(rndGen().scalar01() - 0.5); - } - - if + splitLine ( - magSqr(midPoint - l.start()) > minDistFromSurfaceSqr - && magSqr(midPoint - l.end()) > minDistFromSurfaceSqr - ) - { - initialPoints.append(toPoint(midPoint)); - } + l, + pert, + initialPoints + ); } } } diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/rayShooting/rayShooting.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/rayShooting/rayShooting.H index e66df62b7c..e5bd1936ff 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/rayShooting/rayShooting.H +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/rayShooting/rayShooting.H @@ -61,6 +61,16 @@ private: scalar randomPerturbationCoeff_; + // Private Member Functions + + void splitLine + ( + const line& l, + const scalar& pert, + DynamicList& initialPoints + ) const; + + public: //- Runtime type information