mixtureKEpsilon.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) 2013-2018 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 "mixtureKEpsilon.H"
27 #include "fvOptions.H"
28 #include "bound.H"
29 #include "twoPhaseSystem.H"
30 #include "dragModel.H"
31 #include "virtualMassModel.H"
34 #include "fvmSup.H"
35 
36 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
37 
38 namespace Foam
39 {
40 namespace RASModels
41 {
42 
43 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
44 
45 template<class BasicTurbulenceModel>
47 (
48  const alphaField& alpha,
49  const rhoField& rho,
50  const volVectorField& U,
51  const surfaceScalarField& alphaRhoPhi,
52  const surfaceScalarField& phi,
53  const transportModel& transport,
54  const word& propertiesName,
55  const word& type
56 )
57 :
59  (
60  type,
61  alpha,
62  rho,
63  U,
64  alphaRhoPhi,
65  phi,
66  transport,
67  propertiesName
68  ),
69 
70  liquidTurbulencePtr_(nullptr),
71 
72  Cmu_
73  (
75  (
76  "Cmu",
77  this->coeffDict_,
78  0.09
79  )
80  ),
81  C1_
82  (
84  (
85  "C1",
86  this->coeffDict_,
87  1.44
88  )
89  ),
90  C2_
91  (
93  (
94  "C2",
95  this->coeffDict_,
96  1.92
97  )
98  ),
99  C3_
100  (
102  (
103  "C3",
104  this->coeffDict_,
105  C2_.value()
106  )
107  ),
108  Cp_
109  (
111  (
112  "Cp",
113  this->coeffDict_,
114  0.25
115  )
116  ),
117  sigmak_
118  (
120  (
121  "sigmak",
122  this->coeffDict_,
123  1.0
124  )
125  ),
126  sigmaEps_
127  (
129  (
130  "sigmaEps",
131  this->coeffDict_,
132  1.3
133  )
134  ),
135 
136  k_
137  (
138  IOobject
139  (
140  IOobject::groupName("k", alphaRhoPhi.group()),
141  this->runTime_.timeName(),
142  this->mesh_,
145  ),
146  this->mesh_
147  ),
148  epsilon_
149  (
150  IOobject
151  (
152  IOobject::groupName("epsilon", alphaRhoPhi.group()),
153  this->runTime_.timeName(),
154  this->mesh_,
157  ),
158  this->mesh_
159  )
160 {
161  bound(k_, this->kMin_);
162  bound(epsilon_, this->epsilonMin_);
163 
164  if (type == typeName)
165  {
166  this->printCoeffs(type);
167  }
168 }
169 
170 
171 template<class BasicTurbulenceModel>
173 (
174  const volScalarField& epsilon
175 ) const
176 {
177  const volScalarField::Boundary& ebf = epsilon.boundaryField();
178 
179  wordList ebt = ebf.types();
180 
181  forAll(ebf, patchi)
182  {
183  if (isA<fixedValueFvPatchScalarField>(ebf[patchi]))
184  {
185  ebt[patchi] = fixedValueFvPatchScalarField::typeName;
186  }
187  }
188 
189  return ebt;
190 }
191 
192 
193 template<class BasicTurbulenceModel>
195 (
196  volScalarField& vsf,
197  const volScalarField& refVsf
198 ) const
199 {
200  volScalarField::Boundary& bf = vsf.boundaryFieldRef();
201  const volScalarField::Boundary& refBf =
202  refVsf.boundaryField();
203 
204  forAll(bf, patchi)
205  {
206  if
207  (
208  isA<inletOutletFvPatchScalarField>(bf[patchi])
209  && isA<inletOutletFvPatchScalarField>(refBf[patchi])
210  )
211  {
212  refCast<inletOutletFvPatchScalarField>
213  (bf[patchi]).refValue() =
214  refCast<const inletOutletFvPatchScalarField>
215  (refBf[patchi]).refValue();
216  }
217  }
218 }
219 
220 
221 template<class BasicTurbulenceModel>
223 {
224  if (rhom_.valid()) return;
225 
226  // Local references to gas-phase properties
227  const volScalarField& kg = this->k_;
228  const volScalarField& epsilong = this->epsilon_;
229 
230  // Local references to liquid-phase properties
231  mixtureKEpsilon<BasicTurbulenceModel>& turbc = this->liquidTurbulence();
232  const volScalarField& kl = turbc.k_;
233  const volScalarField& epsilonl = turbc.epsilon_;
234 
235  word startTimeName
236  (
237  this->runTime_.timeName(this->runTime_.startTime().value())
238  );
239 
240  Ct2_.set
241  (
242  new volScalarField
243  (
244  IOobject
245  (
246  "Ct2",
247  startTimeName,
248  this->mesh_,
251  ),
252  Ct2()
253  )
254  );
255 
256  rhom_.set
257  (
258  new volScalarField
259  (
260  IOobject
261  (
262  "rhom",
263  startTimeName,
264  this->mesh_,
267  ),
268  rhom()
269  )
270  );
271 
272  km_.set
273  (
274  new volScalarField
275  (
276  IOobject
277  (
278  "km",
279  startTimeName,
280  this->mesh_,
283  ),
284  mix(kl, kg),
285  kl.boundaryField().types()
286  )
287  );
288  correctInletOutlet(km_(), kl);
289 
290  epsilonm_.set
291  (
292  new volScalarField
293  (
294  IOobject
295  (
296  "epsilonm",
297  startTimeName,
298  this->mesh_,
301  ),
302  mix(epsilonl, epsilong),
303  epsilonBoundaryTypes(epsilonl)
304  )
305  );
306  correctInletOutlet(epsilonm_(), epsilonl);
307 }
308 
309 
310 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
311 
312 template<class BasicTurbulenceModel>
314 {
316  {
317  Cmu_.readIfPresent(this->coeffDict());
318  C1_.readIfPresent(this->coeffDict());
319  C2_.readIfPresent(this->coeffDict());
320  C3_.readIfPresent(this->coeffDict());
321  Cp_.readIfPresent(this->coeffDict());
322  sigmak_.readIfPresent(this->coeffDict());
323  sigmaEps_.readIfPresent(this->coeffDict());
324 
325  return true;
326  }
327  else
328  {
329  return false;
330  }
331 }
332 
333 
334 template<class BasicTurbulenceModel>
336 {
337  this->nut_ = Cmu_*sqr(k_)/epsilon_;
338  this->nut_.correctBoundaryConditions();
339  fv::options::New(this->mesh_).correct(this->nut_);
340 
341  BasicTurbulenceModel::correctNut();
342 }
343 
344 
345 template<class BasicTurbulenceModel>
348 {
349  if (!liquidTurbulencePtr_)
350  {
351  const volVectorField& U = this->U_;
352 
353  const transportModel& gas = this->transport();
354  const twoPhaseSystem& fluid =
355  refCast<const twoPhaseSystem>(gas.fluid());
356  const transportModel& liquid = fluid.otherPhase(gas);
357 
358  liquidTurbulencePtr_ =
360  (
362  (
364  (
366  liquid.name()
367  )
368  )
369  );
370  }
371 
372  return *liquidTurbulencePtr_;
373 }
374 
375 
376 template<class BasicTurbulenceModel>
378 {
379  const mixtureKEpsilon<BasicTurbulenceModel>& liquidTurbulence =
380  this->liquidTurbulence();
381 
382  const transportModel& gas = this->transport();
383  const twoPhaseSystem& fluid = refCast<const twoPhaseSystem>(gas.fluid());
384  const transportModel& liquid = fluid.otherPhase(gas);
385 
386  const volScalarField& alphag = this->alpha_;
387 
388  volScalarField magUr(mag(liquidTurbulence.U() - this->U()));
389 
390  volScalarField beta
391  (
392  (6*this->Cmu_/(4*sqrt(3.0/2.0)))
393  *fluid.Kd()/liquid.rho()
394  *(liquidTurbulence.k_/liquidTurbulence.epsilon_)
395  );
396  volScalarField Ct0((3 + beta)/(1 + beta + 2*gas.rho()/liquid.rho()));
397  volScalarField fAlphad((180 + (-4.71e3 + 4.26e4*alphag)*alphag)*alphag);
398 
399  return sqr(1 + (Ct0 - 1)*exp(-fAlphad));
400 }
401 
402 
403 template<class BasicTurbulenceModel>
405 {
406  const transportModel& gas = this->transport();
407  const twoPhaseSystem& fluid = refCast<const twoPhaseSystem>(gas.fluid());
408  return fluid.otherPhase(gas).rho();
409 }
410 
411 
412 template<class BasicTurbulenceModel>
414 {
415  const transportModel& gas = this->transport();
416  const twoPhaseSystem& fluid = refCast<const twoPhaseSystem>(gas.fluid());
417  const virtualMassModel& virtualMass =
418  fluid.lookupSubModel<virtualMassModel>(gas, fluid.otherPhase(gas));
419  return gas.rho() + virtualMass.Cvm()*fluid.otherPhase(gas).rho();
420 }
421 
422 
423 template<class BasicTurbulenceModel>
425 {
426  const volScalarField& alphag = this->alpha_;
427  const volScalarField& alphal = this->liquidTurbulence().alpha_;
428 
429  return alphal*rholEff() + alphag*rhogEff();
430 }
431 
432 
433 template<class BasicTurbulenceModel>
435 (
436  const volScalarField& fc,
437  const volScalarField& fd
438 ) const
439 {
440  const volScalarField& alphag = this->alpha_;
441  const volScalarField& alphal = this->liquidTurbulence().alpha_;
442 
443  return (alphal*rholEff()*fc + alphag*rhogEff()*fd)/rhom_();
444 }
445 
446 
447 template<class BasicTurbulenceModel>
449 (
450  const volScalarField& fc,
451  const volScalarField& fd
452 ) const
453 {
454  const volScalarField& alphag = this->alpha_;
455  const volScalarField& alphal = this->liquidTurbulence().alpha_;
456 
457  return
458  (alphal*rholEff()*fc + alphag*rhogEff()*Ct2_()*fd)
459  /(alphal*rholEff() + alphag*rhogEff()*Ct2_());
460 }
461 
462 
463 template<class BasicTurbulenceModel>
465 (
466  const surfaceScalarField& fc,
467  const surfaceScalarField& fd
468 ) const
469 {
470  const volScalarField& alphag = this->alpha_;
471  const volScalarField& alphal = this->liquidTurbulence().alpha_;
472 
473  surfaceScalarField alphalf(fvc::interpolate(alphal));
474  surfaceScalarField alphagf(fvc::interpolate(alphag));
475 
476  surfaceScalarField rholEfff(fvc::interpolate(rholEff()));
477  surfaceScalarField rhogEfff(fvc::interpolate(rhogEff()));
478 
479  return
480  (alphalf*rholEfff*fc + alphagf*rhogEfff*fvc::interpolate(Ct2_())*fd)
481  /(alphalf*rholEfff + alphagf*rhogEfff*fvc::interpolate(Ct2_()));
482 }
483 
484 
485 template<class BasicTurbulenceModel>
487 {
488  const mixtureKEpsilon<BasicTurbulenceModel>& liquidTurbulence =
489  this->liquidTurbulence();
490 
491  const transportModel& gas = this->transport();
492  const twoPhaseSystem& fluid = refCast<const twoPhaseSystem>(gas.fluid());
493  const transportModel& liquid = fluid.otherPhase(gas);
494 
495  const dragModel& drag = fluid.lookupSubModel<dragModel>(gas, liquid);
496 
497  volScalarField magUr(mag(liquidTurbulence.U() - this->U()));
498 
499  // Lahey model
500  tmp<volScalarField> bubbleG
501  (
502  Cp_
503  *liquid*liquid.rho()
504  *(
505  pow3(magUr)
506  + pow(drag.CdRe()*liquid.nu()/gas.d(), 4.0/3.0)
507  *pow(magUr, 5.0/3.0)
508  )
509  *gas
510  /gas.d()
511  );
512 
513  // Simple model
514  // tmp<volScalarField> bubbleG
515  // (
516  // Cp_*liquid*drag.K()*sqr(magUr)
517  // );
518 
519  return bubbleG;
520 }
521 
522 
523 template<class BasicTurbulenceModel>
525 {
526  return fvm::Su(bubbleG()/rhom_(), km_());
527 }
528 
529 
530 template<class BasicTurbulenceModel>
532 {
533  return fvm::Su(C3_*epsilonm_()*bubbleG()/(rhom_()*km_()), epsilonm_());
534 }
535 
536 
537 template<class BasicTurbulenceModel>
539 {
540  const transportModel& gas = this->transport();
541  const twoPhaseSystem& fluid = refCast<const twoPhaseSystem>(gas.fluid());
542 
543  // Only solve the mixture turbulence for the gas-phase
544  if (&gas != &fluid.phase1())
545  {
546  // This is the liquid phase but check the model for the gas-phase
547  // is consistent
548  this->liquidTurbulence();
549 
550  return;
551  }
552 
553  if (!this->turbulence_)
554  {
555  return;
556  }
557 
558  // Initialise the mixture fields if they have not yet been constructed
559  initMixtureFields();
560 
561  // Local references to gas-phase properties
562  tmp<surfaceScalarField> phig = this->phi();
563  const volVectorField& Ug = this->U_;
564  const volScalarField& alphag = this->alpha_;
565  volScalarField& kg = this->k_;
566  volScalarField& epsilong = this->epsilon_;
567  volScalarField& nutg = this->nut_;
568 
569  // Local references to liquid-phase properties
570  mixtureKEpsilon<BasicTurbulenceModel>& liquidTurbulence =
571  this->liquidTurbulence();
572  tmp<surfaceScalarField> phil = liquidTurbulence.phi();
573  const volVectorField& Ul = liquidTurbulence.U_;
574  const volScalarField& alphal = liquidTurbulence.alpha_;
575  volScalarField& kl = liquidTurbulence.k_;
576  volScalarField& epsilonl = liquidTurbulence.epsilon_;
577  volScalarField& nutl = liquidTurbulence.nut_;
578 
579  // Local references to mixture properties
580  volScalarField& rhom = rhom_();
581  volScalarField& km = km_();
582  volScalarField& epsilonm = epsilonm_();
583 
584  fv::options& fvOptions(fv::options::New(this->mesh_));
585 
587 
588  // Update the effective mixture density
589  rhom = this->rhom();
590 
591  // Mixture flux
592  surfaceScalarField phim("phim", mixFlux(phil, phig));
593 
594  // Mixture velocity divergence
595  volScalarField divUm
596  (
597  mixU
598  (
599  fvc::div(fvc::absolute(phil, Ul)),
600  fvc::div(fvc::absolute(phig, Ug))
601  )
602  );
603 
605  {
606  tmp<volTensorField> tgradUl = fvc::grad(Ul);
608  (
609  new volScalarField
610  (
611  this->GName(),
612  nutl*(tgradUl() && dev(twoSymm(tgradUl())))
613  )
614  );
615  tgradUl.clear();
616 
617  // Update k, epsilon and G at the wall
618  kl.boundaryFieldRef().updateCoeffs();
619  epsilonl.boundaryFieldRef().updateCoeffs();
620 
621  Gc.ref().checkOut();
622  }
623 
625  {
626  tmp<volTensorField> tgradUg = fvc::grad(Ug);
628  (
629  new volScalarField
630  (
631  this->GName(),
632  nutg*(tgradUg() && dev(twoSymm(tgradUg())))
633  )
634  );
635  tgradUg.clear();
636 
637  // Update k, epsilon and G at the wall
638  kg.boundaryFieldRef().updateCoeffs();
639  epsilong.boundaryFieldRef().updateCoeffs();
640 
641  Gd.ref().checkOut();
642  }
643 
644  // Mixture turbulence generation
645  volScalarField Gm(mix(Gc, Gd));
646 
647  // Mixture turbulence viscosity
648  volScalarField nutm(mixU(nutl, nutg));
649 
650  // Update the mixture k and epsilon boundary conditions
651  km == mix(kl, kg);
652  bound(km, this->kMin_);
653  epsilonm == mix(epsilonl, epsilong);
654  bound(epsilonm, this->epsilonMin_);
655 
656  // Dissipation equation
657  tmp<fvScalarMatrix> epsEqn
658  (
659  fvm::ddt(epsilonm)
660  + fvm::div(phim, epsilonm)
661  - fvm::Sp(fvc::div(phim), epsilonm)
662  - fvm::laplacian(DepsilonEff(nutm), epsilonm)
663  ==
664  C1_*Gm*epsilonm/km
665  - fvm::SuSp(((2.0/3.0)*C1_)*divUm, epsilonm)
666  - fvm::Sp(C2_*epsilonm/km, epsilonm)
667  + epsilonSource()
668  + fvOptions(epsilonm)
669  );
670 
671  epsEqn.ref().relax();
672  fvOptions.constrain(epsEqn.ref());
673  epsEqn.ref().boundaryManipulate(epsilonm.boundaryFieldRef());
674  solve(epsEqn);
675  fvOptions.correct(epsilonm);
676  bound(epsilonm, this->epsilonMin_);
677 
678 
679  // Turbulent kinetic energy equation
680  tmp<fvScalarMatrix> kmEqn
681  (
682  fvm::ddt(km)
683  + fvm::div(phim, km)
684  - fvm::Sp(fvc::div(phim), km)
685  - fvm::laplacian(DkEff(nutm), km)
686  ==
687  Gm
688  - fvm::SuSp((2.0/3.0)*divUm, km)
689  - fvm::Sp(epsilonm/km, km)
690  + kSource()
691  + fvOptions(km)
692  );
693 
694  kmEqn.ref().relax();
695  fvOptions.constrain(kmEqn.ref());
696  solve(kmEqn);
697  fvOptions.correct(km);
698  bound(km, this->kMin_);
700 
701  volScalarField Cc2(rhom/(alphal*rholEff() + alphag*rhogEff()*Ct2_()));
702  kl = Cc2*km;
704  epsilonl = Cc2*epsilonm;
705  epsilonl.correctBoundaryConditions();
706  liquidTurbulence.correctNut();
707 
708  Ct2_() = Ct2();
709  kg = Ct2_()*kl;
711  epsilong = Ct2_()*epsilonl;
712  epsilong.correctBoundaryConditions();
713  nutg = Ct2_()*(liquidTurbulence.nu()/this->nu())*nutl;
714 }
715 
716 
717 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
718 
719 } // End namespace RASModels
720 } // End namespace Foam
721 
722 // ************************************************************************* //
const phaseModel & phase1() const
Return phase model 1.
BasicTurbulenceModel::alphaField alphaField
const modelType & lookupSubModel(const phasePair &key) const
Access a sub model between a phase pair.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
tmp< surfaceScalarField > mixFlux(const surfaceScalarField &fc, const surfaceScalarField &fd) const
tmp< volScalarField > mixU(const volScalarField &fc, const volScalarField &fd) const
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
BasicTurbulenceModel::rhoField rhoField
void clear() const
If object pointer points to valid object:
Definition: tmpI.H:230
fv::options & fvOptions
surfaceScalarField & phi
multiphaseSystem & fluid
Definition: createFields.H:11
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
BasicTurbulenceModel::transportModel transportModel
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
tmp< fvMatrix< Type > > Su(const DimensionedField< Type, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
dimensionedSymmTensor sqr(const dimensionedVector &dv)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
tmp< volScalarField > mix(const volScalarField &fc, const volScalarField &fd) const
surfaceScalarField phig("phig", -rhorAUf *ghf *fvc::snGrad(rho) *mesh.magSf())
dimensionedScalar sqrt(const dimensionedScalar &ds)
wordList epsilonBoundaryTypes(const volScalarField &epsilon) const
Generic dimensioned Type class.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:49
Finite-volume options.
Definition: fvOptions.H:52
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Info<< "Reading strained laminar flame speed field Su\"<< endl;volScalarField Su(IOobject("Su", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading field betav\"<< endl;volScalarField betav(IOobject("betav", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Lobs\"<< endl;volScalarField Lobs(IOobject("Lobs", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field CT\"<< endl;volSymmTensorField CT(IOobject("CT", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Nv\"<< endl;volScalarField Nv(IOobject("Nv", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field nsv\"<< endl;volSymmTensorField nsv(IOobject("nsv", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);IOdictionary PDRProperties(IOobject("PDRProperties", runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED, IOobject::NO_WRITE));autoPtr< PDRDragModel > drag
Definition: createFields.H:182
tmp< volScalarField > rhom() const
mixtureKEpsilon(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName=turbulenceModel::propertiesName, const word &type=typeName)
Construct from components.
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
tmp< volScalarField > bubbleG() const
Class which solves the volume fraction equations for two phases.
tmp< volScalarField > rhogEff() const
Templated abstract base class for RAS turbulence models.
Definition: RASModel.H:49
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
bool read(const char *, int32_t &)
Definition: int32IO.C:85
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
dimensionedScalar exp(const dimensionedScalar &ds)
static const word propertiesName
Default name of the turbulence properties dictionary.
alphal
Definition: alphavPsi.H:12
A class for handling words, derived from string.
Definition: word.H:59
static word groupName(Name name, const word &group)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
virtual bool read()
Re-read model coefficients if they have changed.
virtual tmp< fvScalarMatrix > epsilonSource() const
void constrain(fvMatrix< Type > &eqn)
Apply constraints to equation.
virtual tmp< fvScalarMatrix > kSource() const
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
Bound the given scalar field if it has gone unbounded.
const phaseModel & otherPhase(const phaseModel &phase) const
Return the phase not given as an argument.
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-7/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H" 1{ alphav=max(min((rho - rholSat)/(rhovSat - rholSat), scalar(1)), scalar(0));alphal=1.0 - alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
Definition: pEqn.H:68
virtual tmp< volScalarField > Cvm() const =0
Return the virtual mass coefficient.
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< volScalarField > Ct2() const
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
label patchi
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
Definition: fvcMeshPhi.C:188
U
Definition: pEqn.H:72
const dimensionedScalar & rho() const
Definition: phaseModel.H:165
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
Definition: bound.C:33
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
Mixture k-epsilon turbulence model for two-phase gas-liquid systems.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:46
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
void correctBoundaryConditions()
Correct boundary field.
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
dimensioned< scalar > mag(const dimensioned< Type > &)
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
A class for managing temporary objects.
Definition: PtrList.H:53
virtual tmp< volScalarField > CdRe() const =0
Drag coefficient.
static options & New(const fvMesh &mesh)
Construct fvOptions and register to datbase if not present.
Definition: fvOptions.C:101
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:354
void correctInletOutlet(volScalarField &vsf, const volScalarField &refVsf) const
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
tmp< volScalarField > rholEff() const
volScalarField & nu
Calculate the matrix for implicit and explicit sources.
Namespace for OpenFOAM.