inverseDistance Class Reference

Inverse-distance-weighted interpolation stencil. More...

Inheritance diagram for inverseDistance:
Collaboration diagram for inverseDistance:

Classes

struct  interpolatedDonorSet
 Data for a set of interpolated/donor set. More...
 

Public Member Functions

 TypeName ("inverseDistance")
 Runtime type information. More...
 
 inverseDistance (const fvMesh &, const dictionary &, const bool)
 Construct from fvMesh. More...
 
virtual ~inverseDistance ()
 Destructor. More...
 
virtual bool update ()
 Update stencils. Return false if nothing changed. More...
 
virtual const labelUListcellTypes () const
 Return the cell type list. More...
 
virtual const labelUListinterpolationCells () const
 Indices of interpolated cells. More...
 
virtual const mapDistributecellInterpolationMap () const
 Return a communication schedule. More...
 
virtual const labelListListcellStencil () const
 Per interpolated cell the neighbour cells (in terms of slots as. More...
 
virtual const scalarListListcellInterpolationWeights () const
 Weights for cellStencil. More...
 
virtual const scalarListcellInterpolationWeight () const
 Per interpolated cell the interpolation factor. (0 = use. More...
 
virtual void stencilWeights (const point &sample, const pointList &donorCcs, scalarList &weights) const
 Calculate inverse distance weights for a single acceptor. More...
 
- Public Member Functions inherited from cellCellStencil
 TypeName ("cellCellStencil")
 Runtime type information. More...
 
 declareRunTimeSelectionTable (autoPtr, cellCellStencil, mesh,(const fvMesh &mesh, const dictionary &dict, const bool update),(mesh, dict, update))
 
 cellCellStencil (const fvMesh &)
 Construct from fvMesh. More...
 
virtual ~cellCellStencil ()
 Destructor. More...
 
virtual const wordHashSetnonInterpolatedFields () const
 Return the names of any (stencil or mesh specific) fields that. More...
 
virtual wordHashSetnonInterpolatedFields ()
 Return non-const non-interpolating fields. More...
 
bool localStencil (const labelUList &) const
 Helper: is stencil fully local. More...
 
const labelIOListzoneID () const
 Helper: get reference to registered zoneID. Loads volScalarField. More...
 
template<class T >
void interpolate (const fvMesh &mesh, Field< T > &psi) const
 Explicit interpolation of acceptor cells from donor cells. Looks up cellCellStencil. More...
 
template<class GeoField >
void interpolate (GeoField &psi) const
 Explicit interpolation of acceptor cells from donor cells with boundary condition handling. More...
 
template<class GeoField >
void interpolate (const fvMesh &mesh, const wordHashSet &suppressed) const
 Explicit interpolation of all registered fields. No boundary conditions. Excludes selected fields (and their old-time fields) More...
 
void walkFront (const globalIndex &globalCells, const scalar layerRelax, const labelListList &allStencil, labelList &allCellTypes, scalarField &allWeight, const scalarList &compactCellVol, const labelListList &compactStencil, const labelList &zoneID, const label holeLayers, const label useLayer) const
 Surround holes with layer(s) of interpolated cells. More...
 
void setUpFront (const labelList &allCellTypes, bitSet &isFront) const
 Set up front using allCellTypes. More...
 
void setUpFrontOnOversetPatch (const labelList &allCellTypes, bitSet &isFront) const
 Set up front on overset patches. More...
 
void seedCell (const label cellI, const scalar wantedFraction, bitSet &isFront, scalarField &fraction) const
 Seed faces of cell with wantedFraction (if higher than current) More...
 
InfoProxy< cellCellStencilinfo () const noexcept
 Return info proxy, used to print stencil information to a stream. More...
 
template<class Type >
Foam::tmp< Foam::volScalarFieldcreateField (const fvMesh &mesh, const word &name, const UList< Type > &psi)
 

Protected Member Functions

void markPatchesAsHoles (PstreamBuffers &pBufs, const PtrList< fvMeshSubset > &meshParts, const List< treeBoundBoxList > &patchBb, const List< labelVector > &patchDivisions, const PtrList< PackedList< 2 >> &patchParts, const label srcI, const label tgtI, labelList &allCellTypes) const
 Mark all cells overlapping (a voxel covered by) a src patch. More...
 
bool betterDonor (const label destMesh, const label currentDonorMesh, const label newDonorMesh) const
 If multiple donors meshes: decide which is best. More...
 
void markDonors (const globalIndex &globalCells, PstreamBuffers &pBufs, const PtrList< fvMeshSubset > &meshParts, const List< treeBoundBoxList > &meshBb, const labelList &allCellTypes, const label srcI, const label tgtI, labelListList &allStencil, labelList &allDonor) const
 Determine donors for all tgt cells. More...
 
void uncompactedRegionSplit (const fvMesh &mesh, const globalIndex &globalFaces, const label nZones, const labelList &zoneID, const labelList &cellTypes, const boolList &isBlockedFace, labelList &cellRegion) const
 Replacement of regionSplit. More...
 
autoPtr< globalIndexcompactedRegionSplit (const fvMesh &mesh, const globalIndex &globalFaces, labelList &cellRegion) const
 
void findHoles (const globalIndex &globalCells, const fvMesh &mesh, const labelList &zoneID, const labelListList &stencil, labelList &cellTypes) const
 Do flood filling to detect unreachable (from patches) sections. More...
 
void seedCell (const label cellI, const scalar wantedFraction, bitSet &isFront, scalarField &fraction) const
 Seed faces of cell with wantedFraction (if higher than current) More...
 
virtual void createStencil (const globalIndex &, const bool allowHoleDonors)
 Create stencil starting from the donor containing the acceptor allowHoleDonors : allow donors of type HOLE (otherwise are filtered out) More...
 
void holeExtrapolationStencil (const globalIndex &globalCells)
 Make holes next to live ones type SPECIAL and include in interpolation stencil. More...
 
- Protected Member Functions inherited from cellCellStencil
void suppressMotionFields ()
 Helper: populate nonInterpolatedFields_ with motion solver. More...
 

Static Protected Member Functions

static label index (const labelVector &nDivs, const labelVector &)
 Convert ijk indices into single index. More...
 
static labelVector index3 (const labelVector &nDivs, const label)
 Convert single index into ijk. More...
 
static labelVector index3 (const boundBox &bb, const labelVector &nDivs, const point &pt)
 Convert coordinate into ijk. More...
 
static point position (const boundBox &bb, const labelVector &nDivs, const label boxI)
 Convert index back into coordinate. More...
 
static void fill (PackedList< 2 > &elems, const boundBox &bb, const labelVector &nDivs, const boundBox &subBb, const unsigned int val)
 Fill all elements overlapping subBb with value val. More...
 
static bool overlaps (const boundBox &bb, const labelVector &nDivs, const PackedList< 2 > &voxels, const treeBoundBox &subBb, const unsigned int val)
 Is any voxel inside subBb set to val. More...
 
static void markBoundaries (const fvMesh &mesh, const vector &smallVec, const boundBox &bb, const labelVector &nDivs, PackedList< 2 > &patchTypes, const labelList &cellMap, labelList &patchCellTypes)
 Mark voxels of patchTypes with type of patch face. More...
 
- Static Protected Member Functions inherited from cellCellStencil
template<class Type >
static tmp< volScalarFieldcreateField (const fvMesh &mesh, const word &name, const UList< Type > &)
 Helper: create volScalarField for postprocessing. More...
 
static word baseName (const word &name)
 Helper: strip off trailing _0. More...
 

Protected Attributes

const dictionary dict_
 Dictionary of motion control parameters. More...
 
const bool allowHoleDonors_
 Allow holes as donors. More...
 
const bool allowInterpolatedDonors_
 Allow interpolared as donors. More...
 
vector smallVec_
 Small perturbation vector for geometric tests. More...
 
labelList cellTypes_
 Per cell the cell type. More...
 
labelList interpolationCells_
 Indices of interpolated cells. More...
 
autoPtr< mapDistributecellInterpolationMap_
 Fetch interpolated cells. More...
 
labelListList cellStencil_
 Interpolation stencil. More...
 
scalarListList cellInterpolationWeights_
 Interpolation weights. More...
 
volScalarField cellInterpolationWeight_
 Amount of interpolation. More...
 
- Protected Attributes inherited from cellCellStencil
const fvMeshmesh_
 Reference to the mesh. More...
 
wordHashSet nonInterpolatedFields_
 Set of fields that should not be interpolated. More...
 
const dictionary dict_
 Dictionary of motion control parameters. More...
 

Additional Inherited Members

- Public Types inherited from cellCellStencil
enum  patchCellType { OTHER = 0, PATCH = 1, OVERSET = 2 }
 
enum  cellType {
  CALCULATED = 0, INTERPOLATED = 1, HOLE = 2, SPECIAL = 3,
  POROUS = 4
}
 
- Static Public Member Functions inherited from cellCellStencil
static autoPtr< cellCellStencilNew (const fvMesh &, const dictionary &dict, const bool update=true)
 New function which constructs and returns pointer to a. More...
 
static const labelIOListzoneID (const fvMesh &)
 Helper: get reference to registered zoneID. Loads volScalarField. More...
 
static labelList count (const label size, const labelUList &lst)
 Count occurrences (in parallel) More...
 
static void globalCellCells (const globalIndex &gi, const polyMesh &mesh, const boolList &isValidDonor, const labelList &selectedCells, labelListList &cellCells, pointListList &cellCellCentres)
 Helper: create cell-cell addressing in global numbering. More...
 
template<class T >
static void interpolate (Field< T > &psi, const fvMesh &mesh, const cellCellStencil &overlap, const List< scalarList > &wghts)
 Interpolation of acceptor cells from donor cells. More...
 
template<class GeoField , class SuppressBC >
static void correctBoundaryConditions (GeoField &psi)
 Version of correctBoundaryConditions that excludes 'overset' bcs. More...
 
- Static Protected Attributes inherited from cellCellStencil
static const Enum< cellTypecellTypeNames_
 Mode type names. More...
 

Detailed Description

Inverse-distance-weighted interpolation stencil.

hole finding:

  • mark boundary faces on helper (voxel) mesh
  • mark any cell overlaying these voxels
  • use flood filling to find any unreachable cell Alternative is to use an octree of the boundary faces and determine directly for all cells whether we are outside. Might be slow though.
Source files

Definition at line 63 of file inverseDistanceCellCellStencil.H.

Constructor & Destructor Documentation

◆ inverseDistance()

inverseDistance ( const fvMesh mesh,
const dictionary dict,
const bool  doUpdate 
)

Construct from fvMesh.

Definition at line 1703 of file inverseDistanceCellCellStencil.C.

References Foam::ensightOutput::debug, Foam::endl(), forAll, io(), Foam::Pout, and update().

Here is the call graph for this function:

◆ ~inverseDistance()

~inverseDistance ( )
virtual

Destructor.

Definition at line 1779 of file inverseDistanceCellCellStencil.C.

Member Function Documentation

◆ index()

Foam::label index ( const labelVector nDivs,
const labelVector ijk 
)
staticprotected

Convert ijk indices into single index.

Definition at line 56 of file inverseDistanceCellCellStencil.C.

◆ index3() [1/2]

Foam::labelVector index3 ( const labelVector nDivs,
const label  boxI 
)
staticprotected

Convert single index into ijk.

Definition at line 66 of file inverseDistanceCellCellStencil.C.

References k.

◆ index3() [2/2]

Foam::labelVector index3 ( const boundBox bb,
const labelVector nDivs,
const point pt 
)
staticprotected

Convert coordinate into ijk.

Definition at line 81 of file inverseDistanceCellCellStencil.C.

References boundBox::min(), and boundBox::span().

Here is the call graph for this function:

◆ position()

Foam::point position ( const boundBox bb,
const labelVector nDivs,
const label  boxI 
)
staticprotected

Convert index back into coordinate.

Definition at line 100 of file inverseDistanceCellCellStencil.C.

References boundBox::min(), and boundBox::span().

Here is the call graph for this function:

◆ fill()

void fill ( PackedList< 2 > &  elems,
const boundBox bb,
const labelVector nDivs,
const boundBox subBb,
const unsigned int  val 
)
staticprotected

Fill all elements overlapping subBb with value val.

Definition at line 117 of file inverseDistanceCellCellStencil.C.

References k, boundBox::max(), Foam::max(), boundBox::min(), and Foam::min().

Here is the call graph for this function:

◆ overlaps()

bool overlaps ( const boundBox bb,
const labelVector nDivs,
const PackedList< 2 > &  voxels,
const treeBoundBox subBb,
const unsigned int  val 
)
staticprotected

Is any voxel inside subBb set to val.

Definition at line 231 of file inverseDistanceCellCellStencil.C.

References k, boundBox::max(), Foam::max(), boundBox::min(), and Foam::min().

Here is the call graph for this function:

◆ markBoundaries()

void markBoundaries ( const fvMesh mesh,
const vector smallVec,
const boundBox bb,
const labelVector nDivs,
PackedList< 2 > &  patchTypes,
const labelList cellMap,
labelList patchCellTypes 
)
staticprotected

Mark voxels of patchTypes with type of patch face.

Definition at line 155 of file inverseDistanceCellCellStencil.C.

References fvPatch::faceCells(), forAll, boundBox::grow(), mesh, boundBox::overlaps(), fvPatch::patch(), patchTypes(), pbm, and pp().

Here is the call graph for this function:

◆ markPatchesAsHoles()

void markPatchesAsHoles ( PstreamBuffers pBufs,
const PtrList< fvMeshSubset > &  meshParts,
const List< treeBoundBoxList > &  patchBb,
const List< labelVector > &  patchDivisions,
const PtrList< PackedList< 2 >> &  patchParts,
const label  srcI,
const label  tgtI,
labelList allCellTypes 
) const
protected

Mark all cells overlapping (a voxel covered by) a src patch.

with type HOLE

Definition at line 275 of file inverseDistanceCellCellStencil.C.

References PstreamBuffers::clear(), Foam::exit(), Foam::FatalError, FatalErrorInFunction, PstreamBuffers::finishedSends(), forAll, boundBox::grow(), os(), and treeBoundBox::overlaps().

Here is the call graph for this function:

◆ betterDonor()

bool betterDonor ( const label  destMesh,
const label  currentDonorMesh,
const label  newDonorMesh 
) const
protected

If multiple donors meshes: decide which is best.

Definition at line 399 of file inverseDistanceCellCellStencil.C.

References Foam::mag().

Here is the call graph for this function:

◆ markDonors()

void markDonors ( const globalIndex globalCells,
PstreamBuffers pBufs,
const PtrList< fvMeshSubset > &  meshParts,
const List< treeBoundBoxList > &  meshBb,
const labelList allCellTypes,
const label  srcI,
const label  tgtI,
labelListList allStencil,
labelList allDonor 
) const
protected

◆ uncompactedRegionSplit()

void uncompactedRegionSplit ( const fvMesh mesh,
const globalIndex globalFaces,
const label  nZones,
const labelList zoneID,
const labelList cellTypes,
const boolList isBlockedFace,
labelList cellRegion 
) const
protected

Replacement of regionSplit.

◆ compactedRegionSplit()

autoPtr<globalIndex> compactedRegionSplit ( const fvMesh mesh,
const globalIndex globalFaces,
labelList cellRegion 
) const
protected

◆ findHoles()

void findHoles ( const globalIndex globalCells,
const fvMesh mesh,
const labelList zoneID,
const labelListList stencil,
labelList cellTypes 
) const
protected

Do flood filling to detect unreachable (from patches) sections.

of mesh

Definition at line 894 of file inverseDistanceCellCellStencil.C.

References cellTypes, DebugInfo, mapDistribute::distribute(), Foam::endl(), fvPatch::faceCells(), Foam::findIndices(), forAll, FUNCTION_NAME, mesh, regionSplit::nRegions(), pbm, Foam::reduce(), UList< T >::size(), and Foam::Zero.

Here is the call graph for this function:

◆ seedCell()

void seedCell ( const label  cellI,
const scalar  wantedFraction,
bitSet isFront,
scalarField fraction 
) const
protected

Seed faces of cell with wantedFraction (if higher than current)

◆ createStencil()

void createStencil ( const globalIndex globalCells,
const bool  allowHoleDonors 
)
protectedvirtual

Create stencil starting from the donor containing the acceptor allowHoleDonors : allow donors of type HOLE (otherwise are filtered out)

Note: empty slots can happen for interpolated cells with

bad/insufficient donors. These are handled later on.

Definition at line 1503 of file inverseDistanceCellCellStencil.C.

References Foam::abort(), Foam::FatalError, FatalErrorInFunction, forAll, nSamples(), Foam::returnReduceOr(), samples(), UList< T >::size(), and List< T >::transfer().

Here is the call graph for this function:

◆ holeExtrapolationStencil()

void holeExtrapolationStencil ( const globalIndex globalCells)
protected

Make holes next to live ones type SPECIAL and include in interpolation stencil.

Definition at line 1211 of file inverseDistanceCellCellStencil.C.

References List< T >::append(), DynamicList< T, SizeMin >::append(), Foam::ensightOutput::debug, Foam::endl(), forAll, Foam::identity(), globalIndex::inplaceToGlobal(), Foam::labelMin, Foam::nl, Foam::Pout, Foam::reduce(), and globalIndex::toGlobal().

Here is the call graph for this function:

◆ TypeName()

TypeName ( "inverseDistance"  )

Runtime type information.

◆ update()

◆ cellTypes()

virtual const labelUList& cellTypes ( ) const
inlinevirtual

Return the cell type list.

Implements cellCellStencil.

Definition at line 362 of file inverseDistanceCellCellStencil.H.

References inverseDistance::cellTypes_.

◆ interpolationCells()

virtual const labelUList& interpolationCells ( ) const
inlinevirtual

Indices of interpolated cells.

Implements cellCellStencil.

Definition at line 370 of file inverseDistanceCellCellStencil.H.

References inverseDistance::interpolationCells_.

◆ cellInterpolationMap()

virtual const mapDistribute& cellInterpolationMap ( ) const
inlinevirtual

Return a communication schedule.

Implements cellCellStencil.

Definition at line 378 of file inverseDistanceCellCellStencil.H.

References inverseDistance::cellInterpolationMap_, and inverseDistance::update().

Here is the call graph for this function:

◆ cellStencil()

virtual const labelListList& cellStencil ( ) const
inlinevirtual

Per interpolated cell the neighbour cells (in terms of slots as.

constructed by above cellInterpolationMap) to interpolate

Implements cellCellStencil.

Definition at line 392 of file inverseDistanceCellCellStencil.H.

References inverseDistance::cellStencil_.

◆ cellInterpolationWeights()

virtual const scalarListList& cellInterpolationWeights ( ) const
inlinevirtual

Weights for cellStencil.

Implements cellCellStencil.

Definition at line 400 of file inverseDistanceCellCellStencil.H.

References inverseDistance::cellInterpolationWeights_.

◆ cellInterpolationWeight()

virtual const scalarList& cellInterpolationWeight ( ) const
inlinevirtual

Per interpolated cell the interpolation factor. (0 = use.

calculated, 1 = use interpolated)

Implements cellCellStencil.

Definition at line 410 of file inverseDistanceCellCellStencil.H.

References inverseDistance::cellInterpolationWeight_.

◆ stencilWeights()

void stencilWeights ( const point sample,
const pointList donorCcs,
scalarList weights 
) const
virtual

Calculate inverse distance weights for a single acceptor.

Implements cellCellStencil.

Reimplemented in leastSquares.

Definition at line 1176 of file inverseDistanceCellCellStencil.C.

References forAll, Foam::mag(), List< T >::setSize(), UList< T >::size(), and Foam::sum().

Here is the call graph for this function:

Member Data Documentation

◆ dict_

const dictionary dict_
protected

Dictionary of motion control parameters.

Definition at line 87 of file inverseDistanceCellCellStencil.H.

◆ allowHoleDonors_

const bool allowHoleDonors_
protected

Allow holes as donors.

Definition at line 92 of file inverseDistanceCellCellStencil.H.

◆ allowInterpolatedDonors_

const bool allowInterpolatedDonors_
protected

Allow interpolared as donors.

Definition at line 97 of file inverseDistanceCellCellStencil.H.

◆ smallVec_

vector smallVec_
protected

Small perturbation vector for geometric tests.

Definition at line 102 of file inverseDistanceCellCellStencil.H.

◆ cellTypes_

labelList cellTypes_
protected

Per cell the cell type.

Definition at line 107 of file inverseDistanceCellCellStencil.H.

Referenced by inverseDistance::cellTypes().

◆ interpolationCells_

labelList interpolationCells_
protected

Indices of interpolated cells.

Definition at line 112 of file inverseDistanceCellCellStencil.H.

Referenced by inverseDistance::interpolationCells().

◆ cellInterpolationMap_

autoPtr<mapDistribute> cellInterpolationMap_
protected

Fetch interpolated cells.

Definition at line 117 of file inverseDistanceCellCellStencil.H.

Referenced by inverseDistance::cellInterpolationMap().

◆ cellStencil_

labelListList cellStencil_
protected

Interpolation stencil.

Definition at line 122 of file inverseDistanceCellCellStencil.H.

Referenced by inverseDistance::cellStencil().

◆ cellInterpolationWeights_

scalarListList cellInterpolationWeights_
protected

Interpolation weights.

Definition at line 127 of file inverseDistanceCellCellStencil.H.

Referenced by inverseDistance::cellInterpolationWeights().

◆ cellInterpolationWeight_

volScalarField cellInterpolationWeight_
protected

Amount of interpolation.

Definition at line 132 of file inverseDistanceCellCellStencil.H.

Referenced by inverseDistance::cellInterpolationWeight().


The documentation for this class was generated from the following files: