diff --git a/meshLibrary/utilities/checkMeshDict/checkMeshDict.C b/meshLibrary/utilities/checkMeshDict/checkMeshDict.C index 471a262bb31d5e2729bb113f833bf0b74846d30a..8116ff2256d667e96362031f6b85a7591e19094a 100644 --- a/meshLibrary/utilities/checkMeshDict/checkMeshDict.C +++ b/meshLibrary/utilities/checkMeshDict/checkMeshDict.C @@ -226,6 +226,85 @@ void checkMeshDict::checkObjectRefinements() const } } +void checkMeshDict::checkSurfaceRefinements() const +{ + if( meshDict_.found("surfaceMeshRefinement") ) + { + const dictionary& surfaces = meshDict_.subDict("surfaceMeshRefinement"); + + const wordList surfaceSources = surfaces.toc(); + + forAll(surfaceSources, surfI) + { + if( surfaces.isDict(surfaceSources[surfI]) ) + { + const dictionary& dict = surfaces.subDict(surfaceSources[surfI]); + + if( dict.found("surfaceFile") ) + { + const fileName fName(dict.lookup("surfaceFile")); + + if( !isFile(fName) ) + FatalErrorIn + ( + "void checkMeshDict::checkSurfaceRefinements() const" + ) << "Surface file " << fName + << " does not exist or is not readable!!" + << exit(FatalError); + + if( dict.found("cellSize") ) + { + const scalar cs = readScalar(dict.lookup("cellSize")); + + if( cs < VSMALL ) + FatalErrorIn + ( + "void checkMeshDict::" + "checkSurfaceRefinements() const" + ) << "Cell size for surface " << fName + << " is extremely small or negative!!" + << exit(FatalError); + } + else if( dict.found("additionalRefinementLevels") ) + { + const label nLev = + readLabel(dict.lookup("additionalRefinementLevels")); + + if( nLev < 0 ) + { + FatalErrorIn + ( + "void checkMeshDict::" + "checkSurfaceRefinements() const" + ) << "Number refinement levels for surface " << fName + << " is negative!!" + << exit(FatalError); + } + } + else + { + FatalErrorIn + ( + "void checkMeshDict::checkSurfaceRefinements() const" + ) << "Missing cellSize or additionalRefinementLevels" + << " for surface " << fName + << exit(FatalError); + } + } + } + else + { + FatalErrorIn + ( + "void checkMeshDict::checkSurfaceRefinements() const" + ) << "Dictionary " << surfaceSources[surfI] + << " does not exist!!" + << exit(FatalError); + } + } + } +} + void checkMeshDict::checkBoundaryLayers() const { if( meshDict_.found("boundaryLayers") ) @@ -356,6 +435,8 @@ void checkMeshDict::checkEntries() const checkSubsetCellSize(); + checkSurfaceRefinements(); + checkKeepCellsIntersectingPatches(); checkRemoveCellsIntersectingPatches(); diff --git a/meshLibrary/utilities/checkMeshDict/checkMeshDict.H b/meshLibrary/utilities/checkMeshDict/checkMeshDict.H index c97d2addb6bc4c5bfd06f9ec27a5c6f7c662a50a..04c7c1332a121abcdc30944127e363a2fb097db0 100644 --- a/meshLibrary/utilities/checkMeshDict/checkMeshDict.H +++ b/meshLibrary/utilities/checkMeshDict/checkMeshDict.H @@ -72,6 +72,9 @@ class checkMeshDict //- check objectRefinements entry void checkObjectRefinements() const; + //- check surfaceRefinements entry + void checkSurfaceRefinements() const; + //- check entry for boundary layers void checkBoundaryLayers() const; diff --git a/meshLibrary/utilities/octrees/meshOctree/meshOctreeCreator/meshOctreeCreator.H b/meshLibrary/utilities/octrees/meshOctree/meshOctreeCreator/meshOctreeCreator.H index f2d2ee87bdba8cf341a78ebae8258254f72d8542..231ef332c62974b310891abbd4edecd79385fbfe 100644 --- a/meshLibrary/utilities/octrees/meshOctree/meshOctreeCreator/meshOctreeCreator.H +++ b/meshLibrary/utilities/octrees/meshOctree/meshOctreeCreator/meshOctreeCreator.H @@ -80,6 +80,9 @@ private: //- refine boxes contained inside the objects for refinement void refineBoxesContainedInObjects(); + //- refine boxes intersected by surface meshes used as refinement sources + void refineBoxesIntersectingSurfaces(); + //- refine boxes near DATA boxes to get a nice smooth surface void refineBoxesNearDataBoxes(const direction nLayers = 1); diff --git a/meshLibrary/utilities/octrees/meshOctree/meshOctreeCreator/meshOctreeCreatorAdjustOctreeToSurface.C b/meshLibrary/utilities/octrees/meshOctree/meshOctreeCreator/meshOctreeCreatorAdjustOctreeToSurface.C index a066a72ac89944c000aa388e075bfbf2eebaa044..c05d02a6f7d67ba553dbf25b3c0fdf9f7d53ece5 100644 --- a/meshLibrary/utilities/octrees/meshOctree/meshOctreeCreator/meshOctreeCreatorAdjustOctreeToSurface.C +++ b/meshLibrary/utilities/octrees/meshOctree/meshOctreeCreator/meshOctreeCreatorAdjustOctreeToSurface.C @@ -32,6 +32,7 @@ Description #include "objectRefinementList.H" #include "VRWGraph.H" #include "meshOctreeModifier.H" +#include "helperFunctions.H" #include "HashSet.H" # ifdef USE_OMP @@ -147,7 +148,7 @@ void meshOctreeCreator::refineBoxesContainedInObjects() Info << "Refining boxes inside objects" << endl; objectRefinementList refObjects; - // Read polyPatchList + // Read objects if( meshDictPtr_->isDict("objectRefinements") ) { const dictionary& dict = meshDictPtr_->subDict("objectRefinements"); @@ -309,6 +310,182 @@ void meshOctreeCreator::refineBoxesContainedInObjects() createInsideOutsideInformation(); } +void meshOctreeCreator::refineBoxesIntersectingSurfaces() +{ + if( !meshDictPtr_ || !meshDictPtr_->found("surfaceMeshRefinement") ) + { + return; + } + + Info << "Refining boxes intersecting surface meshes" << endl; + + label nMarked; + + //- read surface meshes and calculate the refinement level for each + //- surface mesh + const dictionary& surfDict = meshDictPtr_->subDict("surfaceMeshRefinement"); + const wordList surfaces = surfDict.toc(); + PtrList<triSurf> surfaceMeshesPtr(surfaces.size()); + List<direction> refLevels(surfaces.size(), globalRefLevel_); + + //- load surface meshes into memory + forAll(surfaceMeshesPtr, surfI) + { + const dictionary& dict = surfDict.subDict(surfaces[surfI]); + + const fileName fName(dict.lookup("surfaceFile")); + + surfaceMeshesPtr.set + ( + surfI, + new triSurf(fName) + ); + + direction addLevel(0); + if( dict.found("cellSize") ) + { + scalar s(readScalar(meshDictPtr_->lookup("maxCellSize"))); + + const scalar cs = readScalar(dict.lookup("cellSize")); + + do + { + nMarked = 0; + if( cs <= s * (1.+SMALL) ) + { + ++nMarked; + ++addLevel; + } + + s /= 2.0; + + } while( nMarked != 0 ); + } + else if( dict.found("additionalRefinementLevels") ) + { + addLevel = + readLabel(dict.lookup("additionalRefinementLevels")); + } + + //- set the refinement level for the current surface + refLevels[surfI] += addLevel; + } + + if( octree_.neiProcs().size() ) + forAll(refLevels, oI) + { + label l = refLevels[oI]; + reduce(l, maxOp<label>()); + refLevels[oI] = l; + } + + //- start refining boxes intersecting triangles in each refinement surface + const boundBox& rootBox = octree_.rootBox(); + const vector tol = SMALL * rootBox.span(); + meshOctreeModifier octreeModifier(octree_); + const LongList<meshOctreeCube*>& leaves = octreeModifier.leavesAccess(); + DynList<label> leavesInBox; + + do + { + # ifdef OCTREETiming + const scalar startIter = omp_get_wtime(); + # endif + + nMarked = 0; + + List<direction> refineCubes(leaves.size(), direction(0)); + + //- select boxes which need to be refined + forAll(surfaceMeshesPtr, surfI) + { + const triSurf& surf = surfaceMeshesPtr[surfI]; + const pointField& points = surf.points(); + + # ifdef USE_OMP + # pragma omp parallel for \ + reduction( + : nMarked) schedule(dynamic, 10) private(leavesInBox) + # endif + forAll(surf, triI) + { + //- find the bounding box of the current triangle + const labelledTri& tri = surf[triI]; + boundBox triBB(points[tri[0]], points[tri[0]]); + for(label pI=1;pI<3;++pI) + { + triBB.min() = Foam::min(triBB.min(), points[tri[pI]]); + triBB.max() = Foam::max(triBB.max(), points[tri[pI]]); + } + + triBB.min() -= tol; + triBB.max() += tol; + + //- find octree leaves inside the bounding box + leavesInBox.clear(); + octree_.findLeavesContainedInBox(triBB, leavesInBox); + + //- check which of the leaves are intersected by the triangle + forAll(leavesInBox, i) + { + const label leafI = leavesInBox[i]; + + if( refineCubes[leafI] ) + continue; + + const meshOctreeCube& oc = *leaves[leafI]; + + if( + (oc.level() < refLevels[surfI]) && + oc.intersectsTriangleExact(surf, rootBox, triI) + ) + { + # ifdef DEBUGSearch + Info << "Marking leaf " << leafI + << " with coordinates " << oc + << " for refinement" << endl; + # endif + + ++nMarked; + refineCubes[leafI] = 1; + } + } + } + } + + //- refine boxes + octreeModifier.refineSelectedBoxes(refineCubes, hexRefinement_); + + # ifdef OCTREETiming + const scalar refTime = omp_get_wtime(); + Info << "Time for refinement " << (refTime-startIter) << endl; + # endif + + if( octree_.neiProcs().size() != 0 ) + { + reduce(nMarked, sumOp<label>()); + if( nMarked ) + { + octreeModifier.distributeLeavesToProcessors(); + + # ifdef OCTREETiming + const scalar distTime = omp_get_wtime(); + Info << "Time for distributing data to processors " + << (distTime-refTime) << endl; + # endif + + loadDistribution(false); + + # ifdef OCTREETiming + Info << "Time for load distribution " + << (omp_get_wtime()-distTime) << endl; + # endif + } + } + } while( nMarked != 0 ); + + Info << "Finished refinement of boxes intersecting surface meshes" << endl; +} + void meshOctreeCreator::refineBoxesNearDataBoxes(const direction nLayers) { # ifdef OCTREETiming diff --git a/meshLibrary/utilities/octrees/meshOctree/meshOctreeCreator/meshOctreeCreatorCreateOctreeBoxes.C b/meshLibrary/utilities/octrees/meshOctree/meshOctreeCreator/meshOctreeCreatorCreateOctreeBoxes.C index 1037b1b35d1bd82048e68cea3388fcdcf964e224..82c146af20065879fa313e9f5da6d1cf61f9b59a 100644 --- a/meshLibrary/utilities/octrees/meshOctree/meshOctreeCreator/meshOctreeCreatorCreateOctreeBoxes.C +++ b/meshLibrary/utilities/octrees/meshOctree/meshOctreeCreator/meshOctreeCreatorCreateOctreeBoxes.C @@ -418,6 +418,9 @@ void meshOctreeCreator::createOctreeBoxes() Info << "Refining boundary" << endl; refineBoundary(); + //- refine parts intersected with surface mesh serving as refinement sources + refineBoxesIntersectingSurfaces(); + //- perform automatic octree refinement if( !Pstream::parRun() ) { diff --git a/meshLibrary/utilities/octrees/meshOctree/meshOctreeCube/meshOctreeCubeRecursiveFunctions.C b/meshLibrary/utilities/octrees/meshOctree/meshOctreeCube/meshOctreeCubeRecursiveFunctions.C index 15ddaa0157b3693a53ba60a1e537c8d298d4257c..a9660869b2f3ec4159daf03af22b210296c1d8b4 100644 --- a/meshLibrary/utilities/octrees/meshOctree/meshOctreeCube/meshOctreeCubeRecursiveFunctions.C +++ b/meshLibrary/utilities/octrees/meshOctree/meshOctreeCube/meshOctreeCubeRecursiveFunctions.C @@ -46,9 +46,8 @@ void meshOctreeCube::leavesInBox DynList<const meshOctreeCube*, 256>& leaves ) const { - point min, max; - this->cubeBox(rootBox, min, max); - const boundBox cubeBox(min, max); + boundBox cubeBox; + this->cubeBox(rootBox, cubeBox.min(), cubeBox.max()); if( cubeBox.overlaps(searchingBox) ) { @@ -74,8 +73,8 @@ void meshOctreeCube::leavesInBox else if( Pstream::parRun() ) { meshOctreeCubeCoordinates cc = refineForPosition(scI); - cc.cubeBox(rootBox, min, max); - const boundBox bb(min, max); + boundBox bb; + cc.cubeBox(rootBox, bb.min(), bb.max()); if( bb.overlaps(searchingBox) ) leaves.append(this); }