34 template<
class CloudType>
40 const fvMesh& mesh = this->owner().mesh();
42 const label nCells = cellZoneCells.size();
43 cachedRandom& rnd = this->owner().rndGen();
45 DynamicList<vector> positions(nCells);
46 DynamicList<label> injectorCells(nCells);
47 DynamicList<label> injectorTetFaces(nCells);
48 DynamicList<label> injectorTetPts(nCells);
50 scalar newParticlesTotal = 0.0;
51 label addParticlesTotal = 0;
55 const label celli = cellZoneCells[i];
58 const scalar newParticles = V[celli]*numberDensity_;
59 newParticlesTotal += newParticles;
60 label addParticles = floor(newParticles);
61 addParticlesTotal += addParticles;
63 const scalar diff = newParticlesTotal - addParticlesTotal;
66 label corr = floor(diff);
68 addParticlesTotal += corr;
72 const List<tetIndices> cellTetIs =
73 polyMeshTetDecomposition::cellTetIndices(mesh, celli);
77 for (
label tetI = 1; tetI < cellTetIs.size() - 1; tetI++)
80 cTetVFrac[tetI-1] + cellTetIs[tetI].tet(mesh).mag()/V[celli];
82 cTetVFrac.last() = 1.0;
85 for (
label pI = 0; pI < addParticles; pI++)
87 const scalar volFrac = rnd.sample01<scalar>();
91 if (cTetVFrac[vfI] > volFrac)
97 positions.append(cellTetIs[tetI].tet(mesh).randomPoint(rnd));
99 injectorCells.append(celli);
100 injectorTetFaces.append(cellTetIs[tetI].face());
101 injectorTetPts.append(cellTetIs[tetI].tetPt());
106 globalIndex globalPositions(positions.size());
107 List<vector> allPositions(globalPositions.size(),
point::max);
108 List<label> allInjectorCells(globalPositions.size(), -1);
109 List<label> allInjectorTetFaces(globalPositions.size(), -1);
110 List<label> allInjectorTetPts(globalPositions.size(), -1);
116 globalPositions.localSize(Pstream::myProcNo()),
117 globalPositions.offset(Pstream::myProcNo())
120 Pstream::listCombineGather(allPositions, minEqOp<point>());
121 Pstream::listCombineScatter(allPositions);
127 globalPositions.localSize(Pstream::myProcNo()),
128 globalPositions.offset(Pstream::myProcNo())
133 globalPositions.localSize(Pstream::myProcNo()),
134 globalPositions.offset(Pstream::myProcNo())
135 ) = injectorTetFaces;
139 globalPositions.localSize(Pstream::myProcNo()),
140 globalPositions.offset(Pstream::myProcNo())
144 positions_.transfer(allPositions);
145 injectorCells_.transfer(allInjectorCells);
146 injectorTetFaces_.transfer(allInjectorTetFaces);
147 injectorTetPts_.transfer(allInjectorTetPts);
152 OFstream
points(
"points.obj");
163 template<
class CloudType>
168 const word& modelName
172 cellZoneName_(this->coeffDict().lookup(
"cellZone")),
179 U0_(this->coeffDict().
lookup(
"U0")),
184 this->coeffDict().subDict(
"sizeDistribution"), owner.rndGen()
192 template<
class CloudType>
199 cellZoneName_(im.cellZoneName_),
200 numberDensity_(im.numberDensity_),
201 positions_(im.positions_),
202 injectorCells_(im.injectorCells_),
203 injectorTetFaces_(im.injectorTetFaces_),
204 injectorTetPts_(im.injectorTetPts_),
205 diameters_(im.diameters_),
207 sizeDistribution_(im.sizeDistribution_().clone().ptr())
213 template<
class CloudType>
220 template<
class CloudType>
224 const fvMesh& mesh = this->owner().mesh();
230 <<
"Unknown cell zone name: " << cellZoneName_
236 const label nCells = cellZoneCells.
size();
240 Info<<
" cell zone size = " << nCellsTotal <<
endl;
241 Info<<
" cell zone volume = " << VCellsTotal <<
endl;
243 if ((nCellsTotal == 0) || (VCellsTotal*numberDensity_ < 1))
246 <<
"Number of particles to be added to cellZone " << cellZoneName_
247 <<
" is zero" <<
endl;
251 setPositions(cellZoneCells);
253 Info<<
" number density = " << numberDensity_ <<
nl 254 <<
" number of particles = " << positions_.size() <<
endl;
257 diameters_.setSize(positions_.size());
260 diameters_[i] = sizeDistribution_->sample();
269 template<
class CloudType>
277 template<
class CloudType>
284 if ((0.0 >= time0) && (0.0 < time1))
286 return positions_.size();
295 template<
class CloudType>
303 if ((0.0 >= time0) && (0.0 < time1))
305 return this->volumeTotal_;
314 template<
class CloudType>
318 const label nParcels,
326 position = positions_[parcelI];
327 cellOwner = injectorCells_[parcelI];
328 tetFacei = injectorTetFaces_[parcelI];
329 tetPtI = injectorTetPts_[parcelI];
333 template<
class CloudType>
346 parcel.d() = diameters_[parcelI];
350 template<
class CloudType>
357 template<
class CloudType>
#define forAll(list, i)
Loop across all elements in list.
virtual ~CellZoneInjection()
Destructor.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
errorManipArg< error, int > exit(error &err, const int errNo=1)
A list of keyword definitions, which are a keyword followed by any number of values (e...
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Templated injection model class.
void size(const label)
Override size to be inconsistent with allocated storage.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
label findZoneID(const word &zoneName) const
Find zone index given a name.
stressControl lookup("compactNormalStress") >> compactNormalStress
Injection positions specified by a particle number density within a cell set.
A class for handling words, derived from string.
wordList names() const
Return a list of zone names.
List< scalar > scalarList
A List of scalars.
List< label > labelList
A List of labels.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
volScalarField scalarField(fieldObject, mesh)
label parcelsToInject(const scalar time0, const scalar time1)
Number of parcels to introduce relative to SOI.
CellZoneInjection(const dictionary &dict, CloudType &owner, const word &modelName)
Construct from dictionary.
dimensionedScalar pow3(const dimensionedScalar &ds)
const cellZoneMesh & cellZones() const
Return cell zone mesh.
#define WarningInFunction
Report a warning using Foam::Warning.
virtual void setProperties(const label parcelI, const label nParcels, const scalar time, typename CloudType::parcelType &parcel)
Set the parcel properties.
Mesh data needed to do the Finite Volume discretisation.
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
virtual bool fullyDescribed() const
Flag to identify whether model fully describes the parcel.
virtual void setPositionAndCell(const label parcelI, const label nParcels, const scalar time, vector &position, label &cellOwner, label &tetFacei, label &tetPtI)
Set the injection position and owner cell, tetFace and tetPt.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
scalar timeEnd() const
Return the end-of-injection time.
Templated base class for dsmc cloud.
virtual void updateMesh()
Set injector locations when mesh is updated.
scalar volumeToInject(const scalar time0, const scalar time1)
Volume of parcels to introduce relative to SOI.
virtual bool validInjection(const label parcelI)
Return flag to identify whether or not injection of parcelI is.