35 template<
class CloudType>
46 particleFluxAccumulators_()
56 if (isType<polyPatch>(patch))
66 this->
coeffDict().subDict(
"numberDensities")
72 particleFluxAccumulators_.setSize(patches_.
size());
92 numberDensitiesDict.
lookup<scalar>(molecules[i]);
96 if (moleculeTypeIds_[i] == -1)
99 <<
"typeId " << molecules[i] <<
"not defined in cloud." <<
nl
104 numberDensities_ /=
cloud.nParticle();
111 template<
class CloudType>
118 template<
class CloudType>
133 pFA[facei].
setSize(patch.size(), 0);
139 template<
class CloudType>
152 label particlesInserted = 0;
156 cloud.boundaryT().boundaryField()
161 cloud.boundaryU().boundaryField()
179 label typeId = moleculeTypeIds_[i];
181 scalar mass =
cloud.constProps(typeId).mass();
186 <<
"Zero boundary temperature detected, check boundaryT "
187 <<
"condition." <<
nl
193 cloud.maxwellianMostProbableSpeed
217 exp(-
sqr(sCosTheta)) + sqrtPi*sCosTheta*(1 +
erf(sCosTheta))
227 const face&
f = patch[pFI];
247 scalar previousCumulativeSum = 0.0;
255 + previousCumulativeSum;
257 previousCumulativeSum = cTriAFracs[triI];
262 cTriAFracs.
last() = 1.0;
279 scalar faceTemperature = boundaryT[
patchi][pFI];
285 scalar& faceAccumulator = pFA[i][pFI];
297 faceAccumulator -= nI;
299 label typeId = moleculeTypeIds_[i];
301 scalar mass =
cloud.constProps(typeId).mass();
303 for (
label i = 0; i < nI; i++)
311 label selectedTriI = -1;
317 if (cTriAFracs[triI] >= triSelection)
325 const tetIndices& faceTetIs = faceTets[selectedTriI];
331 scalar mostProbableSpeed
333 cloud.maxwellianMostProbableSpeed
340 scalar sCosTheta = (faceVelocity &
n)/mostProbableSpeed;
343 scalar uNormProbCoeffA =
344 sCosTheta +
sqrt(
sqr(sCosTheta) + 2.0);
346 scalar uNormProbCoeffB =
350 + sCosTheta*(sCosTheta -
sqrt(
sqr(sCosTheta) + 2.0))
354 scalar randomScaling = 3.0;
358 randomScaling =
mag(sCosTheta) + 1;
366 scalar uNormalThermal;
374 uNormal = uNormalThermal + sCosTheta;
382 P = 2.0*uNormal/uNormProbCoeffA
383 *
exp(uNormProbCoeffB -
sqr(uNormalThermal));
394 + (t1 & faceVelocity)*t1
395 + (t2 & faceVelocity)*t2
396 + mostProbableSpeed*uNormal*
n;
398 scalar Ei =
cloud.equipartitionInternalEnergy
401 cloud.constProps(typeId).internalDegreesOfFreedom()
404 cloud.addNewParcel(
p, celli,
U, Ei, typeId);
414 Info<<
" Particles inserted = "
415 << particlesInserted <<
endl;
#define forAll(list, i)
Loop across all elements in list.
Templated base class for dsmc cloud.
virtual void topoChange()
Update following mesh change.
FreeStream(const dictionary &dict, CloudType &cloud)
Construct from dictionary.
virtual void inflow()
Introduce particles.
virtual ~FreeStream()
Destructor.
Generic GeometricBoundaryField class.
Templated inflow boundary model class.
const dictionary & coeffDict() const
Return the coefficients dictionary.
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
void size(const label)
Override size to be inconsistent with allocated storage.
void setSize(const label)
Reset size of List.
scalar scalar01()
Advance the state and return a scalar sample from a uniform.
scalar scalarNormal()
Advance the state and return a scalar sample from a normal.
scalar deltaTValue() const
Return time step value.
T & last()
Return the last element of the list.
A cloud is a collection of lagrangian particles.
A list of keyword definitions, which are a keyword followed by any number of values (e....
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
wordList toc() const
Return the table of contents.
A face is a list of labels corresponding to mesh vertices.
const Time & time() const
Return time.
static List< tetIndices > faceTetIndices(const polyMesh &mesh, label fI, label cI)
Return the tet decomposition of the given face, with.
Mesh consisting of general polyhedral cells.
virtual const labelList & faceOwner() const
Return face owner.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
virtual const pointField & points() const
Return raw points.
A patch is a list of labels that address the faces in the global face list.
const vectorField::subField faceAreas() const
Return face areas.
label start() const
Return start label of this patch in the polyMesh face list.
const vectorField::subField faceCentres() const
Return face centres.
const scalarField::subField magFaceAreas() const
Return face area magnitudes.
Storage and named access for the indices of a tet which is part of the decomposition of a cell.
triPointRef faceTri(const polyMesh &mesh) const
Return the geometry corresponding to the tri on the face for.
scalar mag() const
Return scalar magnitude.
Point randomPoint(Random &rndGen) const
Return a random point on the triangle from a uniform.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const fvPatchList & patches
dimensionedScalar exp(const dimensionedScalar &ds)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
errorManip< error > abort(error &err)
dimensionedScalar erf(const dimensionedScalar &ds)
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.