diff --git a/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C b/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C index c910b4c1b57a72b61bc583952bf793a07c6f86ef..0f81dff2f5e477815eda5188a35852163057cc9e 100644 --- a/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C +++ b/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C @@ -496,23 +496,63 @@ Foam::labelList Foam::refinementSurfaces::getClosedNamedSurfaces() const } -// Count number of triangles per surface region -Foam::labelList Foam::refinementSurfaces::countRegions(const triSurface& s) -{ - const geometricSurfacePatchList& regions = s.patches(); - - labelList nTris(regions.size(), 0); - - forAll(s, triI) - { - nTris[s[triI].region()]++; - } - return nTris; -} +// // Count number of triangles per surface region +// Foam::labelList Foam::refinementSurfaces::countRegions(const triSurface& s) +// { +// const geometricSurfacePatchList& regions = s.patches(); +// +// labelList nTris(regions.size(), 0); +// +// forAll(s, triI) +// { +// nTris[s[triI].region()]++; +// } +// return nTris; +// } + + +// // Pre-calculate the refinement level for every element +// void Foam::refinementSurfaces::wantedRefinementLevel +// ( +// const shellSurfaces& shells, +// const label surfI, +// const List<pointIndexHit>& info, // Indices +// const pointField& ctrs, // Representative coordinate +// labelList& minLevelField +// ) const +// { +// const searchableSurface& geom = allGeometry_[surfaces_[surfI]]; +// +// // Get per element the region +// labelList region; +// geom.getRegion(info, region); +// +// // Initialise fields to region wise minLevel +// minLevelField.setSize(ctrs.size()); +// minLevelField = -1; +// +// forAll(minLevelField, i) +// { +// if (info[i].hit()) +// { +// minLevelField[i] = minLevel(surfI, region[i]); +// } +// } +// +// // Find out if triangle inside shell with higher level +// // What level does shell want to refine fc to? +// labelList shellLevel; +// shells.findHigherLevel(ctrs, minLevelField, shellLevel); +// +// forAll(minLevelField, i) +// { +// minLevelField[i] = max(minLevelField[i], shellLevel[i]); +// } +// } // Precalculate the refinement level for every element of the searchable -// surface. This is currently hardcoded for triSurfaceMesh only. +// surface. void Foam::refinementSurfaces::setMinLevelFields ( const shellSurfaces& shells @@ -522,58 +562,55 @@ void Foam::refinementSurfaces::setMinLevelFields { const searchableSurface& geom = allGeometry_[surfaces_[surfI]]; - if (isA<triSurfaceMesh>(geom)) + // Precalculation only makes sense if there are different regions + // (so different refinement levels possible) and there are some + // elements. Possibly should have 'enough' elements to have fine + // enough resolution but for now just make sure we don't catch e.g. + // searchableBox (size=6) + if (geom.regions().size() > 1 && geom.globalSize() > 10) { - const triSurfaceMesh& triMesh = refCast<const triSurfaceMesh>(geom); + // Representative local coordinates + const pointField ctrs = geom.coordinates(); - autoPtr<triSurfaceLabelField> minLevelFieldPtr - ( - new triSurfaceLabelField + labelList minLevelField(ctrs.size(), -1); + { + // Get the element index in a roundabout way. Problem is e.g. + // distributed surface where local indices differ from global + // ones (needed for getRegion call) + List<pointIndexHit> info; + geom.findNearest ( - IOobject - ( - "minLevel", - triMesh.objectRegistry::time().timeName(), // instance - "triSurface", // local - triMesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE - ), - triMesh, - dimless - ) - ); - triSurfaceLabelField& minLevelField = minLevelFieldPtr(); + ctrs, + scalarField(ctrs.size(), sqr(GREAT)), + info + ); - const triSurface& s = static_cast<const triSurface&>(triMesh); + // Get per element the region + labelList region; + geom.getRegion(info, region); - // Initialise fields to region wise minLevel - forAll(s, triI) - { - minLevelField[triI] = minLevel(surfI, s[triI].region()); + // From the region get the surface-wise refinement level + forAll(minLevelField, i) + { + if (info[i].hit()) + { + minLevelField[i] = minLevel(surfI, region[i]); + } + } } // Find out if triangle inside shell with higher level - pointField fc(s.size()); - forAll(s, triI) - { - fc[triI] = s[triI].centre(s.points()); - } // What level does shell want to refine fc to? labelList shellLevel; - shells.findHigherLevel(fc, minLevelField, shellLevel); + shells.findHigherLevel(ctrs, minLevelField, shellLevel); - forAll(minLevelField, triI) + forAll(minLevelField, i) { - minLevelField[triI] = max - ( - minLevelField[triI], - shellLevel[triI] - ); + minLevelField[i] = max(minLevelField[i], shellLevel[i]); } - // Store field on triMesh - minLevelFieldPtr.ptr()->store(); + // Store minLevelField on surface + const_cast<searchableSurface&>(geom).setField(minLevelField); } } } @@ -601,44 +638,89 @@ void Foam::refinementSurfaces::findHigherIntersection return; } - // Work arrays - labelList hitMap(identity(start.size())); - pointField p0(start); - pointField p1(end); - List<pointIndexHit> hitInfo(start.size()); - - forAll(surfaces_, surfI) + if (surfaces_.size() == 1) { + // Optimisation: single segmented surface. No need to duplicate + // point storage. + + label surfI = 0; + const searchableSurface& geom = allGeometry_[surfaces_[surfI]]; - geom.findLineAny(p0, p1, hitInfo); + // Do intersection test + List<pointIndexHit> intersectionInfo(start.size()); + geom.findLineAny(start, end, intersectionInfo); + // See if a cached level field available labelList minLevelField; - if (isA<triSurfaceMesh>(geom)) + geom.getField(intersectionInfo, minLevelField); + bool haveLevelField = + ( + returnReduce(minLevelField.size(), sumOp<label>()) + > 0 + ); + + if (!haveLevelField && geom.regions().size() == 1) { - const triSurfaceMesh& triMesh = refCast<const triSurfaceMesh>(geom); - triMesh.getField + minLevelField = labelList ( - "minLevel", - hitInfo, - minLevelField + intersectionInfo.size(), + minLevel(surfI, 0) ); + haveLevelField = true; } - // Copy all hits into arguments, continue with misses - label newI = 0; - forAll(hitInfo, hitI) + if (haveLevelField) + { + forAll(intersectionInfo, i) + { + if + ( + intersectionInfo[i].hit() + && minLevelField[i] > currentLevel[i] + ) + { + surfaces[i] = surfI; // index of surface + surfaceLevel[i] = minLevelField[i]; + } + } + return; + } + } + + + + // Work arrays + pointField p0(start); + pointField p1(end); + labelList intersectionToPoint(identity(start.size())); + List<pointIndexHit> intersectionInfo(start.size()); + + forAll(surfaces_, surfI) + { + const searchableSurface& geom = allGeometry_[surfaces_[surfI]]; + + // Do intersection test + geom.findLineAny(p0, p1, intersectionInfo); + + // See if a cached level field available + labelList minLevelField; + geom.getField(intersectionInfo, minLevelField); + + // Copy all hits into arguments, In-place compact misses. + label missI = 0; + forAll(intersectionInfo, i) { // Get the minLevel for the point label minLocalLevel = -1; - if (hitInfo[hitI].hit()) + if (intersectionInfo[i].hit()) { // Check if minLevelField for this surface. if (minLevelField.size()) { - minLocalLevel = minLevelField[hitI]; + minLocalLevel = minLevelField[i]; } else { @@ -648,36 +730,35 @@ void Foam::refinementSurfaces::findHigherIntersection } } - label pointI = hitMap[hitI]; + + label pointI = intersectionToPoint[i]; if (minLocalLevel > currentLevel[pointI]) { + // Mark point for refinement surfaces[pointI] = surfI; surfaceLevel[pointI] = minLocalLevel; } else { - if (hitI != newI) - { - hitMap[newI] = hitMap[hitI]; - p0[newI] = p0[hitI]; - p1[newI] = p1[hitI]; - } - newI++; + p0[missI] = p0[pointI]; + p1[missI] = p1[pointI]; + intersectionToPoint[missI] = pointI; + missI++; } } // All done? Note that this decision should be synchronised - if (returnReduce(newI, sumOp<label>()) == 0) + if (returnReduce(missI, sumOp<label>()) == 0) { break; } - // Trim and continue - hitMap.setSize(newI); - p0.setSize(newI); - p1.setSize(newI); - hitInfo.setSize(newI); + // Trim misses + p0.setSize(missI); + p1.setSize(missI); + intersectionToPoint.setSize(missI); + intersectionInfo.setSize(missI); } } diff --git a/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.H b/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.H index 86767f57312e470bcf50f183d8a8a73d8e08642d..b1188f1dff79c78754c22f707a00758dda48b1e8 100644 --- a/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.H +++ b/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.H @@ -218,8 +218,8 @@ public: const shellSurfaces& shells ); - //- Helper: count number of triangles per region - static labelList countRegions(const triSurface&); + ////- Helper: count number of triangles per region + //static labelList countRegions(const triSurface&); // Searching