37 template<
class CloudType>
48 particleFluxAccumulators_()
58 if (isType<polyPatch>(
patch))
68 this->
coeffDict().subDict(
"numberDensities")
74 particleFluxAccumulators_.setSize(patches_.
size());
87 moleculeTypeIds_.
setSize(molecules.size());
89 numberDensities_.
setSize(molecules.size());
93 numberDensities_[i] = numberDensitiesDict.
get<scalar>(molecules[i]);
95 moleculeTypeIds_[i] =
cloud.typeIdList().
find(molecules[i]);
97 if (moleculeTypeIds_[i] == -1)
100 <<
"typeId " << molecules[i] <<
"not defined in cloud." <<
nl 105 numberDensities_ /=
cloud.nParticle();
112 template<
class CloudType>
119 template<
class CloudType>
127 label patchi = patches_[
p];
140 template<
class CloudType>
145 const polyMesh&
mesh(cloud.mesh());
149 Random&
rndGen = cloud.rndGen();
153 label particlesInserted = 0;
157 cloud.boundaryT().boundaryField()
162 cloud.boundaryU().boundaryField()
168 label patchi = patches_[
p];
176 List<Field<scalar>>& pFA = particleFluxAccumulators_[
p];
180 label typeId = moleculeTypeIds_[i];
182 scalar mass = cloud.constProps(typeId).mass();
184 if (
min(boundaryT[patchi]) < SMALL)
187 <<
"Zero boundary temperature detected, check boundaryT " 188 <<
"condition." <<
nl 194 cloud.maxwellianMostProbableSpeed
208 (boundaryU[patchi] & -
patch.faceAreas()/
mag(
patch.faceAreas()))
215 mag(
patch.faceAreas())*numberDensities_[i]*deltaT
218 exp(-
sqr(sCosTheta)) + sqrtPi*sCosTheta*(1 +
erf(sCosTheta))
228 const face&
f =
patch[pFI];
230 label globalFaceIndex = pFI +
patch.start();
236 scalar fA =
mag(
patch.faceAreas()[pFI]);
246 List<scalar> cTriAFracs(faceTets.size(),
Zero);
248 scalar previousCummulativeSum = 0.0;
252 const tetIndices& faceTetIs = faceTets[triI];
255 faceTetIs.faceTri(
mesh).mag()/fA
256 + previousCummulativeSum;
258 previousCummulativeSum = cTriAFracs[triI];
263 cTriAFracs.last() = 1.0;
277 scalar faceTemperature = boundaryT[patchi][pFI];
279 const vector& faceVelocity = boundaryU[patchi][pFI];
283 scalar& faceAccumulator = pFA[i][pFI];
286 label nI =
max(label(faceAccumulator), 0);
295 faceAccumulator -= nI;
297 label typeId = moleculeTypeIds_[i];
299 scalar mass = cloud.constProps(typeId).mass();
301 for (label i = 0; i < nI; i++)
309 label selectedTriI = -1;
315 if (cTriAFracs[triI] >= triSelection)
323 const tetIndices& faceTetIs = faceTets[selectedTriI];
329 scalar mostProbableSpeed
331 cloud.maxwellianMostProbableSpeed
338 scalar sCosTheta = (faceVelocity &
n)/mostProbableSpeed;
341 scalar uNormProbCoeffA =
342 sCosTheta +
sqrt(
sqr(sCosTheta) + 2.0);
344 scalar uNormProbCoeffB =
348 + sCosTheta*(sCosTheta -
sqrt(
sqr(sCosTheta) + 2.0))
352 scalar randomScaling = 3.0;
356 randomScaling =
mag(sCosTheta) + 1;
364 scalar uNormalThermal;
372 uNormal = uNormalThermal + sCosTheta;
380 P = 2.0*uNormal/uNormProbCoeffA
381 *
exp(uNormProbCoeffB -
sqr(uNormalThermal));
392 + (t1 & faceVelocity)*t1
393 + (t2 & faceVelocity)*t2
394 + mostProbableSpeed*uNormal*
n;
396 scalar Ei = cloud.equipartitionInternalEnergy
399 cloud.constProps(typeId).internalDegreesOfFreedom()
402 cloud.addNewParcel(
p, celli,
U, Ei, typeId);
410 Info<<
" Particles inserted = "
void size(const label n)
Older name for setAddressableSize.
DSMCCloud< dsmcParcel > CloudType
scalar deltaTValue() const noexcept
Return time step value.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
void transfer(List< T > &list)
Transfer the contents of the argument List into this list and annul the argument list.
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.
const dictionary & coeffDict() const
Return the coefficients dictionary.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Templated inflow boundary model class.
constexpr char nl
The newline '\n' character (0x0a)
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
label k
Boltzmann constant.
virtual void autoMap(const mapPolyMesh &)
Remap the particles to the correct cells following mesh change.
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a T. FatalIOError if not found, or if the number of tokens is incorrect.
const Time & time() const
Return the top-level database.
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Perform reduction on a copy, using specified binary operation.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
wordList toc() const
Return the table of contents.
virtual const pointField & points() const
Return raw points.
#define forAll(list, i)
Loop across all elements in list.
virtual void inflow()
Introduce particles.
Type sample01()
Return a sample whose components lie in the range [0,1].
FreeStream(const dictionary &dict, CloudType &cloud)
Construct from dictionary.
void setSize(const label n)
Alias for resize()
dimensionedScalar exp(const dimensionedScalar &ds)
const polyBoundaryMesh & boundaryMesh() const noexcept
Return boundary mesh.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
iterator find(const Key &key)
Find and return an iterator set at the hashed entry.
A cloud is a registry collection of lagrangian particles.
virtual const labelList & faceOwner() const
Return face owner.
constexpr scalar pi(M_PI)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
errorManip< error > abort(error &err)
dimensionedScalar erf(const dimensionedScalar &ds)
vector point
Point is a vector.
const polyBoundaryMesh & patches
const std::string patch
OpenFOAM patch number as a std::string.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Mesh consisting of general polyhedral cells.
Type GaussNormal()
Return a sample whose components are normally distributed with zero mean and unity variance N(0...
A patch is a list of labels that address the faces in the global face list.
virtual ~FreeStream()
Destructor.
Templated base class for dsmc cloud.
static List< tetIndices > faceTetIndices(const polyMesh &mesh, label fI, label cI)
Return the tet decomposition of the given face, with.
static constexpr const zero Zero
Global zero (0)