34 template<
class CloudType>
40 const fvMesh& mesh = this->owner().mesh();
42 const label nCells = cellZoneCells.size();
43 Random& 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 cTetVFrac[0] = cellTetIs[0].tet(mesh).mag();
78 for (
label tetI = 1; tetI < cellTetIs.size(); tetI++)
81 cTetVFrac[tetI-1] + cellTetIs[tetI].tet(mesh).mag();
83 cTetVFrac /= cTetVFrac.last();
86 for (
label pI = 0; pI < addParticles; pI++)
88 const scalar volFrac = rnd.sample01<scalar>();
92 if (cTetVFrac[vfI] > volFrac)
98 positions.append(cellTetIs[tetI].tet(mesh).randomPoint(rnd));
100 injectorCells.append(celli);
101 injectorTetFaces.append(cellTetIs[tetI].face());
102 injectorTetPts.append(cellTetIs[tetI].tetPt());
107 globalIndex globalPositions(positions.size());
108 List<vector> allPositions(globalPositions.size(),
point::max);
109 List<label> allInjectorCells(globalPositions.size(), -1);
110 List<label> allInjectorTetFaces(globalPositions.size(), -1);
111 List<label> allInjectorTetPts(globalPositions.size(), -1);
117 globalPositions.localSize(Pstream::myProcNo()),
118 globalPositions.offset(Pstream::myProcNo())
121 Pstream::listCombineGather(allPositions, minEqOp<point>());
122 Pstream::listCombineScatter(allPositions);
128 globalPositions.localSize(Pstream::myProcNo()),
129 globalPositions.offset(Pstream::myProcNo())
134 globalPositions.localSize(Pstream::myProcNo()),
135 globalPositions.offset(Pstream::myProcNo())
136 ) = injectorTetFaces;
140 globalPositions.localSize(Pstream::myProcNo()),
141 globalPositions.offset(Pstream::myProcNo())
145 positions_.transfer(allPositions);
146 injectorCells_.transfer(allInjectorCells);
147 injectorTetFaces_.transfer(allInjectorTetFaces);
148 injectorTetPts_.transfer(allInjectorTetPts);
153 OFstream
points(
"points.obj");
164 template<
class CloudType>
169 const word& modelName
173 cellZoneName_(this->coeffDict().lookup(
"cellZone")),
174 numberDensity_(this->coeffDict().template lookup<scalar>(
"numberDensity")),
180 U0_(this->coeffDict().
lookup(
"U0")),
185 this->coeffDict().subDict(
"sizeDistribution"), owner.rndGen()
193 template<
class CloudType>
200 cellZoneName_(im.cellZoneName_),
201 numberDensity_(im.numberDensity_),
202 positions_(im.positions_),
203 injectorCells_(im.injectorCells_),
204 injectorTetFaces_(im.injectorTetFaces_),
205 injectorTetPts_(im.injectorTetPts_),
206 diameters_(im.diameters_),
208 sizeDistribution_(im.sizeDistribution_().clone().ptr())
214 template<
class CloudType>
221 template<
class CloudType>
225 const fvMesh& mesh = this->owner().mesh();
231 <<
"Unknown cell zone name: " << cellZoneName_
237 const label nCells = cellZoneCells.
size();
241 Info<<
" cell zone size = " << nCellsTotal <<
endl;
242 Info<<
" cell zone volume = " << VCellsTotal <<
endl;
244 if ((nCellsTotal == 0) || (VCellsTotal*numberDensity_ < 1))
247 <<
"Number of particles to be added to cellZone " << cellZoneName_
248 <<
" is zero" <<
endl;
252 setPositions(cellZoneCells);
254 Info<<
" number density = " << numberDensity_ <<
nl 255 <<
" number of particles = " << positions_.size() <<
endl;
258 diameters_.setSize(positions_.size());
261 diameters_[i] = sizeDistribution_->sample();
270 template<
class CloudType>
278 template<
class CloudType>
285 if ((0.0 >= time0) && (0.0 < time1))
287 return positions_.size();
296 template<
class CloudType>
304 if ((0.0 >= time0) && (0.0 < time1))
306 return this->volumeTotal_;
315 template<
class CloudType>
319 const label nParcels,
327 position = positions_[parcelI];
328 cellOwner = injectorCells_[parcelI];
329 tetFacei = injectorTetFaces_[parcelI];
330 tetPti = injectorTetPts_[parcelI];
334 template<
class CloudType>
347 parcel.d() = diameters_[parcelI];
351 template<
class CloudType>
358 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.
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
stressControl lookup("compactNormalStress") >> compactNormalStress
virtual bool fullyDescribed() const
Flag to identify whether model fully describes the parcel.
label findZoneID(const word &zoneName) const
Find zone index given a name.
Injection positions specified by a particle number density within a cell set.
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
A class for handling words, derived from string.
const cellZoneMesh & cellZones() const
Return cell zone mesh.
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)
#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.
wordList names() const
Return a list of zone names.
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.
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.
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.