37 template<
class CloudType>
43 const fvMesh&
mesh = this->owner().mesh();
45 const label nCells = cellZoneCells.size();
46 Random& rnd = this->owner().rndGen();
48 DynamicList<vector> positions(nCells);
49 DynamicList<label> injectorCells(nCells);
50 DynamicList<label> injectorTetFaces(nCells);
51 DynamicList<label> injectorTetPts(nCells);
53 scalar newParticlesTotal = 0.0;
54 label addParticlesTotal = 0;
58 const label celli = cellZoneCells[i];
61 const scalar newParticles = V[celli]*numberDensity_;
62 newParticlesTotal += newParticles;
63 label addParticles = floor(newParticles);
64 addParticlesTotal += addParticles;
66 const scalar
diff = newParticlesTotal - addParticlesTotal;
69 label corr = floor(
diff);
71 addParticlesTotal += corr;
75 const List<tetIndices> cellTetIs =
76 polyMeshTetDecomposition::cellTetIndices(
mesh, celli);
80 for (label tetI = 1; tetI < cellTetIs.size() - 1; tetI++)
83 cTetVFrac[tetI-1] + cellTetIs[tetI].tet(
mesh).mag()/V[celli];
85 cTetVFrac.
last() = 1.0;
88 for (label pI = 0; pI < addParticles; pI++)
90 const scalar volFrac = rnd.sample01<scalar>();
94 if (cTetVFrac[vfI] > volFrac)
100 positions.append(cellTetIs[tetI].tet(
mesh).randomPoint(rnd));
102 injectorCells.append(celli);
103 injectorTetFaces.append(cellTetIs[tetI].face());
104 injectorTetPts.append(cellTetIs[tetI].tetPt());
109 const label myProci = UPstream::myProcNo();
110 globalIndex globalPositions(positions.size());
112 List<vector> allPositions(globalPositions.totalSize(),
point::max);
113 List<label> allInjectorCells(globalPositions.totalSize(), -1);
114 List<label> allInjectorTetFaces(globalPositions.totalSize(), -1);
115 List<label> allInjectorTetPts(globalPositions.totalSize(), -1);
121 globalPositions.range(myProci)
124 Pstream::listCombineReduce(allPositions, minEqOp<point>());
130 globalPositions.range(myProci)
135 globalPositions.range(myProci)
136 ) = injectorTetFaces;
140 globalPositions.range(myProci)
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")),
173 numberDensity_(this->coeffDict().getScalar(
"numberDensity")),
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())
213 template<
class CloudType>
220 template<
class CloudType>
225 const label zoneI =
mesh.cellZones().findZoneID(cellZoneName_);
230 <<
"Unknown cell zone name: " << cellZoneName_
231 <<
". Valid cell zones are: " <<
mesh.cellZones().names()
236 const label nCells = cellZoneCells.size();
237 const scalar nCellsTotal =
returnReduce(nCells, sumOp<label>());
239 const scalar VCellsTotal =
returnReduce(VCells, sumOp<scalar>());
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();
293 template<
class CloudType>
301 if ((0.0 >= time0) && (0.0 < time1))
303 return this->volumeTotal_;
310 template<
class CloudType>
314 const label nParcels,
322 position = positions_[parcelI];
323 cellOwner = injectorCells_[parcelI];
324 tetFacei = injectorTetFaces_[parcelI];
325 tetPti = injectorTetPts_[parcelI];
329 template<
class CloudType>
342 parcel.d() = diameters_[parcelI];
346 template<
class CloudType>
353 template<
class CloudType>
List< scalar > scalarList
List of scalar.
scalar diff(const triad &A, const triad &B)
Return a quantity of the difference between two triads.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &f1)
virtual ~CellZoneInjection()
Destructor.
errorManipArg< error, int > exit(error &err, const int errNo=1)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Templated injection model class.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
constexpr char nl
The newline '\n' character (0x0a)
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.
Lookup type of boundary radiation properties.
#define forAll(list, i)
Loop across all elements in list.
virtual bool fullyDescribed() const
Flag to identify whether model fully describes the parcel.
Injection positions specified by a particle number density within a cell set.
A class for handling words, derived from Foam::string.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
constexpr scalar pi(M_PI)
ParcelType parcelType
Type of parcel the cloud was instantiated for.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
int debug
Static debugging option.
label parcelsToInject(const scalar time0, const scalar time1)
Number of parcels to introduce relative to SOI.
T & last()
Access last element of the list, position [size()-1].
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.
A library of runtime-selectable doubly-truncated probability distribution models. Returns random samp...
return returnReduce(nRefine-oldNRefine, sumOp< label >())
messageStream Info
Information stream (stdout output on master, null elsewhere)
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.
List< label > labelList
A List of labels.
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.
static constexpr const zero Zero
Global zero (0)