diff --git a/src/overset/cellCellStencil/trackingInverseDistance/trackingInverseDistanceCellCellStencil.C b/src/overset/cellCellStencil/trackingInverseDistance/trackingInverseDistanceCellCellStencil.C index 988404591f..96fb1a445c 100644 --- a/src/overset/cellCellStencil/trackingInverseDistance/trackingInverseDistanceCellCellStencil.C +++ b/src/overset/cellCellStencil/trackingInverseDistance/trackingInverseDistanceCellCellStencil.C @@ -48,13 +48,13 @@ namespace cellCellStencils // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // -void Foam::cellCellStencils::trackingInverseDistance::markBoundaries +bool Foam::cellCellStencils::trackingInverseDistance::markBoundaries ( const fvMesh& mesh, const vector& smallVec, const boundBox& bb, - const labelVector& nDivs, + labelVector& nDivs, PackedList<2>& patchTypes, const labelList& cellMap, @@ -63,6 +63,8 @@ void Foam::cellCellStencils::trackingInverseDistance::markBoundaries { // Mark all voxels that overlap the bounding box of any patch + static bool hasWarned = false; + const fvBoundaryMesh& pbm = mesh.boundary(); patchTypes = patchCellType::OTHER; @@ -113,6 +115,7 @@ void Foam::cellCellStencils::trackingInverseDistance::markBoundaries { //Info<< "Marking cells on overset patch " << fvp.name() << endl; const polyPatch& pp = fvp.patch(); + const vectorField::subField faceCentres(pp.faceCentres()); forAll(pp, i) { // Mark in overall patch types @@ -122,8 +125,54 @@ void Foam::cellCellStencils::trackingInverseDistance::markBoundaries boundBox faceBb(pp.points(), pp[i], false); faceBb.min() -= smallVec; faceBb.max() += smallVec; - if (bb.overlaps(faceBb)) + if (!bb.contains(faceCentres[i])) { + if (!hasWarned) + { + hasWarned = true; + WarningInFunction << "Patch " << pp.name() + << " : face at " << faceCentres[i] + << " is not inside search bounding box of" + << " voxel mesh " << bb << endl + << " Is your 'searchBox' specification correct?" + << endl; + } + } + else + { + // Test for having voxels already marked as patch + // -> voxel mesh is too coarse + if + ( + voxelMeshSearch::overlaps + ( + bb, + nDivs, + faceBb, + patchTypes, + static_cast(patchCellType::PATCH) + ) + ) + { + // Determine number of voxels from number of cells + // in mesh + const labelVector& dim = mesh.geometricD(); + forAll(dim, i) + { + if (dim[i] != -1) + { + nDivs[i] *= 2; + } + } + + Pout<< "Voxel mesh too coarse. Bounding box " + << faceBb + << " contains both non-overset and overset patches" + << ". Refining voxel mesh to " << nDivs << endl; + + return false; + } + voxelMeshSearch::fill ( patchTypes, @@ -136,6 +185,8 @@ void Foam::cellCellStencils::trackingInverseDistance::markBoundaries } } } + + return true; } @@ -522,64 +573,15 @@ bool Foam::cellCellStencils::trackingInverseDistance::update() DebugInfo<< FUNCTION_NAME << " : Moved zone sub-meshes" << endl; - // Calculate fast search structure for each zone - PtrList meshSearches(nZones); - - List searchBoxDivisions; - if (dict_.readIfPresent("searchBoxDivisions", searchBoxDivisions)) + List searchBoxDivisions(dict_["searchBoxDivisions"]); + if (searchBoxDivisions.size() != nZones) { - forAll(meshParts_, zonei) - { - meshSearches.set - ( - zonei, - new voxelMeshSearch - ( - meshParts_[zonei].subMesh(), - searchBoxDivisions[zonei], - true - ) - ); - } + FatalIOErrorInFunction(dict_) << "Number of searchBoxDivisions " + << searchBoxDivisions.size() + << " should equal the number of zones " << nZones + << exit(FatalIOError); } - else - { - forAll(meshParts_, zonei) - { - meshSearches.set - ( - zonei, - new voxelMeshSearch - ( - meshParts_[zonei].subMesh(), - true - ) - ); - } - } - DebugInfo<< FUNCTION_NAME << " : Constructed cell voxel meshes" << endl; - - PtrList voxelMeshes; - if (debug) - { - voxelMeshes.setSize(meshSearches.size()); - forAll(meshSearches, zonei) - { - IOobject meshIO - ( - word("voxelMesh" + Foam::name(zonei)), - mesh_.time().timeName(), - mesh_, - IOobject::NO_READ - ); - - Pout<< "Constructing voxel mesh " << meshIO.objectPath() << endl; - voxelMeshes.set(zonei, meshSearches[zonei].makeMesh(meshIO)); - } - } - - const boundBox& allBb(mesh_.bounds()); @@ -652,35 +654,78 @@ bool Foam::cellCellStencils::trackingInverseDistance::update() patchBb = meshBb; } } - if (patchDivisions.empty()) - { - patchDivisions.setSize(nZones); - forAll(patchDivisions, zonei) - { - patchDivisions[zonei] = meshSearches[zonei].nDivs(); - } - } forAll(patchParts, zonei) { - const labelVector& g = patchDivisions[zonei]; - patchParts.set(zonei, new PackedList<2>(cmptProduct(g))); + while (true) + { + patchParts.set + ( + zonei, + new PackedList<2>(cmptProduct(patchDivisions[zonei])) + ); - markBoundaries - ( - meshParts_[zonei].subMesh(), - smallVec_, + bool ok = markBoundaries + ( + meshParts_[zonei].subMesh(), + smallVec_, - patchBb[zonei][Pstream::myProcNo()], - patchDivisions[zonei], - patchParts[zonei], + patchBb[zonei][Pstream::myProcNo()], + patchDivisions[zonei], + patchParts[zonei], - meshParts_[zonei].cellMap(), - allPatchTypes - ); + meshParts_[zonei].cellMap(), + allPatchTypes + ); + + //if (returnReduce(ok, andOp())) + if (ok) + { + break; + } + } } DebugInfo<< FUNCTION_NAME << " : Calculated boundary voxel meshes" << endl; + + PtrList meshSearches(meshParts_.size()); + forAll(meshParts_, zonei) + { + meshSearches.set + ( + zonei, + new voxelMeshSearch + ( + meshParts_[zonei].subMesh(), + patchBb[zonei][Pstream::myProcNo()], + patchDivisions[zonei], + true + ) + ); + } + + DebugInfo<< FUNCTION_NAME << " : Constructed cell voxel meshes" << endl; + + PtrList voxelMeshes; + if (debug) + { + voxelMeshes.setSize(meshSearches.size()); + forAll(meshSearches, zonei) + { + IOobject meshIO + ( + word("voxelMesh" + Foam::name(zonei)), + mesh_.time().timeName(), + mesh_, + IOobject::NO_READ + ); + + Pout<< "Constructing voxel mesh " << meshIO.objectPath() << endl; + voxelMeshes.set(zonei, meshSearches[zonei].makeMesh(meshIO)); + } + } + + if (debug) { forAll(patchParts, zonei) diff --git a/src/overset/cellCellStencil/trackingInverseDistance/trackingInverseDistanceCellCellStencil.H b/src/overset/cellCellStencil/trackingInverseDistance/trackingInverseDistanceCellCellStencil.H index d5951c2e22..d5caa9acb6 100644 --- a/src/overset/cellCellStencil/trackingInverseDistance/trackingInverseDistanceCellCellStencil.H +++ b/src/overset/cellCellStencil/trackingInverseDistance/trackingInverseDistanceCellCellStencil.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2017 OpenCFD Ltd. + \\ / A nd | Copyright (C) 2017-2019 OpenCFD Ltd. \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -74,13 +74,13 @@ protected: // Protected Member Functions //- Mark voxels of patchTypes with type of patch face - static void markBoundaries + static bool markBoundaries ( const fvMesh& mesh, const vector& smallVec, const boundBox& bb, - const labelVector& nDivs, + labelVector& nDivs, PackedList<2>& patchTypes, const labelList& cellMap, diff --git a/src/overset/cellCellStencil/trackingInverseDistance/voxelMeshSearch.C b/src/overset/cellCellStencil/trackingInverseDistance/voxelMeshSearch.C index dd33103499..d1e23e278e 100644 --- a/src/overset/cellCellStencil/trackingInverseDistance/voxelMeshSearch.C +++ b/src/overset/cellCellStencil/trackingInverseDistance/voxelMeshSearch.C @@ -281,6 +281,14 @@ Foam::voxelMeshSearch::voxelMeshSearch } } + // Redo the local bounding box + localBb_ = boundBox(mesh_.points(), false); + + const point eps(1e-10, 1e-10, 1e-10); + + localBb_.min() = localBb_.min()-eps; + localBb_.max() = localBb_.max()+eps; + if (debug) { Pout<< "voxelMeshSearch : mesh:" << mesh_.name() @@ -297,11 +305,13 @@ Foam::voxelMeshSearch::voxelMeshSearch Foam::voxelMeshSearch::voxelMeshSearch ( const polyMesh& mesh, + const boundBox& localBb, const labelVector& nDivs, const bool doUpdate ) : mesh_(mesh), + localBb_(localBb), nDivs_(nDivs) { if (doUpdate) @@ -315,14 +325,6 @@ Foam::voxelMeshSearch::voxelMeshSearch bool Foam::voxelMeshSearch::update() { - // Redo the local bounding box - localBb_ = boundBox(mesh_.points(), false); - - const point eps(1e-10, 1e-10, 1e-10); - - localBb_.min() = localBb_.min()-eps; - localBb_.max() = localBb_.max()+eps; - // Initialise seed cell array seedCell_.setSize(cmptProduct(nDivs_)); @@ -373,6 +375,7 @@ Foam::label Foam::voxelMeshSearch::findCell(const point& p) const return -1; } + // Locate the voxel index for this point. Do not clip. label voxeli = index(localBb_, nDivs_, p, false); @@ -397,10 +400,14 @@ Foam::label Foam::voxelMeshSearch::findCell(const point& p) const // found does not have to be the absolute 'correct' one as // long as at least one of the processors finds a cell. - label nextCellOld = -1; - + track_.clear(); while (true) { + if (track_.size() < 5) + { + track_.append(celli); + } + // I am in celli now. How many faces do I have ? label facei = findIntersectedFace(celli, p); @@ -409,6 +416,7 @@ Foam::label Foam::voxelMeshSearch::findCell(const point& p) const return celli; } + const label startOfTrack(max(0, track_.size()-5)); label nextCell; if (mesh_.isInternalFace(facei)) @@ -416,21 +424,24 @@ Foam::label Foam::voxelMeshSearch::findCell(const point& p) const label own = mesh_.faceOwner()[facei]; label nei = mesh_.faceNeighbour()[facei]; nextCell = (own == celli ? nei : own); + + if (findIndex(track_, nextCell, startOfTrack) != -1) + { + return celli; + } } else { nextCell = searchProcPatch(facei, p); - if (nextCell == nextCellOld) - { - return -1; // point is really out - } - if (nextCell == -1 || nextCell == celli) { return nextCell; } - nextCellOld = nextCell; + else if (findIndex(track_, nextCell, startOfTrack) != -1) + { + return -1; // point is really out + } } celli = nextCell; diff --git a/src/overset/cellCellStencil/trackingInverseDistance/voxelMeshSearch.H b/src/overset/cellCellStencil/trackingInverseDistance/voxelMeshSearch.H index 5c929fe44b..d7ce1e543f 100644 --- a/src/overset/cellCellStencil/trackingInverseDistance/voxelMeshSearch.H +++ b/src/overset/cellCellStencil/trackingInverseDistance/voxelMeshSearch.H @@ -68,6 +68,9 @@ class voxelMeshSearch //- Voxel to seed cell labelList seedCell_; + //- Cells visited + mutable DynamicList