diff --git a/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementBaffles.C b/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementBaffles.C index d9ec65959ca7c4bebfa9fda4430b9dba3b39275a..5cdab910ac74b0cc8a6bc406c25d81c3a37dcca9 100644 --- a/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementBaffles.C +++ b/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementBaffles.C @@ -51,9 +51,6 @@ License // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // -// Repatches external face or creates baffle for internal face -// with user specified patches (might be different for both sides). -// Returns label of added face. Foam::label Foam::meshRefinement::createBaffle ( const label faceI, @@ -131,46 +128,6 @@ Foam::label Foam::meshRefinement::createBaffle } -//// Check if we are a boundary face and normal of surface does -//// not align with test vector. In this case there'd probably be -//// a freestanding 'baffle' so we might as well not create it. -//// Note that since it is not a proper baffle we cannot detect it -//// afterwards so this code cannot be merged with the -//// filterDuplicateFaces code. -//bool Foam::meshRefinement::validBaffleTopology -//( -// const label faceI, -// const vector& n1, -// const vector& n2, -// const vector& testDir -//) const -//{ -// -// label patchI = mesh_.boundaryMesh().whichPatch(faceI); -// if (patchI == -1 || mesh_.boundaryMesh()[patchI].coupled()) -// { -// return true; -// } -// else if (mag(n1&n2) > cos(degToRad(30))) -// { -// // Both normals aligned. Check that test vector perpendicularish to -// // surface normal -// scalar magTestDir = mag(testDir); -// if (magTestDir > VSMALL) -// { -// if (mag(n1&(testDir/magTestDir)) < cos(degToRad(45))) -// { -// //Pout<< "** disabling baffling face " -// // << mesh_.faceCentres()[faceI] << endl; -// return false; -// } -// } -// } -// return true; -//} - - -// Determine patches for baffles on all intersected unnamed faces void Foam::meshRefinement::getBafflePatches ( const labelList& globalToMasterPatch, @@ -385,11 +342,11 @@ Foam::Map<Foam::labelPair> Foam::meshRefinement::getZoneBafflePatches } } } + return bafflePatch; } -// Create baffle for every face where ownPatch != -1 Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::createBaffles ( const labelList& ownPatch, @@ -509,6 +466,7 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::createBaffles } } } + // Pick up neighbour side of baffle (added faces) forAll(faceMap, faceI) { @@ -535,7 +493,6 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::createBaffles void Foam::meshRefinement::checkZoneFaces() const { const faceZoneMesh& fZones = mesh_.faceZones(); - const polyBoundaryMesh& pbm = mesh_.boundaryMesh(); forAll(pbm, patchI) @@ -664,27 +621,26 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::createZoneBaffles Info<< "Created " << nZoneFaces << " baffles in = " << mesh_.time().cpuTimeIncrement() << " s\n" << nl << endl; } + return map; } -// Extract those baffles (duplicate) faces that are on the edge of a baffle -// region. These are candidates for merging. -// Done by counting the number of baffles faces per mesh edge. If edge -// has 2 boundary faces and both are baffle faces it is the edge of a baffle -// region. Foam::List<Foam::labelPair> Foam::meshRefinement::freeStandingBaffles ( const List<labelPair>& couples, const scalar planarAngle ) const { + // Done by counting the number of baffles faces per mesh edge. If edge + // has 2 boundary faces and both are baffle faces it is the edge of a baffle + // region. + // All duplicate faces on edge of the patch are to be merged. // So we count for all edges of duplicate faces how many duplicate // faces use them. labelList nBafflesPerEdge(mesh_.nEdges(), 0); - // This algorithm is quite tricky. We don't want to use edgeFaces and // also want it to run in parallel so it is now an algorithm over // all (boundary) faces instead. @@ -698,10 +654,8 @@ Foam::List<Foam::labelPair> Foam::meshRefinement::freeStandingBaffles // So now any edge that is used by baffle faces only will have the // value 2*1000000+2*1. - const label baffleValue = 1000000; - // Count number of boundary faces per edge // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ @@ -729,11 +683,9 @@ Foam::List<Foam::labelPair> Foam::meshRefinement::freeStandingBaffles } } - DynamicList<label> fe0; DynamicList<label> fe1; - // Count number of duplicate boundary faces per edge // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ @@ -927,7 +879,6 @@ Foam::List<Foam::labelPair> Foam::meshRefinement::freeStandingBaffles } -// Merge baffles. Gets pairs of faces. Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::mergeBaffles ( const List<labelPair>& couples @@ -1059,7 +1010,6 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::mergeBaffles } -// Finds region per cell for cells inside closed named surfaces void Foam::meshRefinement::findCellZoneGeometric ( const pointField& neiCc, @@ -1279,7 +1229,6 @@ void Foam::meshRefinement::findCellZoneInsideWalk const labelList& locationSurfaces, // indices of surfaces with inside point const labelList& namedSurfaceIndex, // per face index of named surface const labelList& surfaceToCellZone, // cell zone index per surface - labelList& cellToZone ) const { @@ -1298,6 +1247,7 @@ void Foam::meshRefinement::findCellZoneInsideWalk blockedFace[faceI] = true; } } + // No need to sync since namedSurfaceIndex already is synced // Set region per cell based on walking @@ -1433,14 +1383,11 @@ bool Foam::meshRefinement::calcRegionToZone } } } + return changed; } -// Finds region per cell. Assumes: -// - region containing keepPoint does not go into a cellZone -// - all other regions can be found by crossing faces marked in -// namedSurfaceIndex. void Foam::meshRefinement::findCellZoneTopo ( const point& keepPoint, @@ -1449,6 +1396,11 @@ void Foam::meshRefinement::findCellZoneTopo labelList& cellToZone ) const { + // Assumes: + // - region containing keepPoint does not go into a cellZone + // - all other regions can be found by crossing faces marked in + // namedSurfaceIndex. + // Analyse regions. Reuse regionsplit boolList blockedFace(mesh_.nFaces()); @@ -1652,8 +1604,6 @@ void Foam::meshRefinement::findCellZoneTopo } -// Make namedSurfaceIndex consistent with cellToZone - clear out any blocked -// faces inbetween same cell zone. void Foam::meshRefinement::makeConsistentFaceIndex ( const labelList& cellToZone, @@ -1843,7 +1793,21 @@ Foam::labelList Foam::meshRefinement::freeStandingBaffleFaces const labelList& faceNeighbour = mesh_.faceNeighbour(); - DynamicList<label> faceLabels(mesh_.nFaces()/20); + // We want to pick up the faces to orient. These faces come in + // two variants: + // - faces originating from stand-alone faceZones + // (these will most likely have no cellZone on either side so + // ownZone and neiZone both -1) + // - sticky-up faces originating from a 'bulge' in a outside of + // a cellZone. These will have the same cellZone on either side. + // How to orient these is not really clearly defined so do them + // same as stand-alone faceZone faces for now. (Normally these will + // already have been removed by the 'allowFreeStandingZoneFaces=false' + // default setting) + + // Note that argument neiCellZone will have -1 on uncoupled boundaries. + + DynamicList<label> faceLabels(mesh_.nFaces()/100); for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++) { @@ -1853,7 +1817,7 @@ Foam::labelList Foam::meshRefinement::freeStandingBaffleFaces // Free standing baffle? label ownZone = cellToZone[faceOwner[faceI]]; label neiZone = cellToZone[faceNeighbour[faceI]]; - if (max(ownZone, neiZone) == -1) + if (ownZone == neiZone) { faceLabels.append(faceI); } @@ -1872,13 +1836,14 @@ Foam::labelList Foam::meshRefinement::freeStandingBaffleFaces // Free standing baffle? label ownZone = cellToZone[faceOwner[faceI]]; label neiZone = neiCellZone[faceI-mesh_.nInternalFaces()]; - if (max(ownZone, neiZone) == -1) + if (ownZone == neiZone) { faceLabels.append(faceI); } } } } + return faceLabels.shrink(); } @@ -2135,7 +2100,6 @@ void Foam::meshRefinement::consistentOrientation } - DynamicList<label> changedEdges; DynamicList<patchFaceOrientation> changedInfo; @@ -2224,7 +2188,6 @@ void Foam::meshRefinement::consistentOrientation } - // Walk PatchEdgeFaceWave < @@ -2326,7 +2289,6 @@ void Foam::meshRefinement::consistentOrientation // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // -// Split off unreachable areas of mesh. void Foam::meshRefinement::baffleAndSplitMesh ( const bool doHandleSnapProblems, @@ -2517,7 +2479,6 @@ void Foam::meshRefinement::baffleAndSplitMesh } -// Split off (with optional buffer layers) unreachable areas of mesh. Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::splitMesh ( const label nBufferLayers, @@ -2804,8 +2765,6 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::splitMesh } -// Find boundary points that connect to more than one cell region and -// split them. Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::dupNonManifoldPoints ( const localPointRegion& regionSide @@ -2859,8 +2818,6 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::dupNonManifoldPoints } -// Find boundary points that connect to more than one cell region and -// split them. Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::dupNonManifoldPoints() { // Analyse which points need to be duplicated @@ -2870,7 +2827,6 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::dupNonManifoldPoints() } -// Zoning Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::zonify ( const point& keepPoint, @@ -3172,22 +3128,23 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::zonify // Get coupled neighbour cellZone. Set to -1 on non-coupled patches. - labelList neiCellZone(mesh_.nFaces()-mesh_.nInternalFaces(), -1); + labelList neiCellZone; + syncTools::swapBoundaryCellList(mesh_, cellToZone, neiCellZone); forAll(patches, patchI) { const polyPatch& pp = patches[patchI]; - if (pp.coupled()) + if (!pp.coupled()) { + label bFaceI = pp.start()-mesh_.nInternalFaces(); forAll(pp, i) { - label faceI = pp.start()+i; - neiCellZone[faceI-mesh_.nInternalFaces()] = - cellToZone[mesh_.faceOwner()[faceI]]; + neiCellZone[bFaceI++] = -1; } } } - syncTools::swapBoundaryFaceList(mesh_, neiCellZone); + + // Get per face whether is it master (of a coupled set of faces) const PackedBoolList isMasterFace(syncTools::getMasterFaces(mesh_)); @@ -3228,6 +3185,13 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::zonify Info<< "Detected " << nFreeStanding << " free-standing zone faces" << endl; + if (debug) + { + OBJstream str(mesh_.time().path()/"freeStanding.obj"); + str.write(patch.localFaces(), patch.localPoints(), false); + } + + // Detect non-manifold edges labelList nMasterFacesPerEdge; calcPatchNumMasterFaces(isMasterFace, patch, nMasterFacesPerEdge); @@ -3326,22 +3290,19 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::zonify if (surfI != -1) { // Orient face zone to have slave cells in max cell zone. + // Note: logic to use flipMap should be consistent with logic + // to pick up the freeStandingBaffleFaces! + label ownZone = cellToZone[faceOwner[faceI]]; label neiZone = cellToZone[faceNeighbour[faceI]]; bool flip; - label maxZone = max(ownZone, neiZone); - - if (maxZone == -1) + if (ownZone == neiZone) { // free-standing face. Use geometrically derived orientation flip = meshFlipMap[faceI]; } - else if (ownZone == maxZone) - { - flip = false; - } else { flip = true; @@ -3384,26 +3345,18 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::zonify bool flip; - label maxZone = max(ownZone, neiZone); - - if (maxZone == -1) + if (ownZone == neiZone) { // free-standing face. Use geometrically derived orientation flip = meshFlipMap[faceI]; } - else if (ownZone == neiZone) - { - // Free-standing zone face or coupled boundary. Keep master - // face unflipped. - flip = !isMasterFace[faceI]; - } - else if (neiZone == maxZone) - { - flip = true; - } else { - flip = false; + flip = + ( + ownZone == -1 + || (neiZone != -1 && ownZone > neiZone) + ); } meshMod.setAction