39 void Foam::volPointInterpolation::pushUntransformedData
47 const labelList& meshPoints = cpp.meshPoints();
49 const mapDistribute& slavesMap = gmd.globalCoPointSlavesMap();
52 List<Type> elems(slavesMap.constructSize());
55 elems[i] = pointData[meshPoints[i]];
66 elems[slavePoints[j]] = elems[i];
71 slavesMap.reverseDistribute(elems.size(), elems,
false);
76 pointData[meshPoints[i]] = elems[i];
82 void Foam::volPointInterpolation::addSeparated
84 GeometricField<Type, pointPatchField, pointMesh>& pf
89 Pout<<
"volPointInterpolation::addSeparated" <<
endl;
92 auto& pfi = pf.internalFieldRef();
93 auto& pfbf = pf.boundaryFieldRef();
101 refCast<coupledPointPatchField<Type>>
102 (pfbf[patchi]).initSwapAddSeparated
117 refCast<coupledPointPatchField<Type>>
118 (pfbf[patchi]).swapAddSeparated
131 const GeometricField<Type, fvPatchField, volMesh>& vf,
132 GeometricField<Type, pointPatchField, pointMesh>& pf
137 Pout<<
"volPointInterpolation::interpolateInternalField(" 138 <<
"const GeometricField<Type, fvPatchField, volMesh>&, " 139 <<
"GeometricField<Type, pointPatchField, pointMesh>&) : " 140 <<
"interpolating field " << vf.
name()
141 <<
" from cells to points " << pf.
name() <<
endl;
147 forAll(pointCells, pointi)
149 if (!isPatchPoint_[pointi])
152 const labelList& ppc = pointCells[pointi];
158 pf[pointi] += pw[pointCelli]*vf[ppc[pointCelli]];
168 const DimensionedField<Type, volMesh>& vf,
169 DimensionedField<Type, pointMesh>& pf
174 Pout<<
"volPointInterpolation::interpolateDimensionedInternalField(" 175 <<
"const DimensionedField<Type, volMesh>&, " 176 <<
"DimensionedField<Type, pointMesh>&) : " 177 <<
"interpolating field " << vf.
name() <<
" from cells to points " 181 const fvMesh&
mesh = vf.mesh();
193 forAll(pointCells, pointi)
195 const labelList& ppc = pointCells[pointi];
197 pf[pointi] = Type(
Zero);
201 label celli = ppc[pointCelli];
202 scalar pw = 1.0/
mag(
points[pointi] - cellCentres[celli]);
204 pf[pointi] += pw*vf[celli];
216 scalar
s = sumW[pointi];
227 Foam::volPointInterpolation::flatBoundaryField
229 const GeometricField<Type, fvPatchField, volMesh>& vf
235 auto&
values = tboundaryVals.ref();
237 forAll(vf.boundaryField(), patchi)
239 const auto&
pp =
pbm[patchi];
240 const auto& pfld = vf.boundaryField()[patchi];
245 SubList<Type> slice(
values, pfld.size(),
pp.offset());
249 !isA<emptyPolyPatch>(
pp)
257 return tboundaryVals;
284 label pointi =
mp[i];
286 if (isPatchPoint_[pointi])
289 const scalarList& pWeights = boundaryPointWeights_[i];
291 Type& val = pfi[pointi];
296 if (boundaryIsPatchFace_[
pFaces[j]])
298 val += pWeights[j]*boundaryVals[
pFaces[j]];
311 if (normalisationPtr_)
313 const scalarField& normalisation = normalisationPtr_();
316 pfi[
mp[i]] *= normalisation[i];
324 pushUntransformedData(pfi);
333 const bool overrideFixedValue
336 interpolateBoundaryField(vf, pf);
354 Pout<<
"volPointInterpolation::interpolate(" 355 <<
"const GeometricField<Type, fvPatchField, volMesh>&, " 356 <<
"GeometricField<Type, pointPatchField, pointMesh>&) : " 357 <<
"interpolating field " << vf.
name() <<
" from cells to points " 361 interpolateInternalField(vf, pf);
364 interpolateBoundaryField(vf, pf,
false);
383 "volPointInterpolate(" + vf.
name() +
')',
392 interpolateInternalField(vf, tpf.ref());
395 interpolateBoundaryField(vf, tpf.ref(),
true);
431 PointFieldType* pfPtr =
432 db.objectRegistry::template getObjectPtr<PointFieldType>(
name);
434 if (!cache || vf.
mesh().changing())
437 if (pfPtr && pfPtr->ownedByRegistry())
470 PointFieldType& pf = *pfPtr;
528 PointFieldType* pfPtr =
529 db.objectRegistry::template getObjectPtr<PointFieldType>(
name);
531 if (!cache || vf.
mesh().changing())
534 if (pfPtr && pfPtr->ownedByRegistry())
552 interpolateDimensionedInternalField(vf, tpf.ref());
567 PointFieldType& pf = *pfPtr;
576 interpolateDimensionedInternalField(vf, pf);
List< scalar > scalarList
List of scalar.
const polyBoundaryMesh & pbm
tmp< GeometricField< Type, pointPatchField, pointMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &) const
Interpolate volField using inverse distance weighting.
void size(const label n)
Older name for setAddressableSize.
void interpolateBoundaryField(const GeometricField< Type, fvPatchField, volMesh > &vf, GeometricField< Type, pointPatchField, pointMesh > &pf) const
Interpolate boundary field without applying constraints/boundary.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const word & name() const noexcept
Return the object name.
static const pointConstraints & New(const pointMesh &mesh, Args &&... args)
Get existing or create MeshObject registered with typeName.
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
static label nRequests() noexcept
Number of outstanding requests (on the internal list of requests)
Ostream & endl(Ostream &os)
Add newline and flush stream.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
bool store()
Register object with its registry and transfer ownership to the registry.
Generic GeometricField class.
PrimitivePatch< IndirectList< face >, const pointField & > indirectPrimitivePatch
A PrimitivePatch with an IndirectList for the faces, const reference for the point field...
virtual const fileName & name() const override
Get the name of the output serial stream. (eg, the name of the Fstream file name) ...
List< labelList > labelListList
List of labelList.
static void waitRequests()
Wait for all requests to finish.
virtual const pointField & points() const
Return raw points.
Mesh representing a set of points created from polyMesh.
#define forAll(list, i)
Loop across all elements in list.
Application of (multi-)patch point constraints.
static void cachePrintMessage(const char *message, const word &name, const FieldType &fld)
Helper for printing cache message.
List< T > values(const HashTable< T, Key, Hash > &tbl, const bool doSort=false)
List of values from HashTable, optionally sorted.
A list of faces which address into the list of points.
vectorField pointField
pointField is a vectorField.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
const polyMesh & mesh() const noexcept
Return the mesh reference.
Generic templated field type.
const fvMesh & mesh() const noexcept
Reference to the mesh.
const polyBoundaryMesh & boundaryMesh() const noexcept
Return boundary mesh.
A class for handling words, derived from Foam::string.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static void syncUntransformedData(const polyMesh &mesh, List< Type > &pointData, const CombineOp &cop)
Helper: sync data on collocated points only.
const objectRegistry & thisDb() const
Return database. For now is its polyMesh.
const globalMeshData & globalData() const
Return parallel info (demand-driven)
void interpolateInternalField(const GeometricField< Type, fvPatchField, volMesh > &, GeometricField< Type, pointPatchField, pointMesh > &) const
Interpolate internal field from volField to pointField.
const vectorField & cellCentres() const
const labelListList & pointCells() const
const Mesh & mesh() const noexcept
Return mesh.
static std::string name(const std::string &str)
Return basename (part beyond last /), including its extension.
int debug
Static debugging option.
void interpolateDimensionedInternalField(const DimensionedField< Type, volMesh > &vf, DimensionedField< Type, pointMesh > &pf) const
Interpolate dimensioned internal field from cells to points.
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].reset(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 fileName & instance() const noexcept
Read access to instance path component.
label nFaces() const noexcept
The number of boundary faces in the underlying mesh.
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field values.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
"nonBlocking" (immediate) : (MPI_Isend, MPI_Irecv)
Internal & ref(const bool updateAccessTime=true)
Same as internalFieldRef()
A class representing the concept of 0 (zero) that can be used to avoid manipulating objects known to ...
Field< vector > vectorField
Specialisation of Field<T> for vector.
List< label > labelList
A List of labels.
A class for managing temporary objects.
Registry of regIOobjects.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Defines the attributes of an object for which implicit objectRegistry management is supported...
void constrain(GeometricField< Type, pointPatchField, pointMesh > &pf, const bool overrideValue=false) const
Apply boundary conditions (single-patch constraints) and.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
const dimensionedScalar mp
Proton mass.
uindirectPrimitivePatch pp(UIndirectList< face >(mesh.faces(), faceLabels), mesh.points())
const dimensionSet & dimensions() const noexcept
Return dimensions.
static constexpr const zero Zero
Global zero (0)