Go to the documentation of this file.
46 void Foam::volPointInterpolation::calcBoundaryAddressing()
50 Pout<<
"volPointInterpolation::calcBoundaryAddressing() : "
51 <<
"constructing boundary addressing"
62 mesh().nBoundaryFaces(),
63 mesh().nInternalFaces()
71 boundaryIsPatchFace_ =
false;
74 isPatchPoint_ =
false;
84 const polyPatch& pp = pbm[patchi];
88 !isA<emptyPolyPatch>(pp)
89 && !magSf.boundaryField()[patchi].coupled()
96 boundaryIsPatchFace_[bFacei] =
true;
102 isPatchPoint_[
f[fp]] =
true;
121 forAll(isPatchPoint_, pointi)
123 if (isPatchPoint_[pointi] != oldData[pointi])
125 Pout<<
"volPointInterpolation::calcBoundaryAddressing():"
126 <<
" added dangling mesh point:" << pointi
132 label nPatchFace = 0;
133 forAll(boundaryIsPatchFace_, i)
135 if (boundaryIsPatchFace_[i])
140 label nPatchPoint = 0;
143 if (isPatchPoint_[i])
150 <<
" of which on proper patch:" << nPatchFace <<
nl
152 <<
" of which on proper patch:" << nPatchPoint <<
endl;
157 void Foam::volPointInterpolation::makeInternalWeights(
scalarField& sumWeights)
161 Pout<<
"volPointInterpolation::makeInternalWeights() : "
162 <<
"constructing weighting factors for internal and non-coupled"
163 <<
" points." <<
endl;
171 pointWeights_.clear();
172 pointWeights_.setSize(
points.size());
178 if (!isPatchPoint_[pointi])
180 const labelList& pcp = pointCells[pointi];
188 1.0/
mag(
points[pointi] - cellCentres[pcp[pointCelli]]);
190 sumWeights[pointi] += pw[pointCelli];
197 void Foam::volPointInterpolation::makeBoundaryWeights(
scalarField& sumWeights)
201 Pout<<
"volPointInterpolation::makeBoundaryWeights() : "
202 <<
"constructing weighting factors for boundary points." <<
endl;
210 boundaryPointWeights_.clear();
211 boundaryPointWeights_.setSize(
boundary.meshPoints().size());
217 if (isPatchPoint_[pointi])
224 sumWeights[pointi] = 0.0;
228 if (boundaryIsPatchFace_[
pFaces[i]])
232 pw[i] = 1.0/
mag(
points[pointi] - faceCentres[facei]);
233 sumWeights[pointi] += pw[i];
245 void Foam::volPointInterpolation::makeWeights()
249 Pout<<
"volPointInterpolation::makeWeights() : "
250 <<
"constructing weighting factors"
255 calcBoundaryAddressing();
263 "volPointSumWeights",
273 makeInternalWeights(sumWeights);
277 makeBoundaryWeights(sumWeights);
304 addSeparated(sumWeights);
310 pushUntransformedData(sumWeights);
314 forAll(pointWeights_, pointi)
320 pw[i] /= sumWeights[pointi];
335 pw[i] /= sumWeights[pointi];
342 Pout<<
"volPointInterpolation::makeWeights() : "
343 <<
"finished constructing weighting factors"
351 Foam::volPointInterpolation::volPointInterpolation(
const fvMesh& vm)
387 interpolateInternalField(vf, pf);
390 interpolateBoundaryField(vf, pf);
int debug
Static debugging option.
List< label > labelList
A List of labels.
vectorField pointField
pointField is a vectorField.
virtual const pointField & points() const
Return raw points.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
List< scalar > scalarList
A List of scalars.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
static constexpr const zero Zero
Global zero.
label nInternalFaces() const
Number of internal faces.
Template functions to aid in the implementation of demand driven data.
const fileName & instance() const
static const pointMesh & New(const polyMesh &mesh, Args &&... args)
Get existing or create a new MeshObject.
List< bool > boolList
A List of bools.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
void constrainDisplacement(pointVectorField &displacement, const bool overrideValue=false) const
Apply boundary conditions (single-patch constraints),.
prefixOSstream Pout
An Ostream wrapper for parallel output to std::cout.
void updateMesh(const mapPolyMesh &)
Update mesh topology using the morph engine.
#define forAll(list, i)
Loop across all elements in list.
Field< vector > vectorField
Specialisation of Field<T> for vector.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
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]=cellShape(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]
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
const fvMesh & mesh() const
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Mesh data needed to do the Finite Volume discretisation.
~volPointInterpolation()
Destructor.
GeometricField< scalar, pointPatchField, pointMesh > pointScalarField
Application of (multi-)patch point contraints.
List< labelList > labelListList
A List of labelList.
const vectorField & cellCentres() const
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
static void syncUntransformedData(const polyMesh &mesh, List< Type > &pointData, const CombineOp &cop)
Helper: sync data on collocated points only.
bool movePoints()
Correct weighting factors for moving mesh.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const vectorField & faceCentres() const
const labelListList & pointCells() const
void interpolateDisplacement(const volVectorField &, pointVectorField &) const
Interpolate from volField to pointField.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Templated abstract base-class for optional mesh objects used to automate their allocation to the mesh...
PrimitivePatch< face, SubList, const pointField & > primitivePatch
Addressing for a faceList slice.
void setSize(const label newSize)
Alias for resize(const label)
defineTypeNameAndDebug(combustionModel, 0)
Interpolate from cell centres to points (vertices) using inverse distance weighting.