FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

hexRef8 Class Reference

Refinement of (split) hexes using polyTopoChange. More...

#include <dynamicMesh/hexRef8.H>


Detailed Description

Refinement of (split) hexes using polyTopoChange.

Source files

Definition at line 63 of file hexRef8.H.

Collaboration diagram for hexRef8:

List of all members.

Public Member Functions

 ClassName ("hexRef8")
 Runtime type information.
 hexRef8 (const polyMesh &mesh)
 Construct from mesh, read_if_present refinement data.
 hexRef8 (const polyMesh &mesh, const labelList &cellLevel, const labelList &pointLevel, const refinementHistory &history)
 Construct from mesh and un/refinement data.
 hexRef8 (const polyMesh &mesh, const labelList &cellLevel, const labelList &pointLevel)
 Construct from mesh and refinement data.
const labelIOList &  cellLevel () const
const labelIOList &  pointLevel () const
const refinementHistory &  history () const
scalar  level0EdgeLength () const
 Typical edge length between unrefined points.
label  getAnchorLevel (const label faceI) const
 Gets level such that the face has four points <= level.
labelList  consistentRefinement (const labelList &cellsToRefine, const bool maxSet) const
 Given valid mesh and current cell level and proposed.
labelList  consistentSlowRefinement (const label maxFaceDiff, const labelList &cellsToRefine, const labelList &facesToCheck, const label maxPointDiff, const labelList &pointsToCheck) const
 Like consistentRefinement but slower:
labelList  consistentSlowRefinement2 (const label maxFaceDiff, const labelList &cellsToRefine, const labelList &facesToCheck) const
 Like consistentSlowRefinement but uses different meshWave.
labelListList  setRefinement (const labelList &cells, polyTopoChange &)
 Insert refinement. All selected cells will be split into 8.
void  updateMesh (const mapPolyMesh &)
 Update local numbering for changed mesh.
void  storeData (const labelList &pointsToStore, const labelList &facesToStore, const labelList &cellsToStore)
 Signal points/face/cells for which to store data.
void  updateMesh (const mapPolyMesh &, const Map< label > &pointsToRestore, const Map< label > &facesToRestore, const Map< label > &cellsToRestore)
 Update local numbering + undo.
void  subset (const labelList &pointMap, const labelList &faceMap, const labelList &cellMap)
 Update local numbering for subsetted mesh.
void  distribute (const mapDistributePolyMesh &)
 Update local numbering for mesh redistribution.
void  checkMesh () const
 Debug: Check coupled mesh for correctness.
void  checkRefinementLevels (const label maxPointDiff, const labelList &pointsToCheck) const
 Debug: Check 2:1 consistency across faces.
labelList  getSplitPoints () const
 Return the points at the centre of top-level split cells.
labelList  consistentUnrefinement (const labelList &pointsToUnrefine, const bool maxSet) const
 Given proposed.
void  setUnrefinement (const labelList &splitPointLabels, polyTopoChange &)
 Remove some refinement. Needs to be supplied output of.
void  setInstance (const fileName &inst)
bool  write () const
 Force writing refinement+history to polyMesh directory.

Constructor & Destructor Documentation

hexRef8 ( const polyMesh &   mesh,
const labelList &   cellLevel,
const labelList &   pointLevel,
const refinementHistory &   history  
)

Construct from mesh and un/refinement data.

Definition at line 1849 of file hexRef8.C.

References Foam::abort(), Foam::endl(), Foam::FatalError, and FatalErrorIn.

hexRef8 ( const polyMesh &   mesh,
const labelList &   cellLevel,
const labelList &   pointLevel  
)

Construct from mesh and refinement data.

Definition at line 1946 of file hexRef8.C.

References Foam::abort(), Foam::endl(), Foam::FatalError, and FatalErrorIn.


Member Function Documentation

ClassName ( "hexRef8"    )

Runtime type information.

const labelIOList& cellLevel (  ) const [inline]

Definition at line 347 of file hexRef8.H.

Referenced by dynamicRefineFvMesh::dynamicRefineFvMesh().

const labelIOList& pointLevel (  ) const [inline]

Definition at line 352 of file hexRef8.H.

Referenced by dynamicRefineFvMesh::dynamicRefineFvMesh().

const refinementHistory& history (  ) const [inline]

Definition at line 357 of file hexRef8.H.

Referenced by dynamicRefineFvMesh::update().

scalar level0EdgeLength (  ) const [inline]

Typical edge length between unrefined points.

Definition at line 363 of file hexRef8.H.

Foam::label getAnchorLevel ( const label   faceI  ) const

Gets level such that the face has four points <= level.

Definition at line 798 of file hexRef8.C.

References List< T >::size().

Foam::labelList consistentRefinement ( const labelList &   cellsToRefine,
const bool   maxSet  
) const

Given valid mesh and current cell level and proposed.

cells to refine calculate any clashes (due to 2:1) and return ok list of cells to refine. Either adds cells to refine to set (maxSet = true) or removes cells to refine (maxSet = false)

Definition at line 2031 of file hexRef8.C.

References Foam::endl(), forAll, Foam::Pout, and Foam::reduce().

Foam::labelList consistentSlowRefinement ( const label   maxFaceDiff,
const labelList &   cellsToRefine,
const labelList &   facesToCheck,
const label   maxPointDiff,
const labelList &   pointsToCheck  
) const

Like consistentRefinement but slower:

  • specify number of cells between consecutive refinement levels (consistentRefinement equivalent to 1)

specify max level difference between point-connected cells. (-1 to disable) Note that with normal 2:1 limitation (maxFaceDiff=1) there can be 8:1 size difference across point connected cells so maxPointDiff allows you to make that less. cellsToRefine : cells we're thinking about refining. It will extend this set. All refinement levels will be at least maxFaceDiff layers thick. facesToCheck : additional faces where to implement the maxFaceDiff thickness (usually only boundary faces)

Definition at line 2106 of file hexRef8.C.

References Foam::abort(), List< T >::append(), refinementData::count(), Foam::endl(), Foam::exit(), Foam::FatalError, FatalErrorIn, Foam::findIndices(), forAll, forAllConstIter, refinementData::isRefined(), FaceCellWave< Type >::iterate(), Foam::mag(), Foam::max(), Foam::nl, Foam::Pout, Foam::reduce(), refinementData::refinementCount(), FaceCellWave< Type >::setFaceInfo(), List< T >::size(), syncTools::swapBoundaryFaceList(), syncTools::syncPointList(), and refinementData::updateFace().

Foam::labelList consistentSlowRefinement2 ( const label   maxFaceDiff,
const labelList &   cellsToRefine,
const labelList &   facesToCheck  
) const

Like consistentSlowRefinement but uses different meshWave.

(proper distance instead of toplogical count). No point checks yet.

Definition at line 2599 of file hexRef8.C.

References Foam::abort(), List< T >::append(), Foam::endl(), Foam::exit(), Foam::FatalError, FatalErrorIn, Foam::findIndices(), forAll, PackedList< nBits >::get(), Foam::nl, IOobject::objectPath(), Foam::Pout, Foam::reduce(), HashTable< T, Key, Hash >::size(), List< T >::size(), Ostream::write(), and regIOobject::write().

Foam::labelListList setRefinement ( const labelList &   cells,
polyTopoChange &   meshMod  
)

Insert refinement. All selected cells will be split into 8.

Returns per element in cells the 8 cells they were split into. Guarantees that the 0th element is the original cell label. Mapping: -split cells: 7 new ones get added from original -split faces: original gets modified; 3 new ones get added from original -added internal faces: added from original cell face(if that was internal) or created out-of-nothing (so will not get mapped!). Note: could make this inflate from point but that will allocate interpolation. -points added to split edge: added from edge start() -midpoints added: added from cellPoints[0].

Definition at line 3025 of file hexRef8.C.

References Foam::abort(), List< T >::append(), DynamicList< T, SizeInc, SizeMult, SizeDiv >::append(), edge::centre(), Foam::endl(), Foam::FatalError, FatalErrorIn, forAll, Foam::max(), Foam::min(), OSstream::name(), Foam::nl, Foam::Pout, polyTopoChange::setAction(), List< T >::setSize(), List< T >::size(), syncTools::swapBoundaryFaceList(), syncTools::syncBoundaryFaceList(), syncTools::syncEdgeList(), syncTools::syncFaceList(), List< T >::transfer(), Ostream::write(), and Foam::meshTools::writeOBJ().

void updateMesh ( const mapPolyMesh &   map  )

Update local numbering for changed mesh.

Definition at line 4080 of file hexRef8.C.

void storeData ( const labelList &   pointsToStore,
const labelList &   facesToStore,
const labelList &   cellsToStore  
)

Signal points/face/cells for which to store data.

Definition at line 4055 of file hexRef8.C.

References forAll, List< T >::resize(), and List< T >::size().

void updateMesh ( const mapPolyMesh &   map,
const Map< label > &   pointsToRestore,
const Map< label > &   facesToRestore,
const Map< label > &   cellsToRestore  
)
void subset ( const labelList &   pointMap,
const labelList &   faceMap,
const labelList &   cellMap  
)

Update local numbering for subsetted mesh.

Gets new-to-old maps. Not compatible with unrefinement.

Definition at line 4279 of file hexRef8.C.

References Foam::abort(), Foam::endl(), Foam::FatalError, FatalErrorIn, Foam::findIndex(), forAll, Foam::nl, Foam::Pout, List< T >::size(), and WarningIn.

void distribute ( const mapDistributePolyMesh &   map  )

Update local numbering for mesh redistribution.

Definition at line 4362 of file hexRef8.C.

References mapDistributePolyMesh::distributeCellData(), mapDistributePolyMesh::distributePointData(), Foam::endl(), and Foam::Pout.

void checkRefinementLevels ( const label   maxPointDiff,
const labelList &   pointsToCheck  
) const

Debug: Check 2:1 consistency across faces.

maxPointDiff==-1 : only check 2:1 across faces maxPointDiff!=-1 : check point-connected cells.

Gives problems after first splitting off inside mesher.

// Hanging points { Any patches with points having only two edges.

boolList isHangingPoint(mesh_.nPoints(), false);

const polyBoundaryMesh& patches = mesh_.boundaryMesh();

forAll(patches, patchI) { const polyPatch& pp = patches[patchI];

const labelList& meshPoints = pp.meshPoints();

forAll(meshPoints, i) { label pointI = meshPoints[i];

const labelList& pEdges = mesh_.pointEdges()[pointI];

if (pEdges.size() == 2) { isHangingPoint[pointI] = true; } } }

syncTools::syncPointList ( mesh_, isHangingPoint, andEqOp<bool>(), // only if all decide it is hanging point true, // null false // no separation );

OFstream str(mesh_.time().path()/"hangingPoints.obj");

label nHanging = 0;

forAll(isHangingPoint, pointI) { if (isHangingPoint[pointI]) { nHanging++;

Pout<< "Hanging boundary point " << pointI << " at " << mesh_.points()[pointI] << endl; meshTools::writeOBJ(str, mesh_.points()[pointI]); } }

if (returnReduce(nHanging, sumOp<label>()) > 0) { FatalErrorIn ( "hexRef8::checkRefinementLevels(const label)" ) << "Detected a point used by two edges only (hanging point)" << nl << "This is not allowed" << abort(FatalError); } }

Definition at line 4597 of file hexRef8.C.

References Foam::abort(), Foam::endl(), Foam::FatalError, FatalErrorIn, forAll, Foam::mag(), Foam::max(), Foam::nl, Foam::Pout, syncTools::swapBoundaryFaceList(), and syncTools::syncPointList().

Referenced by hexRef8::hexRef8().

Foam::labelList getSplitPoints (  ) const

Return the points at the centre of top-level split cells.

that can be unsplit.

Definition at line 4868 of file hexRef8.C.

References Foam::abort(), Foam::endl(), Foam::FatalError, FatalErrorIn, forAll, Foam::Pout, and List< T >::size().

Foam::labelList consistentUnrefinement ( const labelList &   pointsToUnrefine,
const bool   maxSet  
) const

Given proposed.

splitPoints to unrefine according to calculate any clashes (due to 2:1) and return ok list of points to unrefine. Either adds points to refine to set (maxSet = true) or removes points to refine (maxSet = false)

Definition at line 5072 of file hexRef8.C.

References Foam::abort(), Foam::endl(), Foam::FatalError, FatalErrorIn, forAll, Foam::Pout, Foam::reduce(), and syncTools::swapBoundaryFaceList().

void setUnrefinement ( const labelList &   splitPointLabels,
polyTopoChange &   meshMod  
)

Remove some refinement. Needs to be supplied output of.

consistentUnrefinement. Only call if undoable set. All 8 pointCells of a split point will be combined into the lowest numbered cell of those 8.

Definition at line 5301 of file hexRef8.C.

References Foam::abort(), Foam::endl(), Foam::FatalError, FatalErrorIn, forAll, Foam::min(), Foam::nl, IOobject::objectPath(), pFaces, Foam::Pout, HashTable< T, Key, Hash >::size(), List< T >::size(), and regIOobject::write().

void setInstance ( const fileName &   inst  )

Definition at line 1739 of file hexRef8.C.

References Foam::endl(), and Foam::Pout.

bool write (  ) const

Force writing refinement+history to polyMesh directory.

Definition at line 5498 of file hexRef8.C.


The documentation for this class was generated from the following files:
  • src/dynamicMesh/polyTopoChange/polyTopoChange/hexRef8.H
  • src/dynamicMesh/polyTopoChange/polyTopoChange/hexRef8.C