leastSquares Class Reference

Least-squares-weighted interpolation stencil. More...

Inheritance diagram for leastSquares:
Collaboration diagram for leastSquares:

Public Member Functions

 TypeName ("leastSquares")
 Runtime type information. More...
 
 leastSquares (const fvMesh &, const dictionary &, const bool)
 Construct from fvMesh. More...
 
virtual ~leastSquares ()
 Destructor. More...
 
virtual void stencilWeights (const point &sample, const pointList &donorCcs, scalarList &weights) const
 Calculate lsq weights for single acceptor. More...
 
- Public Member Functions inherited from inverseDistance
 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...
 
- 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)
 

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...
 
- Protected Member Functions inherited from inverseDistance
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 inherited from inverseDistance
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 inherited from inverseDistance
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...
 
- Static Protected Attributes inherited from cellCellStencil
static const Enum< cellTypecellTypeNames_
 Mode type names. More...
 

Detailed Description

Least-squares-weighted interpolation stencil.

Base machinery is similar to inverse distance interpolation stencil but weights minimize error in LSQ sense recovering exact solution for linear solution problems. Gradient and values are found simultaneously.

Source files

Definition at line 54 of file leastSquaresCellCellStencil.H.

Constructor & Destructor Documentation

◆ leastSquares()

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

Construct from fvMesh.

Definition at line 151 of file leastSquaresCellCellStencil.C.

References update().

Here is the call graph for this function:

◆ ~leastSquares()

~leastSquares ( )
virtual

Destructor.

Definition at line 168 of file leastSquaresCellCellStencil.C.

Member Function Documentation

◆ TypeName()

TypeName ( "leastSquares"  )

Runtime type information.

◆ stencilWeights()

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

Calculate lsq weights for single acceptor.

Reimplemented from inverseDistance.

Definition at line 39 of file leastSquaresCellCellStencil.C.

References A, forAll, Foam::mag(), Foam::magSqr(), List< T >::setSize(), UList< T >::size(), dimensioned< Type >::T(), and Foam::Zero.

Here is the call graph for this function:

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