Go to the documentation of this file.
48 bool Foam::volPointInterpolation::hasSeparated(
const pointMesh& pMesh)
50 const pointBoundaryMesh& pbm = pMesh.boundary();
52 bool hasSpecial =
false;
53 for (
const auto& pp : pbm)
55 if (isA<coupledFacePointPatch>(pp) && !isType<processorPointPatch>(pp))
62 reduce(hasSpecial, orOp<bool>());
67 void Foam::volPointInterpolation::calcBoundaryAddressing()
71 Pout<<
"volPointInterpolation::calcBoundaryAddressing() : "
72 <<
"constructing boundary addressing"
83 mesh().nBoundaryFaces(),
84 mesh().nInternalFaces()
91 boundaryIsPatchFace_.setSize(
boundary.size());
92 boundaryIsPatchFace_ =
false;
107 const polyPatch& pp = pbm[patchi];
111 !isA<emptyPolyPatch>(pp)
112 && !magSf.boundaryField()[patchi].coupled()
119 boundaryIsPatchFace_[bFacei] =
true;
125 isPatchPoint[
f[fp]] =
true;
142 isPatchPoint_.assign(isPatchPoint);
146 label nPatchFace = 0;
147 forAll(boundaryIsPatchFace_, i)
149 if (boundaryIsPatchFace_[i])
154 label nPatchPoint = 0;
157 if (isPatchPoint_[i])
164 <<
" of which on proper patch:" << nPatchFace <<
nl
166 <<
" of which on proper patch:" << nPatchPoint <<
endl;
171 void Foam::volPointInterpolation::makeInternalWeights(
scalarField& sumWeights)
175 Pout<<
"volPointInterpolation::makeInternalWeights() : "
176 <<
"constructing weighting factors for internal and non-coupled"
177 <<
" points." <<
endl;
185 pointWeights_.clear();
186 pointWeights_.setSize(
points.size());
192 if (!isPatchPoint_[pointi])
194 const labelList& pcp = pointCells[pointi];
202 1.0/
mag(
points[pointi] - cellCentres[pcp[pointCelli]]);
204 sumWeights[pointi] += pw[pointCelli];
211 void Foam::volPointInterpolation::makeBoundaryWeights(
scalarField& sumWeights)
215 Pout<<
"volPointInterpolation::makeBoundaryWeights() : "
216 <<
"constructing weighting factors for boundary points." <<
endl;
224 boundaryPointWeights_.clear();
225 boundaryPointWeights_.setSize(
boundary.meshPoints().size());
229 label pointi =
boundary.meshPoints()[i];
231 if (isPatchPoint_[pointi])
238 sumWeights[pointi] = 0.0;
242 if (boundaryIsPatchFace_[
pFaces[i]])
246 pw[i] = 1.0/
mag(
points[pointi] - faceCentres[facei]);
247 sumWeights[pointi] += pw[i];
259 void Foam::volPointInterpolation::interpolateOne
261 const tmp<scalarField>& tnormalisation,
267 Pout<<
"volPointInterpolation::interpolateOne("
268 <<
"pointScalarField&) : "
269 <<
"interpolating oneField"
270 <<
" from cells to BOUNDARY points "
276 Field<scalar>& pfi = pf.primitiveFieldRef();
283 const label pointi =
mp[i];
284 if (!isPatchPoint_[pointi])
288 scalar& val = pfi[pointi];
293 val += pw[pointCelli];
304 const label pointi =
mp[i];
306 if (isPatchPoint_[pointi])
309 const scalarList& pWeights = boundaryPointWeights_[i];
311 scalar& val = pfi[pointi];
316 if (boundaryIsPatchFace_[
pFaces[j]])
338 const scalarField& normalisation = tnormalisation();
341 pfi[
mp[i]] *= normalisation[i];
351 void Foam::volPointInterpolation::makeWeights()
355 Pout<<
"volPointInterpolation::makeWeights() : "
356 <<
"constructing weighting factors"
363 calcBoundaryAddressing();
367 tmp<pointScalarField> tsumWeights
373 "volPointSumWeights",
385 makeInternalWeights(sumWeights);
389 makeBoundaryWeights(sumWeights);
406 addSeparated(sumWeights);
413 pushUntransformedData(sumWeights);
417 forAll(pointWeights_, pointi)
423 pw[i] /= sumWeights[pointi];
430 const label pointi =
mp[i];
436 pw[i] /= sumWeights[pointi];
446 Pout<<
"volPointInterpolation::makeWeights() : "
447 <<
"detected separated coupled patches"
448 <<
" - allocating normalisation" <<
endl;
452 interpolateOne(tmp<scalarField>(), sumWeights);
456 normalisationPtr_.ref() = scalar(1.0);
461 normalisationPtr_.clear();
467 Pout<<
"volPointInterpolation::makeWeights() : "
468 <<
"finished constructing weighting factors"
476 Foam::volPointInterpolation::volPointInterpolation(
const fvMesh& vm)
516 interpolateInternalField(vf, pf);
519 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 dimensionedScalar mp
Proton mass.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
static constexpr const zero Zero
Global zero (0)
label nInternalFaces() const
Number of internal faces.
Template functions to aid in the implementation of demand driven data.
const fileName & instance() const
static const pointConstraints & New(const pointMesh &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
OSstream wrapped stdout (std::cout) with parallel prefix.
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.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
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.
virtual const fileName & name() const
Return the name of the stream.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Mesh data needed to do the Finite Volume discretisation.
Mesh representing a set of points created from polyMesh.
~volPointInterpolation()
Destructor.
Application of (multi-)patch point constraints.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
List< labelList > labelListList
A List of labelList.
const vectorField & cellCentres() const
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
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< SubList< face >, const pointField & > primitivePatch
A PrimitivePatch with a SubList addressing for the faces, const reference for the point field.
void setSize(const label newSize)
Alias for resize(const label)
defineTypeNameAndDebug(combustionModel, 0)
Interpolate from cell centres to points (vertices) using inverse distance weighting.
GeometricField< scalar, pointPatchField, pointMesh > pointScalarField