48 void Foam::meshRefinement::markBoundaryFace
56 isBoundaryFace[facei] =
true;
62 isBoundaryEdge[fEdges[fp]] =
true;
65 const face& f = mesh_.
faces()[facei];
69 isBoundaryPoint[f[fp]] =
true;
74 void Foam::meshRefinement::findNearest
77 List<pointIndexHit>& nearestInfo,
101 nearestNormal.setSize(nearestInfo.size());
102 nearestRegion.setSize(nearestInfo.size());
104 forAll(allSurfaces, surfI)
106 DynamicList<pointIndexHit> localHits;
110 if (nearestSurface[i] == surfI)
112 localHits.append(nearestInfo[i]);
119 surfaces_.
geometry()[geomI].getNormal(localHits, localNormals);
122 surfaces_.
geometry()[geomI].getRegion(localHits, localRegion);
127 if (nearestSurface[i] == surfI)
129 nearestNormal[i] = localNormals[localI];
130 nearestRegion[i] = localRegion[localI];
145 autoPtr<indirectPrimitivePatch> ppPtr
159 DynamicList<label> candidateFaces(pp.size()/20);
167 const labelList& eFaces = edgeFaces[edgeI];
169 if (eFaces.size() == 2)
171 label face0 = pp.addressing()[eFaces[0]];
172 label face1 = pp.addressing()[eFaces[1]];
177 if (cellLevel[cell0] > cellLevel[cell1])
180 const vector& n0 = pp.faceNormals()[eFaces[0]];
181 const vector& n1 = pp.faceNormals()[eFaces[1]];
183 if (
mag(n0 & n1) < 0.1)
185 candidateFaces.append(face0);
188 else if (cellLevel[cell1] > cellLevel[cell0])
191 const vector& n0 = pp.faceNormals()[eFaces[0]];
192 const vector& n1 = pp.faceNormals()[eFaces[1]];
194 if (
mag(n0 & n1) < 0.1)
196 candidateFaces.append(face1);
201 candidateFaces.shrink();
204 <<
" faces on edge-connected cells of differing level." 207 if (debug&meshRefinement::MESH)
209 faceSet fSet(mesh_,
"edgeConnectedFaces", candidateFaces);
211 Pout<<
"Writing " << fSet.size()
212 <<
" with problematic topology to faceSet " 213 << fSet.objectPath() <<
endl;
221 List<pointIndexHit> nearestInfo;
238 Map<label> candidateCells(candidateFaces.size());
240 faceSet perpFaces(mesh_,
"perpendicularFaces", pp.size()/100);
244 label facei = candidateFaces[i];
255 scalar angle =
degToRad(perpendicularAngle[region]);
261 perpFaces.insert(facei);
262 candidateCells.insert
265 globalToPatch[region]
271 if (debug&meshRefinement::MESH)
274 Pout<<
"Writing " << perpFaces.size()
275 <<
" faces that are perpendicular to the surface to set " 276 << perpFaces.objectPath() <<
endl;
279 return candidateCells;
286 bool Foam::meshRefinement::isCollapsedFace
290 const scalar minFaceArea,
291 const scalar maxNonOrtho,
296 const scalar severeNonorthogonalityThreshold =
299 const vector s = mesh_.
faces()[facei].area(points);
300 const scalar magS =
mag(s);
303 if (magS < minFaceArea)
316 scalar dDotS = (d &
s)/(
mag(d)*magS + vSmall);
318 if (dDotS < severeNonorthogonalityThreshold)
335 scalar dDotS = (d &
s)/(
mag(d)*magS + vSmall);
337 if (dDotS < severeNonorthogonalityThreshold)
357 bool Foam::meshRefinement::isCollapsedCell
360 const scalar volFraction,
364 scalar vol = mesh_.
cells()[celli].mag(points, mesh_.
faces());
386 if (adaptPatchIDs.size())
395 const polyPatch& pp = patches[adaptPatchIDs[i]];
400 List<topoDistanceData> cellData(mesh_.
nCells());
401 List<topoDistanceData> faceData(mesh_.
nFaces());
405 List<topoDistanceData> patchData(nFaces);
409 label patchi = adaptPatchIDs[i];
410 const polyPatch& pp = patches[
patchi];
414 patchFaces[nFaces] = pp.start()+i;
415 patchData[nFaces] = topoDistanceData(patchi, 0);
421 FaceCellWave<topoDistanceData> deltaCalc
433 bool haveWarned =
false;
436 if (!faceData[facei].valid(deltaCalc.data()))
441 <<
"Did not visit some faces, e.g. face " << facei
443 <<
"Assigning these cells to patch " 451 nearestAdaptPatch[facei] = faceData[facei].data();
458 nearestAdaptPatch.setSize(mesh_.
nFaces(), 0);
461 return nearestAdaptPatch;
473 const dictionary& motionDict,
474 const bool removeEdgeConnectedCells,
496 const polyPatch& pp = patches[adaptPatchIDs[i]];
498 label facei = pp.start();
516 const labelList nearestAdaptPatch(nearestPatch(adaptPatchIDs));
526 calcNeighbourData(neiLevel, neiCc);
530 label nBaffleFaces = 0;
534 label nPrevented = 0;
536 if (removeEdgeConnectedCells &&
max(perpendicularAngle) >= 0)
538 Info<<
"markFacesOnProblemCells :" 539 <<
" Checking for edge-connected cells of highly differing sizes." 544 Map<label> problemCells
546 findEdgeConnectedProblemCells
556 const cell& cFaces = mesh_.
cells()[iter.key()];
560 label facei = cFaces[i];
564 facePatch[facei] = nearestAdaptPatch[facei];
578 Info<<
"markFacesOnProblemCells : Marked " 580 <<
" additional internal faces to be converted into baffles" 583 <<
" cells edge-connected to lower level cells." <<
endl;
585 if (debug&meshRefinement::MESH)
587 cellSet problemCellSet(mesh_,
"problemCells", problemCells.toc());
588 problemCellSet.instance() =
timeName();
589 Pout<<
"Writing " << problemCellSet.size()
590 <<
" cells that are edge connected to coarser cell to set " 591 << problemCellSet.objectPath() <<
endl;
592 problemCellSet.
write();
624 const scalar volFraction =
625 motionDict.lookupOrDefault<scalar>(
"minVolCollapseRatio", -1);
627 const bool checkCollapse = (volFraction > 0);
629 scalar maxNonOrtho = -1;
637 minArea = motionDict.lookup<scalar>(
"minArea");
638 maxNonOrtho = motionDict.lookup<scalar>(
"maxNonOrtho");
640 Info<<
"markFacesOnProblemCells :" 641 <<
" Deleting all-anchor surface cells only if" 642 <<
" snapping them violates mesh quality constraints:" <<
nl 643 <<
" snapped/original cell volume < " << volFraction <<
nl 644 <<
" face area < " << minArea <<
nl 645 <<
" non-orthogonality > " << maxNonOrtho <<
nl 649 autoPtr<indirectPrimitivePatch> ppPtr
658 const pointField& localPoints = pp.localPoints();
659 const labelList& meshPoints = pp.meshPoints();
661 List<pointIndexHit> hitInfo;
673 newPoints = mesh_.
points();
677 if (hitInfo[i].hit())
679 newPoints[meshPoints[i]] = hitInfo[i].hitPoint();
683 if (debug&meshRefinement::MESH)
685 const_cast<Time&
>(mesh_.
time())++;
714 PackedBoolList hasSevenBoundaryAnchorPoints(mesh_.
nCells());
719 DynamicList<label> dynFEdges;
720 DynamicList<label> dynCPoints;
728 label nBoundaryAnchors = 0;
729 label nNonAnchorBoundary = 0;
730 label nonBoundaryAnchor = -1;
734 label pointi = cPoints[i];
736 if (pointLevel[pointi] <= cellLevel[celli])
739 if (isBoundaryPoint[pointi])
746 nonBoundaryAnchor = pointi;
749 else if (isBoundaryPoint[pointi])
751 nNonAnchorBoundary++;
755 if (nBoundaryAnchors == 8)
757 const cell& cFaces = mesh_.
cells()[celli];
764 if (isBoundaryFace[cFaces[cFacei]])
780 && !isCollapsedCell(newPoints, volFraction, celli)
797 label facei = cFaces[cf];
801 facePatch[facei] == -1
805 facePatch[facei] = nearestAdaptPatch[facei];
821 else if (nBoundaryAnchors == 7)
824 hasSevenBoundaryAnchorPoints.set(celli, 1u);
825 nonBoundaryAnchors.insert(nonBoundaryAnchor);
834 DynamicList<label> dynPCells;
838 label pointi = iter.key();
847 if (hasSevenBoundaryAnchorPoints.get(pCells[i]) == 1u)
858 label celli = pCells[i];
860 if (hasSevenBoundaryAnchorPoints.get(celli) == 1u)
865 && !isCollapsedCell(newPoints, volFraction, celli)
880 const cell& cFaces = mesh_.
cells()[celli];
884 label facei = cFaces[cf];
888 facePatch[facei] == -1
892 facePatch[facei] = nearestAdaptPatch[facei];
942 if (facePatch[facei] == -1)
945 label nFaceBoundaryEdges = 0;
949 if (isBoundaryEdge[fEdges[fe]])
951 nFaceBoundaryEdges++;
955 if (nFaceBoundaryEdges == fEdges.size())
979 facePatch[facei] = nearestAdaptPatch[facei];
991 const polyPatch& pp = patches[
patchi];
995 label facei = pp.start();
999 if (facePatch[facei] == -1)
1002 label nFaceBoundaryEdges = 0;
1006 if (isBoundaryEdge[fEdges[fe]])
1008 nFaceBoundaryEdges++;
1012 if (nFaceBoundaryEdges == fEdges.size())
1036 facePatch[facei] = nearestAdaptPatch[facei];
1037 if (isMasterFace[facei])
1066 Info<<
"markFacesOnProblemCells : marked " 1068 <<
" additional internal faces to be converted into baffles." 1073 Info<<
"markFacesOnProblemCells : prevented " 1075 <<
" internal faces from getting converted into baffles." 1085 Foam::labelList Foam::meshRefinement::markFacesOnProblemCellsGeometric
1087 const snapParameters& snapParams,
1088 const dictionary& motionDict
1098 autoPtr<indirectPrimitivePatch> ppPtr
1116 Info<<
"Constructing mesh displacer ..." <<
endl;
1117 Info<<
"Using mesh parameters " << motionDict <<
nl <<
endl;
1121 motionSmoother meshMover
1132 Info<<
"Checking initial mesh ..." <<
endl;
1141 Info<<
"Detected " << nInitErrors <<
" illegal faces" 1142 <<
" (concave, zero area or negative cell pyramid volume)" 1146 Info<<
"Checked initial mesh in = " 1173 const labelList& meshPoints = pp.meshPoints();
1178 newPoints[meshPoints[i]] += disp[i];
1185 minMagSqrEqOp<point>(),
1186 vector(great, great, great)
1200 label nBaffleFaces = 0;
1203 faceSet wrongFaces(mesh_,
"wrongFaces", 100);
1211 label nWrongFaces = 0;
1242 scalar minArea(motionDict.lookup<scalar>(
"minArea"));
1243 if (minArea > -small)
1261 Info<<
" faces with area < " 1262 <<
setw(5) << minArea
1264 << nNewWrongFaces-nWrongFaces <<
endl;
1266 nWrongFaces = nNewWrongFaces;
1269 scalar minDet(motionDict.lookup<scalar>(
"minDeterminant"));
1288 Info<<
" faces on cells with determinant < " 1289 <<
setw(5) << minDet <<
" : " 1290 << nNewWrongFaces-nWrongFaces <<
endl;
1292 nWrongFaces = nNewWrongFaces;
1303 facePatch[iter.key()] = nearestAdaptPatch[iter.key()];
1319 Info<<
"markFacesOnProblemCellsGeometric : marked " 1321 <<
" additional internal and coupled faces" 1322 <<
" to be converted into baffles." <<
endl;
static bool checkMesh(const bool report, const polyMesh &mesh, const dictionary &dict, labelHashSet &wrongFaces)
Check mesh with mesh settings in dict. Collects incorrect faces.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
List< labelList > labelListList
A List of labelList.
#define forAll(list, i)
Loop across all elements in list.
fileName path() const
Return path.
virtual Ostream & write(const char)=0
Write character.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
const labelList & surfaces() const
bool checkFaceArea(const bool report, const scalar minArea, const polyMesh &, const vectorField &faceAreas, const labelList &checkFaces, labelHashSet *setPtr)
Check for small faces.
const labelIOList & pointLevel() const
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
const labelListList & faceEdges() const
label nInternalFaces() const
virtual const labelList & faceNeighbour() const
Return face neighbour.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Unit conversion functions.
void size(const label)
Override size to be inconsistent with allocated storage.
Ostream & endl(Ostream &os)
Add newline and flush stream.
label nTotalCells() const
Return total number of cells in decomposed mesh.
static writeType writeLevel()
Get/set write level.
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
const cellList & cells() const
Vector< scalar > vector
A scalar version of the templated Vector.
const Time & time() const
Return the top-level database.
PrimitivePatch< IndirectList< face >, const pointField & > indirectPrimitivePatch
Foam::indirectPrimitivePatch.
void findNearest(const labelList &surfacesToTest, const pointField &samples, const scalarField &nearestDistSqr, labelList &surfaces, List< pointIndexHit > &) const
Find nearest point on surfaces.
scalar degToRad(const scalar deg)
Conversion from degrees to radians.
static const pointMesh & New(const polyMesh &mesh)
virtual const pointField & points() const
Return raw points.
List< bool > boolList
Bool container classes.
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
const searchableSurfaces & geometry() const
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
vectorField pointField
pointField is a vectorField.
dimensionedScalar cos(const dimensionedScalar &ds)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
double cpuTimeIncrement() const
Return CPU time (in seconds) since last call to cpuTimeIncrement()
PolyMesh checking routines. Checks various criteria for a mesh and supplied geometry, with the mesh only used for topology. Coupled patch aware (i.e., coupled faces are treated as internal).
virtual const labelList & faceOwner() const
Return face owner.
const globalMeshData & globalData() const
Return parallel info.
List< label > labelList
A List of labels.
static void preSmoothPatch(const meshRefinement &meshRefiner, const snapParameters &snapParams, const label nInitErrors, const List< labelPair > &baffles, motionSmoother &)
Smooth the mesh (patch and internal) to increase visibility.
virtual const faceList & faces() const
Return raw faces.
const PtrList< surfaceZonesInfo > & surfZones() const
const vectorField & cellCentres() const
static labelList getUnnamedSurfaces(const PtrList< surfaceZonesInfo > &surfList)
Get indices of unnamed surfaces (surfaces without faceZoneName)
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
const labelListList & pointCells() const
word timeName() const
Replacement for Time::timeName() : return oldInstance (if.
static scalarField calcSnapDistance(const fvMesh &mesh, const snapParameters &snapParams, const indirectPrimitivePatch &)
Calculate edge length per patch point.
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
dimensionedScalar sin(const dimensionedScalar &ds)
Istream and Ostream manipulators taking arguments.
labelList meshedPatches() const
Get patchIDs for patches added in addMeshedPatch.
static autoPtr< indirectPrimitivePatch > makePatch(const polyMesh &, const labelList &)
Create patch from set of patches.
virtual tmp< scalarField > movePoints(const pointField &)
Move points, returns volumes swept by faces in motion.
const vectorField & faceCentres() const
void setSize(const label)
Reset size of List.
bool checkCellDeterminant(const bool report, const scalar minDet, const polyMesh &, const vectorField &faceAreas, const labelList &checkFaces, labelHashSet *setPtr)
Check the area of internal faces v.s. boundary faces.
vector point
Point is a vector.
#define WarningInFunction
Report a warning using Foam::Warning.
prefixOSstream Pout(cout, "Pout")
const vectorField & faceAreas() const
label globalRegion(const label surfI, const label regionI) const
From surface and region on surface to global region.
dimensioned< scalar > mag(const dimensioned< Type > &)
Field< vector > vectorField
Specialisation of Field<T> for vector.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
Omanip< int > setw(const int i)
bool write() const
Write mesh and all data.
static tmp< pointVectorField > makeDisplacementField(const pointMesh &pMesh, const labelList &adaptPatchIDs)
Helper function to make a pointVectorField with correct.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
const labelIOList & cellLevel() const
const labelListList & cellPoints() const
const scalarField & cellVolumes() const
A HashTable to objects of type <T> with a label key.
static vectorField calcNearestSurface(const meshRefinement &meshRefiner, const scalarField &snapDist, const indirectPrimitivePatch &, pointField &nearestPoint, vectorField &nearestNormal)
Per patch point calculate point on nearest surface. Set as.