Go to the documentation of this file.
33 #include "surfaceInterpolate.H"
44 namespace functionObjects
50 extractEulerianParticles,
72 const label nFaces = fz.size();
79 <<
": Number of faceZone faces (" << allFaces
80 <<
") is less than the number of requested locations ("
87 <<
" faces : " << allFaces <<
nl
100 if (!nInjectorLocations_)
105 const faceZone& fz = mesh_.faceZones()[zoneID_];
114 const label nFaces = fz.size();
115 label nLocations = nInjectorLocations_;
120 scalar fraction = scalar(nFaces)/scalar(nGlobalFaces);
121 nLocations = ceil(fraction*nInjectorLocations_);
124 Pout<<
"nFaces:" << nFaces
125 <<
", nGlobalFaces:" << nGlobalFaces
126 <<
", fraction:" << fraction
127 <<
", nLocations:" << nLocations
145 label nCoarseFaces = 0;
148 fineToCoarseAddr_ = ppa.restrictTopBottomAddressing();
149 nCoarseFaces =
max(fineToCoarseAddr_) + 1;
155 <<
" coarse faces" <<
endl;
192 patchIDs_.setSize(fz.size(), -1);
193 patchFaceIDs_.setSize(fz.size(), -1);
195 label nBlockedFaces = 0;
198 const label facei = fz[localFacei];
200 if (mesh_.isInternalFace(facei))
202 if (alphaf[facei] > alphaThreshold_)
204 blockedFaces[localFacei] =
true;
209 label patchi = mesh_.boundaryMesh().whichPatch(facei);
210 label patchFacei = -1;
215 if (isA<coupledPolyPatch>(pp))
218 refCast<const coupledPolyPatch>(pp);
225 else if (!isA<emptyPolyPatch>(pp))
230 if (patchFacei == -1)
234 else if (alphafp[patchFacei] > alphaThreshold_)
236 blockedFaces[localFacei] =
true;
239 patchIDs_[localFacei] = patchi;
240 patchFaceIDs_[localFacei] = patchFacei;
256 const label particlei = regionToParticleMap_[regioni];
259 if (
p.faceIHit != -1 && nInjectorLocations_)
262 label coarseFacei = fineToCoarseAddr_[
p.faceIHit];
263 p.faceIHit = globalCoarseFaces_.toGlobal(coarseFacei);
270 if ((pDiameter > minDiameter_) && (pDiameter < maxDiameter_))
275 const point position =
p.VC/(
p.V + ROOTVSMALL);
276 const vector U =
p.VU/(
p.V + ROOTVSMALL);
278 if (nInjectorLocations_)
294 cloud_.addParticle(ip);
296 collectedVolume_ +=
p.V;
299 ++nCollectedParticles_;
304 ++nDiscardedParticles_;
305 discardedVolume_ +=
p.V;
312 const label nNewRegions,
321 labelList oldToNewRegion(particles_.size(), -1);
324 forAll(regionFaceIDs, facei)
326 label newRegioni = regionFaceIDs[facei];
327 label oldRegioni = regions0_[facei];
329 if (newRegioni != -1 && oldRegioni != -1)
333 newToNewRegion[newRegioni] =
334 max(newRegioni, oldToNewRegion[oldRegioni]);
335 oldToNewRegion[oldRegioni] = newRegioni;
344 label nParticle = -1;
347 forAll(newToNewRegion, newRegioni0)
349 label newRegioni = newToNewRegion[newRegioni0];
350 if (newRegions.
insert(newRegioni))
356 newRegionToParticleMap.insert(newRegioni0, nParticle);
364 forAll(oldToNewRegion, oldRegioni)
366 label newRegioni = oldToNewRegion[oldRegioni];
367 if (newRegioni == -1)
371 <<
"Collecting particle from oldRegion:" << oldRegioni
374 collectParticle(time, oldRegioni);
379 label newParticlei = newRegionToParticleMap[newRegioni];
380 label oldParticlei = regionToParticleMap_[oldRegioni];
383 <<
"Combining newRegioni: " << newRegioni
384 <<
"(p:" << newParticlei <<
") and "
385 <<
"oldRegioni: " << oldRegioni
386 <<
"(p:" << oldParticlei <<
")"
389 newParticles[newParticlei] =
392 newParticles[newParticlei],
393 particles_[oldParticlei]
399 particles_.transfer(newParticles);
400 regionToParticleMap_ = newRegionToParticleMap;
404 regions0_ = regionFaceIDs;
421 const scalar deltaT = mesh_.time().deltaTValue();
422 const pointField& faceCentres = mesh_.faceCentres();
424 forAll(regionFaceIDs, localFacei)
426 const label newRegioni = regionFaceIDs[localFacei];
427 if (newRegioni != -1)
429 const label particlei = regionToParticleMap_[newRegioni];
430 const label meshFacei = fz[localFacei];
436 p.faceIHit = localFacei;
443 scalar magPhii =
mag(faceValue(
phi, localFacei, meshFacei));
444 vector Ufi = faceValue(
Uf, localFacei, meshFacei);
445 scalar dV = magPhii*deltaT;
447 p.VC += dV*faceCentres[meshFacei];
465 cloud_(mesh_,
"eulerianParticleCloud"),
471 alphaThreshold_(0.1),
475 nInjectorLocations_(0),
477 globalCoarseFaces_(),
480 regionToParticleMap_(),
481 minDiameter_(ROOTVSMALL),
483 nCollectedParticles_(getProperty<label>(
"nCollectedParticles", 0)),
484 collectedVolume_(getProperty<scalar>(
"collectedVolume", 0)),
485 nDiscardedParticles_(getProperty<label>(
"nDiscardedParticles", 0)),
486 discardedVolume_(getProperty<scalar>(
"discardedVolume", 0))
488 if (mesh_.nSolutionD() != 3)
491 <<
name <<
" function object only applicable to 3-D cases"
510 dict.readEntry(
"faceZone", faceZoneName_);
511 dict.readEntry(
"alpha", alphaName_);
513 dict.readIfPresent(
"alphaThreshold", alphaThreshold_);
514 dict.readIfPresent(
"U", UName_);
515 dict.readIfPresent(
"rho", rhoName_);
516 dict.readIfPresent(
"phi", phiName_);
517 dict.readIfPresent(
"nLocations", nInjectorLocations_);
518 dict.readIfPresent(
"minDiameter", minDiameter_);
519 dict.readIfPresent(
"maxDiameter", maxDiameter_);
523 if (nInjectorLocations_)
546 typeName +
":alphaf",
550 const faceZone& fz = mesh_.faceZones()[zoneID_];
558 boolList blockedFaces(fz.size(),
false);
559 setBlockedFaces(alphaf, fz, blockedFaces);
566 const label nRegionsNew = regionFaceIDs.
nRegions();
574 mesh_.time().value(),
580 accumulateParticleInfo(alphaf, tphi(), regionFaceIDs, fz);
582 Log <<
" Collected particles : " << nCollectedParticles_ <<
nl
583 <<
" Collected volume : " << collectedVolume_ <<
nl
584 <<
" Discarded particles : " << nDiscardedParticles_ <<
nl
585 <<
" Discarded volume : " << discardedVolume_ <<
nl
586 <<
" Particles in progress : " << particles_.size() <<
nl
599 setProperty(
"nCollectedParticles", nCollectedParticles_);
600 setProperty(
"collectedVolume", collectedVolume_);
601 setProperty(
"nDiscardedParticles", nDiscardedParticles_);
602 setProperty(
"discardedVolume", discardedVolume_);
int debug
Static debugging option.
extractEulerianParticles(const word &name, const Time &runTime, const dictionary &dict)
Construct from components.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
label nInjectorLocations_
Number of sample locations to generate.
A class for handling words, derived from Foam::string.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
virtual bool read(const dictionary &)
Read the field min/max data.
A class for managing temporary objects.
Splits a patch into regions based on a mask field. Result is a globally consistent label list of regi...
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
virtual void initialiseBins()
Initialise the particle collection bins.
bool read(const char *buf, int32_t &val)
Same as readInt32.
static bool & parRun()
Is this a parallel run?
autoPtr< surfaceVectorField > Uf
Primitive patch pair agglomerate method.
Ostream & endl(Ostream &os)
Add newline and flush stream.
The coupledPolyPatch is an abstract base class for patches that couple regions of the computational d...
prefixOSstream Pout
An Ostream wrapper for parallel output to std::cout.
virtual bool owner() const =0
Does this side own the patch ?
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
#define forAll(list, i)
Loop across all elements in list.
label nRegions() const
Return the global number of regions.
const faceZoneMesh & faceZones() const
Return face zone mesh.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
virtual bool read(const dictionary &dict)
Read.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
virtual void collectParticle(const scalar time, const label regioni)
Collect particles that have passed through the faceZone.
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
virtual void checkFaceZone()
Check that the faceZone is valid.
A subset of mesh faces organised as a primitive patch.
messageStream Info
Information stream (uses stdout - output is on the master only)
#define DebugInFunction
Report an information message using Foam::Info.
A patch is a list of labels that address the faces in the global face list.
word name(const complex &c)
Return string representation of complex.
virtual void calculateAddressing(const label nRegionsNew, const scalar time, labelList ®ionFaceIDs)
A List with indirect addressing.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
labelList regions0_
Region indices in faceZone faces from last iteration.
Primarily stores particle properties so that it can be injected at a later time. Note that this store...
virtual bool write()
Write.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
virtual tmp< surfaceScalarField > phiU() const
Return the volumetric flux.
virtual bool read(const dictionary &dict)
Read optional controls.
Macros for easy insertion into run-time selection tables.
virtual bool write(const token &tok)=0
Write token to stream or otherwise handle it.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
label findZoneID(const word &zoneName) const
Find zone index given a name, return -1 if not found.
wordList names() const
A list of the zone names.
errorManipArg< error, int > exit(error &err, const int errNo=1)
static bool master(const label communicator=0)
Am I the master process.
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
#define DebugInfo
Report an information message using Foam::Info.
label whichFace(const label l) const
Return label of face in patch from global face label.
constexpr scalar pi(M_PI)
const word & name() const
Return the name of this functionObject.
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
const std::string patch
OpenFOAM patch number as a std::string.
void agglomerate()
Agglomerate patch.
virtual const word & type() const =0
Runtime type information.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
const fvMesh & mesh_
Reference to the fvMesh.
labelList identity(const label len, label start=0)
Create identity map of the given length with (map[i] == i)
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
defineTypeNameAndDebug(ObukhovLength, 0)
static const word null
An empty word.
static const Vector< scalar > zero
Base class for writing single files from the function objects.
word faceZoneName_
Name of faceZone to sample.
virtual void setBlockedFaces(const surfaceScalarField &alphaf, const faceZone &fz, boolList &blockedFaces)
Set the blocked faces, i.e. where alpha > alpha threshold value.
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
dimensionedScalar cbrt(const dimensionedScalar &ds)
void setSize(const label newSize)
Alias for resize(const label)
label zoneID_
Index of the faceZone.
virtual bool execute()
Execute.
virtual void accumulateParticleInfo(const surfaceScalarField &alphaf, const surfaceScalarField &phi, const labelList ®ionFaceIDs, const faceZone &fz)
Process latest region information.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
A list of faces which address into the list of points.