Go to the documentation of this file.
38 template<
class CloudType>
41 const scalar a = 0.147;
43 scalar
h =
log(1.0 -
y*
y)/a;
57 template<
class CloudType>
61 const objectRegistry& obr = this->owner().
mesh();
77 <<
"Turbulence model not found in mesh database" <<
nl
78 <<
"Database objects include: " << obr.sortedToc()
87 template<
class CloudType>
97 lambda_(this->coeffs().getScalar(
"lambda")),
98 turbulence_(this->coeffs().getBool(
"turbulence")),
104 template<
class CloudType>
111 rndGen_(bmf.rndGen_),
112 lambda_(bmf.lambda_),
113 turbulence_(bmf.turbulence_),
121 template<
class CloudType>
128 template<
class CloudType>
159 template<
class CloudType>
163 const typename CloudType::parcelType::trackingData& td,
172 const scalar dp =
p.d();
173 const scalar Tc = td.Tc();
175 const scalar
alpha = 2.0*lambda_/dp;
184 const label celli =
p.cell();
186 const scalar kc =
k[celli];
213 const scalar u = 2*rnd.
sample01<scalar>() - 1;
215 const scalar a =
sqrt(1 -
sqr(u));
BrownianMotionForce(CloudType &owner, const fvMesh &mesh, const dictionary &dict)
Construct from mesh.
constexpr const char *const group
Group name for atomic constants.
const dimensionedScalar k
Boltzmann constant.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
virtual forceSuSp calcCoupled(const typename CloudType::parcelType &p, const typename CloudType::parcelType::trackingData &td, const scalar dt, const scalar mass, const scalar Re, const scalar muc) const
Calculate the coupled force.
Different types of constants.
Template functions to aid in the implementation of demand driven data.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
dimensionedScalar sin(const dimensionedScalar &ds)
virtual ~BrownianMotionForce()
Destructor.
Type sample01()
Return a sample whose components lie in the range [0,1].
static const word propertiesName
Default name of the turbulence properties dictionary.
virtual void cacheFields(const bool store)
Cache fields.
bool isTmp() const noexcept
Identical to is_pointer()
const Type & value() const
Return const reference to value.
dimensionedScalar exp(const dimensionedScalar &ds)
const vector & Su() const
Return const access to the explicit contribution [kg.m/s2].
Random & rndGen()
Return reference to the random object.
void deleteDemandDrivenData(DataPtr &dataPtr)
ThermalDiffusivity< CompressibleTurbulenceModel< fluidThermo > > turbulenceModel
const dimensionedScalar h
Planck constant.
const fvMesh & mesh() const
Return reference to the mesh.
Helper container for force Su and Sp terms.
Abstract base class for particle forces.
Templated base class for dsmc cloud.
constexpr scalar twoPi(2 *M_PI)
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
dimensionedScalar pow5(const dimensionedScalar &ds)
dimensionedScalar log(const dimensionedScalar &ds)
Mesh data needed to do the Finite Volume discretisation.
T * ptr() const
Return managed pointer for reuse, or clone() the object reference.
errorManip< error > abort(error &err)
Calculates particle Brownian motion force.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
constexpr scalar pi(M_PI)
Fundamental dimensioned constants.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label k
Boltzmann constant.
scalarField Re(const UList< complex > &cf)
Extract real component.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
const fvMesh & mesh() const
Return the mesh database.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
compressible::turbulenceModel & turb
dimensionedScalar cos(const dimensionedScalar &ds)