ThermoSurfaceFilm.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2021 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "ThermoSurfaceFilm.H"
27 #include "thermoSingleLayer.H"
29 #include "mathematicalConstants.H"
30 #include "Pstream.H"
31 
32 using namespace Foam::constant::mathematical;
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 template<class CloudType>
38 (
39  IStringStream
40  (
41  "(absorb bounce splashBai)"
42  )()
43 );
44 
45 
46 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
47 
48 template<class CloudType>
51 {
52  forAll(interactionTypeNames_, i)
53  {
54  if (interactionTypeNames_[i] == it)
55  {
56  return interactionType(i);
57  }
58  }
59 
61  << "Unknown interaction type " << it
62  << ". Valid interaction types include: " << interactionTypeNames_
63  << abort(FatalError);
64 
65  return interactionType(0);
66 }
67 
68 
69 template<class CloudType>
71 (
72  const interactionType& it
73 ) const
74 {
75  if (it >= interactionTypeNames_.size())
76  {
78  << "Unknown interaction type enumeration" << abort(FatalError);
79  }
80 
81  return interactionTypeNames_[it];
82 }
83 
84 
85 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
86 
87 template<class CloudType>
89 (
90  const vector& tanVec1,
91  const vector& tanVec2,
92  const vector& nf
93 ) const
94 {
95  // Azimuthal angle [rad]
96  const scalar phiSi = twoPi*rndGen_.sample01<scalar>();
97 
98  // Ejection angle [rad]
99  const scalar thetaSi = pi/180.0*(rndGen_.sample01<scalar>()*(50 - 5) + 5);
100 
101  // Direction vector of new parcel
102  const scalar alpha = sin(thetaSi);
103  const scalar dcorr = cos(thetaSi);
104  const vector normal = alpha*(tanVec1*cos(phiSi) + tanVec2*sin(phiSi));
105  vector dirVec = dcorr*nf;
106  dirVec += normal;
107 
108  return dirVec/mag(dirVec);
109 }
110 
111 
112 template<class CloudType>
114 (
116  const parcelType& p,
117  const polyPatch& pp,
118  const label facei,
119  const scalar mass,
120  bool& keepParticle
121 )
122 {
123  if (debug)
124  {
125  Info<< "Parcel " << p.origId() << " absorbInteraction" << endl;
126  }
127 
128  const fluidThermo& carrierThermo =
129  static_cast<const ThermoCloud<CloudType>&>(this->owner())
130  .carrierThermo();
131 
132  const parcelThermo& thermo =
133  static_cast<const ThermoCloud<CloudType>&>(this->owner()).thermo();
134 
135  // Patch face normal
136  const vector& nf = pp.faceNormals()[facei];
137 
138  // Patch velocity
139  const vector& Up = this->owner().U().boundaryField()[pp.index()][facei];
140 
141  // Relative parcel velocity
142  const vector Urel = p.U() - Up;
143 
144  // Parcel normal velocity
145  const vector Un = nf*(Urel & nf);
146 
147  // Parcel tangential velocity
148  const vector Ut = Urel - Un;
149 
150  const liquidProperties& liq = thermo.liquids().properties()[0];
151 
152  // Local pressure
153  const scalar pc = carrierThermo.p()[p.cell()];
154 
155  filmModel.addSources
156  (
157  pp.index(),
158  facei,
159  mass, // mass
160  mass*Ut, // tangential momentum
161  mass*mag(Un), // impingement pressure
162  mass*liq.Hs(pc, p.T()) // energy
163  );
164 
165  this->nParcelsTransferred()++;
166 
167  keepParticle = false;
168 }
169 
170 
171 template<class CloudType>
173 (
174  parcelType& p,
175  const polyPatch& pp,
176  const label facei,
177  bool& keepParticle
178 ) const
179 {
180  if (debug)
181  {
182  Info<< "Parcel " << p.origId() << " bounceInteraction" << endl;
183  }
184 
185  // Patch face normal
186  const vector& nf = pp.faceNormals()[facei];
187 
188  // Patch velocity
189  const vector& Up = this->owner().U().boundaryField()[pp.index()][facei];
190 
191  // Relative parcel velocity
192  const vector Urel = p.U() - Up;
193 
194  // Flip parcel normal velocity component
195  p.U() -= 2.0*nf*(Urel & nf);
196 
197  keepParticle = true;
198 }
199 
200 
201 template<class CloudType>
203 (
205  const parcelType& p,
206  const polyPatch& pp,
207  const label facei,
208  bool& keepParticle
209 )
210 {
211  if (debug)
212  {
213  Info<< "Parcel " << p.origId() << " drySplashInteraction" << endl;
214  }
215 
216  const fluidThermo& carrierThermo =
217  static_cast<const ThermoCloud<CloudType>&>(this->owner())
218  .carrierThermo();
219 
220  const parcelThermo& thermo =
221  static_cast<const ThermoCloud<CloudType>&>(this->owner()).thermo();
222 
223  const liquidProperties& liq = thermo.liquids().properties()[0];
224 
225  // Patch face velocity and normal
226  const vector& Up = this->owner().U().boundaryField()[pp.index()][facei];
227  const vector& nf = pp.faceNormals()[facei];
228 
229  // Local pressure
230  const scalar pc = carrierThermo.p()[p.cell()];
231 
232  // Retrieve parcel properties
233  const scalar m = p.mass()*p.nParticle();
234  const scalar rho = p.rho();
235  const scalar d = p.d();
236  const scalar sigma = liq.sigma(pc, p.T());
237  const scalar mu = liq.mu(pc, p.T());
238  const vector Urel = p.U() - Up;
239  const vector Un = nf*(Urel & nf);
240 
241  // Laplace number
242  const scalar La = rho*sigma*d/sqr(mu);
243 
244  // Weber number
245  const scalar We = rho*magSqr(Un)*d/sigma;
246 
247  // Critical Weber number
248  const scalar Wec = Adry_*pow(La, -0.183);
249 
250  if (We < Wec) // Adhesion - assume absorb
251  {
252  absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
253  }
254  else // Splash
255  {
256  // Ratio of incident mass to splashing mass
257  const scalar mRatio = 0.2 + 0.6*rndGen_.sample01<scalar>();
258  splashInteraction
259  (filmModel, p, pp, facei, mRatio, We, Wec, sigma, keepParticle);
260  }
261 }
262 
263 
264 template<class CloudType>
266 (
268  parcelType& p,
269  const polyPatch& pp,
270  const label facei,
271  bool& keepParticle
272 )
273 {
274  if (debug)
275  {
276  Info<< "Parcel " << p.origId() << " wetSplashInteraction" << endl;
277  }
278 
279  const fluidThermo& carrierThermo =
280  static_cast<const ThermoCloud<CloudType>&>(this->owner())
281  .carrierThermo();
282 
283  const parcelThermo& thermo =
284  static_cast<const ThermoCloud<CloudType>&>(this->owner()).thermo();
285 
286  const liquidProperties& liq = thermo.liquids().properties()[0];
287 
288  // Patch face velocity and normal
289  const vector& Up = this->owner().U().boundaryField()[pp.index()][facei];
290  const vector& nf = pp.faceNormals()[facei];
291 
292  // Local pressure
293  const scalar pc = carrierThermo.p()[p.cell()];
294 
295  // Retrieve parcel properties
296  const scalar m = p.mass()*p.nParticle();
297  const scalar rho = p.rho();
298  const scalar d = p.d();
299  vector& U = p.U();
300  const scalar sigma = liq.sigma(pc, p.T());
301  const scalar mu = liq.mu(pc, p.T());
302  const vector Urel = p.U() - Up;
303  const vector Un = nf*(Urel & nf);
304  const vector Ut = Urel - Un;
305 
306  // Laplace number
307  const scalar La = rho*sigma*d/sqr(mu);
308 
309  // Weber number
310  const scalar We = rho*magSqr(Un)*d/sigma;
311 
312  // Critical Weber number
313  const scalar Wec = Awet_*pow(La, -0.183);
314 
315  if (We < 2) // Adhesion - assume absorb
316  {
317  absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
318  }
319  else if ((We >= 2) && (We < 20)) // Bounce
320  {
321  // Incident angle of impingement
322  const scalar theta = pi/2 - acos(U/mag(U) & nf);
323 
324  // Restitution coefficient
325  const scalar epsilon = 0.993 - theta*(1.76 - theta*(1.56 - theta*0.49));
326 
327  // Update parcel velocity
328  U = -epsilon*(Un) + 5.0/7.0*(Ut);
329 
330  keepParticle = true;
331  return;
332  }
333  else if ((We >= 20) && (We < Wec)) // Spread - assume absorb
334  {
335  absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
336  }
337  else // Splash
338  {
339  // Ratio of incident mass to splashing mass
340  // splash mass can be > incident mass due to film entrainment
341  const scalar mRatio = 0.2 + 0.9*rndGen_.sample01<scalar>();
342  splashInteraction
343  (filmModel, p, pp, facei, mRatio, We, Wec, sigma, keepParticle);
344  }
345 }
346 
347 
348 template<class CloudType>
350 (
352  const parcelType& p,
353  const polyPatch& pp,
354  const label facei,
355  const scalar mRatio,
356  const scalar We,
357  const scalar Wec,
358  const scalar sigma,
359  bool& keepParticle
360 )
361 {
362  // Patch face velocity and normal
363  const fvMesh& mesh = this->owner().mesh();
364  const vector& Up = this->owner().U().boundaryField()[pp.index()][facei];
365  const vector& nf = pp.faceNormals()[facei];
366 
367  // Determine direction vectors tangential to patch normal
368  const vector tanVec1 = normalised(perpendicular(nf));
369  const vector tanVec2 = nf^tanVec1;
370 
371  // Retrieve parcel properties
372  const scalar np = p.nParticle();
373  const scalar m = p.mass()*np;
374  const scalar d = p.d();
375  const vector Urel = p.U() - Up;
376  const vector Un = nf*(Urel & nf);
377  const vector Ut = Urel - Un;
378  const vector& posC = mesh.C()[p.cell()];
379  const vector& posCf = mesh.Cf().boundaryField()[pp.index()][facei];
380 
381  // Total mass of (all) splashed parcels
382  const scalar mSplash = m*mRatio;
383 
384  // Number of splashed particles per incoming particle
385  const scalar Ns = 5.0*(We/Wec - 1.0);
386 
387  // Average diameter of splashed particles
388  const scalar dBarSplash = 1/cbrt(6.0)*cbrt(mRatio/Ns)*d + rootVSmall;
389 
390  // Cumulative diameter splash distribution
391  const scalar dMax = 0.9*cbrt(mRatio)*d;
392  const scalar dMin = 0.1*dMax;
393  const scalar K = exp(-dMin/dBarSplash) - exp(-dMax/dBarSplash);
394 
395  // Surface energy of secondary parcels [J]
396  scalar ESigmaSec = 0;
397 
398  // Sample splash distribution to determine secondary parcel diameters
399  scalarList dNew(parcelsPerSplash_);
400  scalarList npNew(parcelsPerSplash_);
401  forAll(dNew, i)
402  {
403  const scalar y = rndGen_.sample01<scalar>();
404  dNew[i] = -dBarSplash*log(exp(-dMin/dBarSplash) - y*K);
405  npNew[i] = mRatio*np*pow3(d)/pow3(dNew[i])/parcelsPerSplash_;
406  ESigmaSec += npNew[i]*sigma*p.areaS(dNew[i]);
407  }
408 
409  // Incident kinetic energy [J]
410  const scalar EKIn = 0.5*m*magSqr(Un);
411 
412  // Incident surface energy [J]
413  const scalar ESigmaIn = np*sigma*p.areaS(d);
414 
415  // Dissipative energy
416  const scalar Ed = max(0.8*EKIn, np*Wec/12*pi*sigma*sqr(d));
417 
418  // Total energy [J]
419  const scalar EKs = EKIn + ESigmaIn - ESigmaSec - Ed;
420 
421  // Switch to absorb if insufficient energy for splash
422  if (EKs <= 0)
423  {
424  absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
425  return;
426  }
427 
428  // Helper variables to calculate magUns0
429  const scalar logD = log(d);
430  const scalar coeff2 = log(dNew[0]) - logD + rootVSmall;
431  scalar coeff1 = 0.0;
432  forAll(dNew, i)
433  {
434  coeff1 += sqr(log(dNew[i]) - logD);
435  }
436 
437  // Magnitude of the normal velocity of the first splashed parcel
438  const scalar magUns0 =
439  sqrt(2.0*parcelsPerSplash_*EKs/mSplash/(1.0 + coeff1/sqr(coeff2)));
440 
441  // Set splashed parcel properties
442  forAll(dNew, i)
443  {
444  const vector dirVec = splashDirection(tanVec1, tanVec2, -nf);
445 
446  // Create a new parcel by copying source parcel
447  parcelType* pPtr = new parcelType(p);
448 
449  pPtr->origId() = pPtr->getNewParticleID();
450 
451  pPtr->origProc() = Pstream::myProcNo();
452 
453  if (splashParcelType_ >= 0)
454  {
455  pPtr->typeId() = splashParcelType_;
456  }
457 
458  // Perturb new parcels towards the owner cell centre
459  pPtr->track(0.5*rndGen_.sample01<scalar>()*(posC - posCf), 0);
460 
461  pPtr->nParticle() = npNew[i];
462 
463  pPtr->d() = dNew[i];
464 
465  pPtr->U() = dirVec*(mag(Cf_*Ut) + magUns0*(log(dNew[i]) - logD)/coeff2);
466 
467  // Apply correction to velocity for 2-D cases
468  meshTools::constrainDirection(mesh, mesh.solutionD(), pPtr->U());
469 
470  // Add the new parcel
471  this->owner().addParticle(pPtr);
472 
473  nParcelsSplashed_++;
474  }
475 
476  // Transfer remaining part of parcel to film 0 - splashMass can be -ve
477  // if entraining from the film
478  const scalar mDash = m - mSplash;
479  absorbInteraction(filmModel, p, pp, facei, mDash, keepParticle);
480 }
481 
482 
483 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
484 
485 template<class CloudType>
487 (
488  const dictionary& dict,
489  CloudType& owner
490 )
491 :
492  SurfaceFilmModel<CloudType>(dict, owner, typeName),
493  rndGen_(owner.rndGen()),
494  TFilmPatch_(0),
495  CpFilmPatch_(0),
496  interactionType_
497  (
498  interactionTypeEnum(this->coeffDict().lookup("interactionType"))
499  ),
500  deltaWet_(0.0),
501  splashParcelType_(0),
502  parcelsPerSplash_(0),
503  Adry_(0.0),
504  Awet_(0.0),
505  Cf_(0.0),
506  nParcelsSplashed_(0)
507 {
508  Info<< " Applying " << interactionTypeStr(interactionType_)
509  << " interaction model" << endl;
510 
511  if (interactionType_ == itSplashBai)
512  {
513  this->coeffDict().lookup("deltaWet") >> deltaWet_;
514  splashParcelType_ =
515  this->coeffDict().lookupOrDefault("splashParcelType", -1);
516  parcelsPerSplash_ =
517  this->coeffDict().lookupOrDefault("parcelsPerSplash", 2);
518  this->coeffDict().lookup("Adry") >> Adry_;
519  this->coeffDict().lookup("Awet") >> Awet_;
520  this->coeffDict().lookup("Cf") >> Cf_;
521  }
522 }
523 
524 
525 template<class CloudType>
527 (
529 )
530 :
532  rndGen_(sfm.rndGen_),
533  TFilmPatch_(sfm.TFilmPatch_),
534  CpFilmPatch_(sfm.CpFilmPatch_),
535  interactionType_(sfm.interactionType_),
536  deltaWet_(sfm.deltaWet_),
537  splashParcelType_(sfm.splashParcelType_),
538  parcelsPerSplash_(sfm.parcelsPerSplash_),
539  Adry_(sfm.Adry_),
540  Awet_(sfm.Awet_),
541  Cf_(sfm.Cf_),
542  nParcelsSplashed_(sfm.nParcelsSplashed_)
543 {}
544 
545 
546 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
547 
548 template<class CloudType>
550 {}
551 
552 
553 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
554 
555 template<class CloudType>
557 (
558  parcelType& p,
559  const polyPatch& pp,
560  bool& keepParticle
561 )
562 {
563  // Retrieve the film model from the owner database
566  (
567  this->owner().db().time().objectRegistry::template
568  lookupObject
570  (
571  "surfaceFilmProperties"
572  )
573  );
574 
575  const label patchi = pp.index();
576 
577  if (filmModel.isRegionPatch(patchi))
578  {
579  const label facei = pp.whichFace(p.face());
580 
581  switch (interactionType_)
582  {
583  case itBounce:
584  {
585  bounceInteraction(p, pp, facei, keepParticle);
586 
587  break;
588  }
589  case itAbsorb:
590  {
591  const scalar m = p.nParticle()*p.mass();
592  absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
593 
594  break;
595  }
596  case itSplashBai:
597  {
598  bool dry = this->deltaFilmPatch_[patchi][facei] < deltaWet_;
599 
600  if (dry)
601  {
602  drySplashInteraction(filmModel, p, pp, facei, keepParticle);
603  }
604  else
605  {
606  wetSplashInteraction(filmModel, p, pp, facei, keepParticle);
607  }
608 
609  break;
610  }
611  default:
612  {
614  << "Unknown interaction type enumeration"
615  << abort(FatalError);
616  }
617  }
618 
619  // Transfer parcel/parcel interactions complete
620  return true;
621  }
622 
623  // Parcel not interacting with film
624  return false;
625 }
626 
627 
628 template<class CloudType>
630 (
631  const label filmPatchi,
632  const label primaryPatchi,
634 )
635 {
637  (
638  filmPatchi,
639  primaryPatchi,
640  filmModel
641  );
642 
644  refCast<const regionModels::surfaceFilmModels::thermoSingleLayer>
645  (
646  filmModel
647  );
648 
649  TFilmPatch_ = thermalFilmModel.thermo().T().boundaryField()[filmPatchi];
650  filmModel.toPrimary(filmPatchi, TFilmPatch_);
651 
652  CpFilmPatch_ =
653  thermalFilmModel.thermo().Cpv()().boundaryField()[filmPatchi];
654  filmModel.toPrimary(filmPatchi, CpFilmPatch_);
655 }
656 
657 
658 template<class CloudType>
660 (
661  parcelType& p,
662  const label filmFacei
663 ) const
664 {
666 
667  // Set parcel properties
668  p.T() = TFilmPatch_[filmFacei];
669  p.Cp() = CpFilmPatch_[filmFacei];
670 }
671 
672 
673 template<class CloudType>
675 {
677 
678  label nSplash0 = this->template getModelProperty<label>("nParcelsSplashed");
679  label nSplashTotal =
680  nSplash0 + returnReduce(nParcelsSplashed_, sumOp<label>());
681 
682  os << " New film splash parcels = " << nSplashTotal << endl;
683 
684  if (this->writeTime())
685  {
686  this->setModelProperty("nParcelsSplashed", nSplashTotal);
687  nParcelsSplashed_ = 0;
688  }
689 }
690 
691 
692 // ************************************************************************* //
scalar Awet_
Wet surface roughness coefficient.
virtual void info(Ostream &os)
Write surface film info to stream.
Thermo parcel surface film model.
virtual scalar Hs(scalar p, scalar T) const =0
Liquid sensible enthalpy [J/kg].
dictionary dict
dimensionedScalar acos(const dimensionedScalar &ds)
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
virtual void addSources(const label patchi, const label facei, const scalar massSource, const vector &momentumSource, const scalar pressureSource, const scalar energySource)=0
External hook to add sources to the film.
scalar Cf_
Skin friction typically in the range 0.6 < Cf < 0.8.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
fluidReactionThermo & thermo
Definition: createFields.H:28
dimensionedScalar log(const dimensionedScalar &ds)
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual void setParcelProperties(parcelType &p, const label filmFacei) const
Set the individual parcel properties.
bool isRegionPatch(const label primaryPatchi) const
Return true if patchi on the primary region is a coupled.
Definition: regionModelI.H:149
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
const surfaceVectorField & Cf() const
Return face centres as surfaceVectorField.
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package, and provides:
Definition: parcelThermo.H:58
const Boundary & boundaryField() const
Return const-reference to the boundary field.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const Vector< label > & solutionD() const
Return the vector of solved-for directions in mesh.
Definition: polyMesh.C:831
ThermoSurfaceFilm(const dictionary &dict, CloudType &owner)
Construct from components.
scalar Adry_
Dry surface roughness coefficient.
virtual ~ThermoSurfaceFilm()
Destructor.
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
scalarList TFilmPatch_
Film temperature / patch face.
CGAL::Exact_predicates_exact_constructions_kernel K
virtual volScalarField & p()=0
Pressure [Pa].
virtual tmp< volScalarField > Cpv() const =0
Heat capacity at constant pressure/volume [J/kg/K].
void drySplashInteraction(regionModels::surfaceFilmModels::surfaceFilmRegionModel &, const parcelType &p, const polyPatch &pp, const label facei, bool &keepParticle)
Parcel interaction with dry surface.
Macros for easy insertion into run-time selection tables.
Templated base class for thermodynamic cloud.
Definition: ThermoCloud.H:76
virtual void info(Ostream &os)
Write surface film info to stream.
virtual scalar sigma(scalar p, scalar T) const =0
Surface tension [N/m].
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m^2/K^4].
scalar y
void toPrimary(const label regionPatchi, List< Type > &regionField) const
Convert a local region field to the primary region.
interactionType interactionType_
Interaction type enumeration.
void absorbInteraction(regionModels::surfaceFilmModels::surfaceFilmRegionModel &, const parcelType &p, const polyPatch &pp, const label facei, const scalar mass, bool &keepParticle)
Absorb parcel into film.
stressControl lookup("compactNormalStress") >> compactNormalStress
dynamicFvMesh & mesh
Form normalised(const VectorSpace< Form, Cmpt, Ncmpts > &vs)
Definition: VectorSpaceI.H:413
dimensionedScalar cos(const dimensionedScalar &ds)
dimensionedScalar exp(const dimensionedScalar &ds)
interactionType interactionTypeEnum(const word &it) const
static wordList interactionTypeNames_
Word descriptions of interaction type names.
mathematical constants.
A class for handling words, derived from string.
Definition: word.H:59
Base-class for fluid thermodynamic properties.
Definition: fluidThermo.H:53
dimensionedScalar cbrt(const dimensionedScalar &ds)
The thermophysical properties of a liquid.
virtual void cacheFilmFields(const label filmPatchi, const label primaryPatchi, const regionModels::surfaceFilmModels::surfaceFilmRegionModel &filmModel)
Cache the film fields in preparation for injection.
const Field< PointType > & faceNormals() const
Return face normals for patch.
Urel
Definition: pEqn.H:60
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const scalar twoPi(2 *pi)
vector splashDirection(const vector &tanVec1, const vector &tanVec2, const vector &nf) const
Return splashed parcel direction.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
dimensionedScalar sin(const dimensionedScalar &ds)
const dimensionedScalar mu
Atomic mass unit.
virtual void setParcelProperties(parcelType &p, const label filmFacei) const
Set the individual parcel properties.
void bounceInteraction(parcelType &p, const polyPatch &pp, const label facei, bool &keepParticle) const
Bounce parcel (flip parcel normal velocity)
Vector< Cmpt > perpendicular(const Vector< Cmpt > &v)
Definition: VectorI.H:166
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
label patchi
const liquidMixtureProperties & liquids() const
Return reference to the global (additional) liquids.
Definition: parcelThermo.C:82
U
Definition: pEqn.H:72
virtual const volScalarField & T() const =0
Temperature [K].
CloudType::parcelType parcelType
Convenience typedef to the cloud&#39;s parcel type.
void splashInteraction(regionModels::surfaceFilmModels::surfaceFilmRegionModel &, const parcelType &p, const polyPatch &pp, const label facei, const scalar mRatio, const scalar We, const scalar Wec, const scalar sigma, bool &keepParticle)
Bai parcel splash interaction model.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
label index() const
Return the index of this patch in the boundaryMesh.
label parcelsPerSplash_
Number of new parcels resulting from splash event.
scalar epsilon
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual void cacheFilmFields(const label filmPatchi, const label primaryPatchi, const regionModels::surfaceFilmModels::surfaceFilmRegionModel &)
Cache the film fields in preparation for injection.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
label nParcelsSplashed_
Counter for number of new splash parcels.
Random & rndGen_
Reference to the cloud random number generator.
void constrainDirection(const polyMesh &mesh, const Vector< label > &dirs, vector &d)
Set the constrained components of directions/velocity to zero.
Definition: meshTools.C:671
const volVectorField & C() const
Return cell centres as volVectorField.
const PtrList< liquidProperties > & properties() const
Return the liquid properties.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
word interactionTypeStr(const interactionType &it) const
virtual scalar mu(scalar p, scalar T) const =0
Liquid viscosity [Pa s].
virtual bool transferParcel(parcelType &p, const polyPatch &pp, bool &keepParticle)
Transfer parcel from cloud to surface film.
void wetSplashInteraction(regionModels::surfaceFilmModels::surfaceFilmRegionModel &, parcelType &p, const polyPatch &pp, const label facei, bool &keepParticle)
Parcel interaction with wetted surface.
scalar deltaWet_
Film thickness beyond which patch is assumed to be wet.
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:75
label splashParcelType_
Splash parcel type label - id assigned to identify parcel for.
scalarList CpFilmPatch_
Film specific heat capacity / patch face.
Thermodynamic form of single-cell layer surface film model.
label whichFace(const label l) const
Return label of face in patch from global face label.
Definition: polyPatch.H:389