Go to the documentation of this file.
31 #include "surfaceInterpolate.H"
57 if (protectedCell_.empty())
59 unrefineableCell.
clear();
63 const labelList& cellLevel = meshCutter_.cellLevel();
65 unrefineableCell = protectedCell_;
68 labelList neiLevel(nFaces()-nInternalFaces());
70 for (label facei = nInternalFaces(); facei < nFaces(); ++facei)
72 neiLevel[facei-nInternalFaces()] = cellLevel[faceOwner()[facei]];
74 syncTools::swapBoundaryFaceList(*
this, neiLevel);
85 for (label facei = 0; facei < nInternalFaces(); ++facei)
87 const label own = faceOwner()[facei];
88 const label nei = faceNeighbour()[facei];
94 unrefineableCell.
test(own)
95 && (cellLevel[nei] > cellLevel[own])
100 unrefineableCell.
test(nei)
101 && (cellLevel[own] > cellLevel[nei])
108 for (label facei = nInternalFaces(); facei < nFaces(); facei++)
110 const label own = faceOwner()[facei];
116 unrefineableCell.
test(own)
117 && (neiLevel[facei-nInternalFaces()] > cellLevel[own])
129 bool hasExtended =
false;
131 for (label facei = 0; facei < nInternalFaces(); ++facei)
133 if (seedFace.
test(facei))
135 if (unrefineableCell.
set(faceOwner()[facei]))
139 if (unrefineableCell.
set(faceNeighbour()[facei]))
145 for (label facei = nInternalFaces(); facei < nFaces(); ++facei)
147 if (seedFace.
test(facei))
149 const label own = faceOwner()[facei];
151 if (unrefineableCell.
set(own))
188 for (
const auto& pr : fluxVelocities)
212 bitSet masterFaces(nFaces());
216 const label oldFacei =
faceMap[facei];
220 const label masterFacei = reverseFaceMap[oldFacei];
225 <<
"Problem: should not have removed faces"
227 <<
nl <<
"face:" << facei <<
endl
230 else if (masterFacei != facei)
232 masterFaces.
set(masterFacei);
239 Pout<<
"Found " << masterFaces.
count() <<
" split faces " <<
endl;
244 lookupClass<surfaceScalarField>()
254 if (!correctFluxes_.found(iter.key()))
257 <<
"Cannot find surfaceScalarField " << iter.key()
258 <<
" in user-provided flux mapping table "
259 << correctFluxes_ <<
endl
260 <<
" The flux mapping table is used to recreate the"
261 <<
" flux on newly created faces." <<
endl
262 <<
" Either add the entry if it is a flux or use ("
263 << iter.key() <<
" none) to suppress this warning."
268 const word&
UName = correctFluxes_[iter.key()];
279 Pout<<
"Setting surfaceScalarField " << iter.key()
280 <<
" to NaN" <<
endl;
289 Pout<<
"Mapping flux " << iter.key()
290 <<
" using interpolated flux " <<
UName
298 lookupObject<volVectorField>(
UName)
304 for (label facei = 0; facei < nInternalFaces(); ++facei)
306 const label oldFacei =
faceMap[facei];
311 phi[facei] = phiU[facei];
313 else if (reverseFaceMap[oldFacei] != facei)
316 phi[facei] = phiU[facei];
333 const label oldFacei =
faceMap[facei];
338 patchPhi[i] = patchPhiU[i];
340 else if (reverseFaceMap[oldFacei] != facei)
343 patchPhi[i] = patchPhiU[i];
351 for (
const label facei : masterFaces)
353 if (isInternalFace(facei))
355 phi[facei] = phiU[facei];
367 patchPhi[i] = patchPhiU[i];
379 mapNewInternalFaces<scalar>(this->Sf(), this->magSf(),
faceMap);
380 mapNewInternalFaces<vector>(this->Sf(), this->magSf(),
faceMap);
383 mapNewInternalFaces<sphericalTensor>(
faceMap);
384 mapNewInternalFaces<symmTensor>(
faceMap);
385 mapNewInternalFaces<tensor>(
faceMap);
401 meshCutter_.setRefinement(cellsToRefine, meshMod);
407 Info<<
"Refined from "
409 <<
" to " << globalData().nTotalCells() <<
" cells." <<
endl;
414 for (label facei = 0; facei < nInternalFaces(); ++facei)
416 const label oldFacei = map().faceMap()[facei];
418 if (oldFacei >= nInternalFaces())
421 <<
"New internal face:" << facei
422 <<
" fc:" << faceCentres()[facei]
423 <<
" originates from boundary oldFace:" << oldFacei
455 meshCutter_.updateMesh(*map);
458 if (protectedCell_.size())
460 bitSet newProtectedCell(nCells());
462 forAll(newProtectedCell, celli)
464 const label oldCelli = map().cellMap()[celli];
465 if (protectedCell_.test(oldCelli))
467 newProtectedCell.
set(celli);
470 protectedCell_.transfer(newProtectedCell);
474 meshCutter_.checkRefinementLevels(-1,
labelList());
489 meshCutter_.setUnrefinement(splitPoints, meshMod);
498 Map<label> faceToSplitPoint(3*splitPoints.size());
501 for (
const label pointi : splitPoints)
503 const labelList& pEdges = pointEdges()[pointi];
505 for (
const label edgei : pEdges)
507 const label otherPointi = edges()[edgei].otherVertex(pointi);
511 for (
const label facei :
pFaces)
513 faceToSplitPoint.insert(facei, otherPointi);
524 Info<<
"Unrefined from "
526 <<
" to " << globalData().nTotalCells() <<
" cells."
549 const labelList& reversePointMap = map().reversePointMap();
550 const labelList& reverseFaceMap = map().reverseFaceMap();
554 lookupClass<surfaceScalarField>()
558 if (!correctFluxes_.found(iter.key()))
561 <<
"Cannot find surfaceScalarField " << iter.key()
562 <<
" in user-provided flux mapping table "
563 << correctFluxes_ <<
endl
564 <<
" The flux mapping table is used to recreate the"
565 <<
" flux on newly created faces." <<
endl
566 <<
" Either add the entry if it is a flux or use ("
567 << iter.key() <<
" none) to suppress this warning."
572 const word&
UName = correctFluxes_[iter.key()];
580 <<
"Mapping flux " << iter.key()
581 <<
" using interpolated flux " <<
UName
586 surfaceScalarField::Boundary& phiBf =
593 lookupObject<volVectorField>(
UName)
601 const label oldFacei = iter.key();
602 const label oldPointi = iter.val();
604 if (reversePointMap[oldPointi] < 0)
607 const label facei = reverseFaceMap[oldFacei];
611 if (isInternalFace(facei))
613 phi[facei] = phiU[facei];
623 patchPhi[i] = patchPhiU[i];
633 meshCutter_.updateMesh(*map);
636 if (protectedCell_.size())
638 bitSet newProtectedCell(nCells());
640 forAll(newProtectedCell, celli)
642 const label oldCelli = map().cellMap()[celli];
643 if (protectedCell_.test(oldCelli))
645 newProtectedCell.
set(celli);
648 protectedCell_.transfer(newProtectedCell);
652 meshCutter_.checkRefinementLevels(-1,
labelList());
667 for (
const label celli : pCells)
669 vFld[celli] =
max(vFld[celli], pFld[pointi]);
685 for (
const label celli : pCells)
687 pFld[pointi] =
max(pFld[pointi], vFld[celli]);
704 for (
const label celli : pCells)
708 pFld[pointi] =
sum/pCells.size();
717 const scalar minLevel,
718 const scalar maxLevel
725 scalar err =
min(
fld[i]-minLevel, maxLevel-
fld[i]);
738 const scalar lowerRefineLevel,
739 const scalar upperRefineLevel,
762 if (cellError[celli] > 0)
764 candidateCell.
set(celli);
772 const label maxCells,
773 const label maxRefinement,
774 const bitSet& candidateCell
778 label nTotToRefine = (maxCells - globalData().nTotalCells()) / 7;
780 const labelList& cellLevel = meshCutter_.cellLevel();
785 calculateProtectedCells(unrefineableCell);
788 label nLocalCandidates = candidateCell.
count();
794 if (nCandidates < nTotToRefine)
796 for (
const label celli : candidateCell)
800 (!unrefineableCell.
test(celli))
801 && cellLevel[celli] < maxRefinement
811 for (label level = 0; level < maxRefinement; ++level)
813 for (
const label celli : candidateCell)
817 (!unrefineableCell.
test(celli))
818 && cellLevel[celli] == level
835 meshCutter_.consistentRefinement
843 <<
" cells for refinement out of " << globalData().nTotalCells()
846 return consistentSet;
852 const scalar unrefineLevel,
858 const labelList splitPoints(meshCutter_.getSplitPoints());
868 if (protectedCell_.size())
875 if (protectedCell_.test(celli))
877 protectedPoint.
set(pointi);
893 <<
" protected cells found "
895 <<
" protected points." <<
endl;
901 for (
const label pointi : splitPoints)
903 if (!protectedPoint[pointi] && pFld[pointi] < unrefineLevel)
906 bool hasMarked =
false;
910 if (markedCell.
test(celli))
919 newSplitPoints.
append(pointi);
925 newSplitPoints.shrink();
930 meshCutter_.consistentUnrefinement
937 <<
" split points out of a possible "
941 return consistentSet;
951 bitSet markedFace(nFaces());
953 for (
const label celli : markedCell)
961 for (label facei = 0; facei < nInternalFaces(); ++facei)
963 if (markedFace.
test(facei))
965 markedCell.set(faceOwner()[facei]);
966 markedCell.set(faceNeighbour()[facei]);
969 for (label facei = nInternalFaces(); facei < nFaces(); ++facei)
971 if (markedFace.
test(facei))
973 markedCell.set(faceOwner()[facei]);
984 const labelList& cellLevel = meshCutter_.cellLevel();
985 const labelList& pointLevel = meshCutter_.pointLevel();
989 forAll(pointLevel, pointi)
993 for (
const label celli : pCells)
995 if (pointLevel[pointi] <= cellLevel[celli])
998 if (nAnchorPoints[celli] == 8)
1000 protectedCell.
set(celli);
1003 if (!protectedCell.
test(celli))
1005 ++nAnchorPoints[celli];
1012 forAll(protectedCell, celli)
1014 if (nAnchorPoints[celli] != 8)
1016 protectedCell.
set(celli);
1024 Foam::dynamicRefineFvMesh::dynamicRefineFvMesh
1047 protectedCell_.setSize(nCells());
1048 nRefinementIterations_ = 0;
1055 const labelList& cellLevel = meshCutter_.cellLevel();
1056 const labelList& pointLevel = meshCutter_.pointLevel();
1072 for (
const label celli : pCells)
1074 if (!protectedCell_.test(celli))
1076 if (pointLevel[pointi] <= cellLevel[celli])
1080 if (nAnchors[celli] > 8)
1082 protectedCell_.
set(celli);
1098 for (label facei = 0; facei < nInternalFaces(); ++facei)
1100 neiLevel[facei] = cellLevel[faceNeighbour()[facei]];
1102 for (label facei = nInternalFaces(); facei < nFaces(); ++facei)
1104 neiLevel[facei] = cellLevel[faceOwner()[facei]];
1109 bitSet protectedFace(nFaces());
1111 forAll(faceOwner(), facei)
1113 const label faceLevel =
max
1115 cellLevel[faceOwner()[facei]],
1119 const face&
f = faces()[facei];
1123 for (
const label pointi :
f)
1125 if (pointLevel[pointi] <= faceLevel)
1131 protectedFace.
set(facei);
1140 for (label facei = 0; facei < nInternalFaces(); ++facei)
1142 if (protectedFace.
test(facei))
1144 protectedCell_.set(faceOwner()[facei]);
1145 protectedCell_.set(faceNeighbour()[facei]);
1148 for (label facei = nInternalFaces(); facei < nFaces(); ++facei)
1150 if (protectedFace.
test(facei))
1152 protectedCell_.set(faceOwner()[facei]);
1161 if (cFaces.size() < 6)
1163 protectedCell_.
set(celli);
1167 for (
const label cfacei : cFaces)
1169 if (faces()[cfacei].size() < 4)
1171 protectedCell_.set(celli);
1179 checkEightAnchorPoints(protectedCell_);
1184 protectedCell_.clear();
1197 <<
" cells that are protected from refinement."
1198 <<
" Writing these to cellSet "
1199 << protectedCells.
name()
1202 protectedCells.
write();
1229 ).optionalSubDict(typeName +
"Coeffs")
1232 const label refineInterval = refineDict.
get<label>(
"refineInterval");
1234 bool hasChanged =
false;
1236 if (refineInterval == 0)
1238 topoChanging(hasChanged);
1242 else if (refineInterval < 0)
1245 <<
"Illegal refineInterval " << refineInterval <<
nl
1246 <<
"The refineInterval setting in the dynamicMeshDict should"
1247 <<
" be >= 1." <<
nl
1257 const label maxCells = refineDict.
get<label>(
"maxCells");
1262 <<
"Illegal maximum number of cells " << maxCells <<
nl
1263 <<
"The maxCells setting in the dynamicMeshDict should"
1268 const label maxRefinement = refineDict.
get<label>(
"maxRefinement");
1270 if (maxRefinement <= 0)
1273 <<
"Illegal maximum refinement level " << maxRefinement <<
nl
1274 <<
"The maxCells setting in the dynamicMeshDict should"
1279 const word fieldName(refineDict.
get<
word>(
"field"));
1281 const volScalarField& vFld = lookupObject<volScalarField>(fieldName);
1283 const scalar lowerRefineLevel =
1284 refineDict.
get<scalar>(
"lowerRefineLevel");
1285 const scalar upperRefineLevel =
1286 refineDict.
get<scalar>(
"upperRefineLevel");
1287 const scalar unrefineLevel = refineDict.
getOrDefault<scalar>
1292 const label nBufferLayers = refineDict.
get<label>(
"nBufferLayers");
1298 selectRefineCandidates
1306 if (globalData().nTotalCells() < maxCells)
1325 if (nCellsToRefine > 0)
1333 const labelList& cellMap = map().cellMap();
1334 const labelList& reverseCellMap = map().reverseCellMap();
1336 bitSet newRefineCell(cellMap.size());
1340 const label oldCelli = cellMap[celli];
1345 || (reverseCellMap[oldCelli] != celli)
1349 newRefineCell.
set(celli);
1357 for (label i = 0; i < nBufferLayers; ++i)
1371 selectUnrefinePoints
1381 pointsToUnrefine.size(),
1385 if (nSplitPoints > 0)
1388 unrefine(pointsToUnrefine);
1395 if ((nRefinementIterations_ % 10) == 0)
1401 nRefinementIterations_++;
1404 topoChanging(hasChanged);
1423 const_cast<hexRef8&
>(meshCutter_).setInstance(time().timeName());
1428 && meshCutter_.write(valid)
1448 const labelList& cellLevel = meshCutter_.cellLevel();
1452 scalarCellLevel[celli] = cellLevel[celli];
1455 writeOk = writeOk && scalarCellLevel.write();
int debug
Static debugging option.
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
void clearOut()
Clear all geometry and addressing.
List< label > labelList
A List of labels.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
void reset()
Clear all bits but do not adjust the addressable size.
void resize(const label numElem, const unsigned int val=0u)
Reset addressable list size, does not shrink the allocated size.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
scalarField error(const scalarField &fld, const scalar minLevel, const scalar maxLevel) const
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
void extendMarkedCells(bitSet &markedCell) const
Extend markedCell with cell-face-cell.
A class for handling words, derived from Foam::string.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
static constexpr const zero Zero
Global zero (0)
virtual label start() const
Return start label of this patch in the polyMesh face list.
A topoSetPointSource to select all the points from given cellSet(s).
const word UName(propsDict.getOrDefault< word >("U", "U"))
virtual bool init(const bool doInit)
Initialise all non-demand-driven data.
virtual void mapFields(const mapPolyMesh &mpm)
Map all fields in time using given map.
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Direct mesh changes based on v1.3 polyTopoChange syntax.
virtual bool init(const bool doInit)
Initialise all non-demand-driven data.
Abstract base class for geometry and/or topology changing fvMesh.
virtual void mapFields(const mapPolyMesh &mpm)
Map all fields in time using given map.
virtual bool update()
Update the mesh for both mesh motion and topology change.
virtual labelList selectRefineCells(const label maxCells, const label maxRefinement, const bitSet &candidateCell) const
Subset candidate cells for refinement.
void set(const bitSet &bitset)
Set specified bits from another bitset.
Ostream & endl(Ostream &os)
Add newline and flush stream.
HashTable< word > correctFluxes_
Fluxes to map.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
bool test(const label pos) const
Test value at specified position, never auto-vivify entries.
DynamicList< T, SizeMin > & shrink()
Shrink the allocated space to the number of elements used.
unsigned int count(const bool on=true) const
Count number of bits set.
void readDict()
Read the projection parameters from dictionary.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
#define forAll(list, i)
Loop across all elements in list.
virtual bool writeObject(IOstreamOption streamOpt, const bool valid) const
Write using stream options.
All refinement history. Used in unrefinement.
void checkEightAnchorPoints(bitSet &protectedCell) const
Check all cells have 8 anchor points.
virtual bool write(const bool valid=true) const
Write using setting from DB.
static void fillNan(UList< scalar > &list)
Fill data block with NaN values.
messageStream Info
Information stream (stdout output on master, null elsewhere)
DynamicList< T, SizeMin > & append(const T &val)
Append an element to the end of this list.
virtual void selectRefineCandidates(const scalar lowerRefineLevel, const scalar upperRefineLevel, const scalarField &vFld, bitSet &candidateCell) const
Select candidate cells for refinement.
The IOstreamOption is a simple container for options an IOstream can normally have.
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
scalarField maxCellField(const volScalarField &) const
Get point max of connected cell.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< ' ';}gmvFile<< nl;for(const word &name :lagrangianScalarNames){ IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
autoPtr< mapPolyMesh > changeMesh(polyMesh &mesh, const labelUList &patchMap, const bool inflate, const bool syncParallel=true, const bool orderCells=false, const bool orderPoints=false)
Inplace changes mesh without change of patches.
#define forAllIters(container, iter)
Iterate across all elements in the container object.
virtual autoPtr< mapPolyMesh > refine(const labelList &)
Refine cells. Update mesh and fields.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Macros for easy insertion into run-time selection tables.
errorManip< error > abort(error &err)
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
A collection of cell labels.
virtual labelList selectUnrefinePoints(const scalar unrefineLevel, const bitSet &markedCell, const scalarField &pFld) const
Select points that can be unrefined.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
A HashTable similar to std::unordered_map.
const labelList & reverseFaceMap() const
Reverse face map.
const word & name() const noexcept
Return name.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
void calculateProtectedCells(bitSet &unrefineableCell) const
Calculate cells that cannot be refined since would trigger.
Container with cells to refine. Refinement given as single direction.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
#define DebugInfo
Report an information message using Foam::Info.
forAllConstIters(mixture.phases(), phase)
const fvPatch & patch() const
Return patch.
Refinement of (split) hexes using polyTopoChange.
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
labelHashSet used(const bitSet &select)
Convert a bitset to a labelHashSet of the indices used.
bool dumpLevel_
Dump cellLevel for post-processing.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const labelList & faceMap() const
Old face map.
const dictionary & optionalSubDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary, otherwise return this dictionary.
std::enable_if< std::is_same< bool, TypeT >::value, bool >::type set(const label i, bool val=true)
A bitSet::set() method for a list of bool.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
const dimensionedScalar c
Speed of light in a vacuum.
const Time & time() const
Return the top-level database.
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface....
virtual autoPtr< mapPolyMesh > unrefine(const labelList &)
Unrefine cells. Gets passed in centre points of cells to combine.
scalarField cellToPoint(const scalarField &vFld) const
A face is a list of labels corresponding to mesh vertices.
Smooth ATC in cells having a point to a set of patches supplied by type.
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
constant condensation/saturation model.
defineTypeNameAndDebug(combustionModel, 0)
scalarField maxPointField(const scalarField &) const
Get per cell max of connected point.
#define WarningInFunction
Report a warning using Foam::Warning.
A cell is defined as a list of faces with extra functionality.
Class to handle errors and exceptions in a simple, consistent stream-based manner.
Info<< "Finished reading KIVA file"<< endl;cellShapeList cellShapes(nPoints);labelList cellZoning(nPoints, -1);const cellModel &hex=cellModel::ref(cellModel::HEX);labelList hexLabels(8);label activeCells=0;labelList pointMap(nPoints);forAll(pointMap, i){ pointMap[i]=i;}for(label i=0;i< nPoints;i++){ if(f[i] > 0.0) { hexLabels[0]=i;hexLabels[1]=i1tab[i];hexLabels[2]=i3tab[i1tab[i]];hexLabels[3]=i3tab[i];hexLabels[4]=i8tab[i];hexLabels[5]=i1tab[i8tab[i]];hexLabels[6]=i3tab[i1tab[i8tab[i]]];hexLabels[7]=i3tab[i8tab[i]];cellShapes[activeCells].reset(hex, hexLabels);edgeList edges=cellShapes[activeCells].edges();forAll(edges, ei) { if(edges[ei].mag(points)< SMALL) { label start=pointMap[edges[ei].start()];while(start !=pointMap[start]) { start=pointMap[start];} label end=pointMap[edges[ei].end()];while(end !=pointMap[end]) { end=pointMap[end];} label minLabel=min(start, end);pointMap[start]=pointMap[end]=minLabel;} } cellZoning[activeCells]=idreg[i];activeCells++;}}cellShapes.setSize(activeCells);cellZoning.setSize(activeCells);forAll(cellShapes, celli){ cellShape &cs=cellShapes[celli];forAll(cs, i) { cs[i]=pointMap[cs[i]];} cs.collapse();}label bcIDs[11]={-1, 0, 2, 4, -1, 5, -1, 6, 7, 8, 9};const label nBCs=12;const word *kivaPatchTypes[nBCs]={ &wallPolyPatch::typeName, &wallPolyPatch::typeName, &wallPolyPatch::typeName, &wallPolyPatch::typeName, &symmetryPolyPatch::typeName, &wedgePolyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &symmetryPolyPatch::typeName, &oldCyclicPolyPatch::typeName};enum patchTypeNames{ PISTON, VALVE, LINER, CYLINDERHEAD, AXIS, WEDGE, INFLOW, OUTFLOW, PRESIN, PRESOUT, SYMMETRYPLANE, CYCLIC};const char *kivaPatchNames[nBCs]={ "piston", "valve", "liner", "cylinderHead", "axis", "wedge", "inflow", "outflow", "presin", "presout", "symmetryPlane", "cyclic"};List< SLList< face > > pFaces[nBCs]
void clear()
Clear the list, i.e. set addressable size to zero.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const dimensionSet dimless
Dimensionless.
virtual bool writeObject(IOstreamOption streamOpt, const bool valid) const
Write the underlying polyMesh and other data.