30 template<
class CloudType>
44 const typename CloudType::parcelType& p = iter();
50 if (useEquivalentSize_)
52 dEff *=
cbrt(p.nParticle()*volumeFactor_);
55 rMin =
min(dEff, rMin);
57 rhoMax =
max(p.rho(), rhoMax);
61 mag(p.U()) +
mag(p.omega())*dEff/2,
73 template<
class CloudType>
76 typename CloudType::parcelType& p,
78 const WallSiteData<vector>& data,
84 vector r_PW = p.position() - site;
86 vector U_PW = p.U() - data.wallData();
88 scalar r_PW_mag =
mag(r_PW);
90 scalar normalOverlapMag =
max(pREff - r_PW_mag, 0.0);
92 vector rHat_PW = r_PW/(r_PW_mag + vSmall);
94 scalar etaN = alpha_*
sqrt(p.mass()*kN)*
pow025(normalOverlapMag);
98 *(kN*
pow(normalOverlapMag, b_) - etaN*(U_PW & rHat_PW));
105 -cohesionEnergyDensity_
113 U_PW - (U_PW & rHat_PW)*rHat_PW
114 + (p.omega() ^ (pREff*-rHat_PW));
116 scalar deltaT = this->owner().mesh().time().deltaTValue();
118 vector& tangentialOverlap_PW =
119 p.collisionRecords().matchWallRecord(-r_PW, pREff).collisionData();
121 tangentialOverlap_PW += USlip_PW*deltaT;
123 scalar tangentialOverlapMag =
mag(tangentialOverlap_PW);
125 if (tangentialOverlapMag > vSmall)
127 scalar kT = 8.0*
sqrt(pREff*normalOverlapMag)*Gstar_;
134 if (kT*tangentialOverlapMag > mu_*
mag(fN_PW))
139 fT_PW = -mu_*
mag(fN_PW)*USlip_PW/
mag(USlip_PW);
141 tangentialOverlap_PW =
Zero;
145 fT_PW = - kT*tangentialOverlap_PW - etaT*USlip_PW;
150 p.torque() += (pREff*-rHat_PW) ^ fT_PW;
157 template<
class CloudType>
167 alpha_(this->coeffDict().
template lookup<scalar>(
"alpha")),
168 b_(this->coeffDict().
template lookup<scalar>(
"b")),
169 mu_(this->coeffDict().
template lookup<scalar>(
"mu")),
170 cohesionEnergyDensity_
172 this->coeffDict().
template lookup<scalar>(
"cohesionEnergyDensity")
175 collisionResolutionSteps_
177 this->coeffDict().
template lookup<scalar>(
"collisionResolutionSteps")
180 useEquivalentSize_(
Switch(this->coeffDict().
lookup(
"useEquivalentSize")))
182 if (useEquivalentSize_)
185 this->coeffDict().template lookup<scalar>(
"volumeFactor");
188 scalar nu = this->coeffDict().template lookup<scalar>(
"poissonsRatio");
190 scalar E = this->coeffDict().template lookup<scalar>(
"youngsModulus");
192 scalar pNu = this->owner().constProps().poissonsRatio();
194 scalar pE = this->owner().constProps().youngsModulus();
196 Estar_ = 1/((1 -
sqr(pNu))/pE + (1 -
sqr(nu))/E);
198 Gstar_ = 1/(2*((2 + pNu -
sqr(pNu))/pE + (2 + nu -
sqr(nu))/E));
200 cohesion_ = (
mag(cohesionEnergyDensity_) > vSmall);
206 template<
class CloudType>
213 template<
class CloudType>
219 if (useEquivalentSize_)
221 return p.d()/2*
cbrt(p.nParticle()*volumeFactor_);
230 template<
class CloudType>
237 template<
class CloudType>
240 if (!(this->owner().size()))
249 findMinMaxProperties(rMin, rhoMax, UMagMax);
252 scalar minCollisionDeltaT =
255 *
pow(rhoMax/(Estar_*
sqrt(UMagMax) + vSmall), 0.4)
256 /collisionResolutionSteps_;
258 return ceil(this->owner().time().deltaTValue()/minCollisionDeltaT);
262 template<
class CloudType>
272 scalar pREff = this->pREff(p);
274 scalar kN = (4.0/3.0)*
sqrt(pREff)*Estar_;
276 forAll(flatSitePoints, siteI)
281 flatSitePoints[siteI],
289 forAll(sharpSitePoints, siteI)
296 sharpSitePoints[siteI],
297 sharpSiteData[siteI],
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
DSMCCloud< dsmcParcel > CloudType
A list of keyword definitions, which are a keyword followed by any number of values (e...
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensionedScalar pow025(const dimensionedScalar &ds)
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none/any.
Forces between particles and walls, interacting with a spring, slider, damper model.
virtual label nSubCycles() const
For WallModels that control the timestep, calculate the.
Vector< scalar > vector
A scalar version of the templated Vector.
virtual ~WallSpringSliderDashpot()
Destructor.
stressControl lookup("compactNormalStress") >> compactNormalStress
dimensionedScalar cbrt(const dimensionedScalar &ds)
Templated wall interaction class.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Stores the patch ID and templated data to represent a collision with a wall to be passed to the wall ...
ParcelType parcelType
Type of parcel the cloud was instantiated for.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
WallSpringSliderDashpot(const dictionary &dict, CloudType &cloud)
Construct from dictionary.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
vector point
Point is a vector.
dimensioned< scalar > mag(const dimensioned< Type > &)
Templated base class for dsmc cloud.
virtual scalar pREff(const typename CloudType::parcelType &p) const
Return the effective radius for a particle for the model.
virtual bool controlsTimestep() const
Whether the WallModel has a timestep limit that will.