36 template<
class CloudType>
39 const scalar a = 0.147;
41 scalar
h =
log(1.0 -
y*
y)/a;
55 template<
class CloudType>
59 const objectRegistry& obr = this->owner().
mesh();
71 <<
"Turbulence model not found in mesh database" <<
nl
72 <<
"Database objects include: " << obr.sortedToc()
75 return tmp<volScalarField>(
nullptr);
82 template<
class CloudType>
92 lambda_(this->coeffs().template lookup<scalar>(
"lambda")),
93 turbulence_(
readBool(this->coeffs().lookup(
"turbulence"))),
99 template<
class CloudType>
106 rndGen_(bmf.rndGen_),
107 lambda_(bmf.lambda_),
108 turbulence_(bmf.turbulence_),
116 template<
class CloudType>
123 template<
class CloudType>
154 template<
class CloudType>
158 const typename CloudType::parcelType::trackingData& td,
167 const scalar dp =
p.d();
168 const scalar Tc = td.Tc();
170 const scalar
alpha = 2.0*lambda_/dp;
179 const label celli =
p.cell();
181 const scalar kc =
k[celli];
205 Random& rnd = this->owner().rndGen();
208 const scalar u = 2*rnd.
scalar01() - 1;
210 const scalar a =
sqrt(1 -
sqr(u));
Calculates particle Brownian motion force.
BrownianMotionForce(CloudType &owner, const fvMesh &mesh, const dictionary &dict)
Construct from mesh.
virtual void cacheFields(const bool store)
Cache fields.
virtual ~BrownianMotionForce()
Destructor.
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.
Templated base class for dsmc cloud.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Generic GeometricField class.
Abstract base class for particle forces.
const fvMesh & mesh() const
Return the mesh database.
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.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Helper container for force Su and Sp terms.
const vector & Su() const
Return const access to the explicit contribution [kg m/s^2].
Mesh data needed to do the Finite Volume discretisation.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
A class for managing temporary objects.
bool isTmp() const
Return true if this is really a temporary object.
T * ptr() const
Return tmp pointer for reuse.
Template functions to aid in the implementation of demand driven data.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Fundamental dimensioned constants.
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
compressibleMomentumTransportModel momentumTransportModel
const scalar twoPi(2 *pi)
const dimensionedScalar h
Planck constant.
dimensionedScalar pow5(const dimensionedScalar &ds)
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)
errorManip< error > abort(error &err)
dimensionedScalar sin(const dimensionedScalar &ds)
void deleteDemandDrivenData(DataPtr &dataPtr)
dimensionedScalar log(const dimensionedScalar &ds)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
dimensionedScalar cos(const dimensionedScalar &ds)
scalarField Re(const UList< complex > &cf)