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-2023 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 "fvModels.H"
28 #include "bound.H"
29 #include "phaseSystem.H"
30 #include "dispersedDragModel.H"
36 #include "fvmSup.H"
37 
38 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
39 
40 namespace Foam
41 {
42 namespace RASModels
43 {
44 
45 
46 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
47 
48 template<class BasicMomentumTransportModel>
50 (
51  const volScalarField& Cc2
52 )
53 {
54  epsilonm_() = max
55  (
56  epsilonm_(),
57  Cmu_*sqr(k_)/(this->nutMaxCoeff_*Cc2*this->nu())
58  );
59 }
60 
61 
62 template<class BasicMomentumTransportModel>
64 {
65  this->nut_ = Cmu_*sqr(k_)/epsilon_;
66  this->nut_.correctBoundaryConditions();
67  fvConstraints::New(this->mesh_).constrain(this->nut_);
68 }
69 
70 
71 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
72 
73 template<class BasicMomentumTransportModel>
75 (
76  const alphaField& alpha,
77  const rhoField& rho,
78  const volVectorField& U,
79  const surfaceScalarField& alphaRhoPhi,
80  const surfaceScalarField& phi,
81  const viscosity& viscosity,
82  const word& type
83 )
84 :
85  eddyViscosity<RASModel<BasicMomentumTransportModel>>
86  (
87  type,
88  alpha,
89  rho,
90  U,
91  alphaRhoPhi,
92  phi,
93  viscosity
94  ),
95 
96  gasTurbulencePtr_(nullptr),
97 
98  Cmu_
99  (
100  dimensioned<scalar>::lookupOrAddToDict
101  (
102  "Cmu",
103  this->coeffDict_,
104  0.09
105  )
106  ),
107  C1_
108  (
109  dimensioned<scalar>::lookupOrAddToDict
110  (
111  "C1",
112  this->coeffDict_,
113  1.44
114  )
115  ),
116  C2_
117  (
118  dimensioned<scalar>::lookupOrAddToDict
119  (
120  "C2",
121  this->coeffDict_,
122  1.92
123  )
124  ),
125  C3_
126  (
127  dimensioned<scalar>::lookupOrAddToDict
128  (
129  "C3",
130  this->coeffDict_,
131  C2_.value()
132  )
133  ),
134  Cp_
135  (
136  dimensioned<scalar>::lookupOrAddToDict
137  (
138  "Cp",
139  this->coeffDict_,
140  0.25
141  )
142  ),
143  alphap_
144  (
145  dimensioned<scalar>::lookupOrAddToDict
146  (
147  "alphap",
148  this->coeffDict_,
149  1
150  )
151  ),
152  sigmak_
153  (
154  dimensioned<scalar>::lookupOrAddToDict
155  (
156  "sigmak",
157  this->coeffDict_,
158  1.0
159  )
160  ),
161  sigmaEps_
162  (
163  dimensioned<scalar>::lookupOrAddToDict
164  (
165  "sigmaEps",
166  this->coeffDict_,
167  1.3
168  )
169  ),
170 
171  k_
172  (
173  IOobject
174  (
175  this->groupName("k"),
176  this->runTime_.name(),
177  this->mesh_,
178  IOobject::MUST_READ,
179  IOobject::AUTO_WRITE
180  ),
181  this->mesh_
182  ),
183  epsilon_
184  (
185  IOobject
186  (
187  this->groupName("epsilon"),
188  this->runTime_.name(),
189  this->mesh_,
190  IOobject::MUST_READ,
191  IOobject::AUTO_WRITE
192  ),
193  this->mesh_
194  )
195 {
196  bound(k_, this->kMin_);
197 
198  if (type == typeName)
199  {
200  this->printCoeffs(type);
201  }
202 
203  const phaseModel& phase = refCast<const phaseModel>(this->properties());
204 
205  // Construct mixture properties only for liquid phase (phase 1)
206  if (phase.index() == 1)
207  {
208  km_.set
209  (
210  new volScalarField
211  (
212  IOobject
213  (
214  "km",
215  this->runTime_.name(),
216  this->mesh_,
219  ),
220  this->mesh_
221  )
222  );
223 
224  epsilonm_.set
225  (
226  new volScalarField
227  (
228  IOobject
229  (
230  "epsilonm",
231  this->runTime_.name(),
232  this->mesh_,
235  ),
236  this->mesh_
237  )
238  );
239 
240  Ct2_.set
241  (
242  new volScalarField
243  (
244  IOobject
245  (
246  "Ct2",
247  this->runTime_.name(),
248  this->mesh_,
251  ),
252  this->mesh_,
254  )
255  );
256 
257  rhom_.set
258  (
259  new volScalarField
260  (
261  IOobject
262  (
263  "rhom",
264  this->runTime_.name(),
265  this->mesh_,
268  ),
269  this->mesh_,
271  )
272  );
273  }
274 }
275 
276 
277 template<class BasicMomentumTransportModel>
279 {
280  static bool initialised = false;
281 
282  if (!initialised)
283  {
284  Ct2_() = Ct2();
285  rhom_() = rhom();
286  initialised = true;
287  }
288 }
289 
290 
291 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
292 
293 template<class BasicMomentumTransportModel>
295 {
297  {
298  Cmu_.readIfPresent(this->coeffDict());
299  C1_.readIfPresent(this->coeffDict());
300  C2_.readIfPresent(this->coeffDict());
301  C3_.readIfPresent(this->coeffDict());
302  Cp_.readIfPresent(this->coeffDict());
303  sigmak_.readIfPresent(this->coeffDict());
304  sigmaEps_.readIfPresent(this->coeffDict());
305 
306  return true;
307  }
308  else
309  {
310  return false;
311  }
312 }
313 
314 
315 template<class BasicMomentumTransportModel>
318 {
319  if (!gasTurbulencePtr_)
320  {
321  const volVectorField& U = this->U_;
322 
323  const phaseModel& liquid =
324  refCast<const phaseModel>(this->properties());
325  const phaseSystem& fluid = liquid.fluid();
326  const phaseModel& gas = fluid.phases()[0];
327 
328  gasTurbulencePtr_ =
330  (
331  U.db().lookupObject
332  <
334  >
335  (
337  (
338  momentumTransportModel::typeName,
339  gas.name()
340  )
341  )
342  );
343  }
344 
345  return *gasTurbulencePtr_;
346 }
347 
348 
349 template<class BasicMomentumTransportModel>
351 {
352  const mixtureKEpsilon<BasicMomentumTransportModel>& gasTurbulence =
353  this->gasTurbulence();
354 
355  const phaseModel& liquid = refCast<const phaseModel>(this->properties());
356  const phaseSystem& fluid = liquid.fluid();
357  const phaseModel& gas = fluid.phases()[0];
358 
362 
363  const volScalarField& alphag = gasTurbulence.alpha();
364  const volScalarField magUr(mag(this->U() - gasTurbulence.U()));
365 
366  volScalarField beta
367  (
368  (6*this->Cmu_/(4*sqrt(3.0/2.0)))
369  *drag.K()/liquid.rho()
370  *(k_/epsilon_)
371  );
372  volScalarField Ct0((3 + beta)/(1 + beta + 2*gas.rho()/liquid.rho()));
373  volScalarField fAlphad((180 + (-4.71e3 + 4.26e4*alphag)*alphag)*alphag);
374 
375  return sqr(1 + (Ct0 - 1)*exp(-fAlphad));
376 }
377 
378 
379 template<class BasicMomentumTransportModel>
382 {
383  const phaseModel& liquid = refCast<const phaseModel>(this->properties());
384  return liquid.rho();
385 }
386 
387 
388 template<class BasicMomentumTransportModel>
391 {
392  const phaseModel& liquid = refCast<const phaseModel>(this->properties());
393  const phaseSystem& fluid = liquid.fluid();
394  const phaseModel& gas = fluid.phases()[0];
395 
400 
401  return gas.rho() + virtualMass.Cvm()*liquid.rho();
402 }
403 
404 
405 template<class BasicMomentumTransportModel>
407 {
408  const volScalarField& alphal = this->alpha_;
409  const volScalarField& alphag = this->gasTurbulence().alpha_;
410 
411  return alphal*rholEff() + alphag*rhogEff();
412 }
413 
414 
415 template<class BasicMomentumTransportModel>
417 (
418  const volScalarField& fc,
419  const volScalarField& fd
420 ) const
421 {
422  const volScalarField& alphal = this->alpha_;
423  const volScalarField& alphag = this->gasTurbulence().alpha_;
424 
425  return (alphal*rholEff()*fc + alphag*rhogEff()*fd)/rhom_();
426 }
427 
428 
429 template<class BasicMomentumTransportModel>
431 (
432  const volScalarField& fc,
433  const volScalarField& fd
434 ) const
435 {
436  const volScalarField& alphal = this->alpha_;
437  const volScalarField& alphag = this->gasTurbulence().alpha_;
438 
439  return
440  (alphal*rholEff()*fc + alphag*rhogEff()*Ct2_()*fd)
441  /(alphal*rholEff() + alphag*rhogEff()*Ct2_());
442 }
443 
444 
445 template<class BasicMomentumTransportModel>
447 (
448  const surfaceScalarField& fc,
449  const surfaceScalarField& fd
450 ) const
451 {
452  const volScalarField& alphal = this->alpha_;
453  const volScalarField& alphag = this->gasTurbulence().alpha_;
454 
455  surfaceScalarField alphalf(fvc::interpolate(alphal));
456  surfaceScalarField alphagf(fvc::interpolate(alphag));
457 
458  surfaceScalarField rholEfff(fvc::interpolate(rholEff()));
459  surfaceScalarField rhogEfff(fvc::interpolate(rhogEff()));
460 
461  return
462  (alphalf*rholEfff*fc + alphagf*rhogEfff*fvc::interpolate(Ct2_())*fd)
463  /(alphalf*rholEfff + alphagf*rhogEfff*fvc::interpolate(Ct2_()));
464 }
465 
466 
467 template<class BasicMomentumTransportModel>
470 {
471  const mixtureKEpsilon<BasicMomentumTransportModel>& gasTurbulence =
472  this->gasTurbulence();
473 
474  const phaseModel& liquid = refCast<const phaseModel>(this->properties());
475  const phaseSystem& fluid = liquid.fluid();
476  const phaseModel& gas = fluid.phases()[0];
477 
481 
482  volScalarField magUr(mag(this->U() - gasTurbulence.U()));
483 
484  // Lahey model
485  tmp<volScalarField> bubbleG
486  (
487  Cp_
488  *pos(alphap_ - gas)*liquid*liquid.rho()
489  *(
490  pow3(magUr)
491  + pow(drag.CdRe()*liquid.fluidThermo().nu()/gas.d(), 4.0/3.0)
492  *pow(magUr, 5.0/3.0)
493  )
494  *gas
495  /gas.d()
496  );
497 
498  // Simple model
499  // tmp<volScalarField> bubbleG
500  // (
501  // Cp_*liquid*drag.K()*sqr(magUr)
502  // );
503 
504  return bubbleG;
505 }
506 
507 
508 template<class BasicMomentumTransportModel>
511 {
512  return fvm::Su(bubbleG()/rhom_(), km_());
513 }
514 
515 
516 template<class BasicMomentumTransportModel>
519 {
520  return fvm::Su(C3_*epsilonm_()*bubbleG()/(rhom_()*km_()), epsilonm_());
521 }
522 
523 
524 template<class BasicMomentumTransportModel>
526 {
527  const phaseModel& phase = refCast<const phaseModel>(this->properties());
528 
529  // Only solve the mixture turbulence for the liquid phase (phase 1)
530  if (phase.index() == 0)
531  {
532  // This is the liquid phase but check the model for the gas-phase
533  // is consistent
534  this->gasTurbulence();
535 
536  return;
537  }
538  else
539  {
540  initMixtureFields();
541  }
542 
543  if (!this->turbulence_)
544  {
545  return;
546  }
547 
548  // Local references to liquid-phase properties
549  tmp<surfaceScalarField> phil = this->phi();
550  const volVectorField& Ul = this->U_;
551  const volScalarField& alphal = this->alpha_;
552  volScalarField& kl = this->k_;
553  volScalarField& epsilonl = this->epsilon_;
554  volScalarField& nutl = this->nut_;
555 
556  // Local references to gas-phase properties
558  this->gasTurbulence();
559  tmp<surfaceScalarField> phig = gasTurbulence.phi();
560  const volVectorField& Ug = gasTurbulence.U_;
561  const volScalarField& alphag = gasTurbulence.alpha_;
562  volScalarField& kg = gasTurbulence.k_;
563  volScalarField& epsilong = gasTurbulence.epsilon_;
564  volScalarField& nutg = gasTurbulence.nut_;
565 
566  // Local references to mixture properties
567  volScalarField& rhom = rhom_();
568  volScalarField& km = km_();
569  volScalarField& epsilonm = epsilonm_();
570 
571  const Foam::fvModels& fvModels(Foam::fvModels::New(this->mesh_));
573  (
574  Foam::fvConstraints::New(this->mesh_)
575  );
576 
578 
579  // Update the effective mixture density
580  rhom = this->rhom();
581 
582  // Mixture flux
583  const surfaceScalarField phim("phim", mixFlux(phil, phig));
584 
585  // Mixture velocity divergence
586  const volScalarField divUm
587  (
588  mixU
589  (
590  fvc::div(fvc::absolute(phil, Ul)),
591  fvc::div(fvc::absolute(phig, Ug))
592  )
593  );
594 
596  {
597  tmp<volTensorField> tgradUl = fvc::grad(Ul);
599  (
600  new volScalarField
601  (
602  this->GName(),
603  nutl*(tgradUl() && dev(twoSymm(tgradUl())))
604  )
605  );
606  tgradUl.clear();
607 
608  // Update k, epsilon and G at the wall
610  epsilonl.boundaryFieldRef().updateCoeffs();
611 
612  Gc.ref().checkOut();
613  }
614 
616  {
617  tmp<volTensorField> tgradUg = fvc::grad(Ug);
619  (
620  new volScalarField
621  (
622  this->GName(),
623  nutg*(tgradUg() && dev(twoSymm(tgradUg())))
624  )
625  );
626  tgradUg.clear();
627 
628  // Update k, epsilon and G at the wall
630  epsilong.boundaryFieldRef().updateCoeffs();
631 
632  Gd.ref().checkOut();
633  }
634 
635  // Mixture turbulence generation
636  const volScalarField Gm(mix(Gc, Gd));
637 
638  // Mixture turbulence viscosity
639  const volScalarField nutm(mixU(nutl, nutg));
640 
641  // Update the mixture k and epsilon boundary conditions
642  km == mix(kl, kg);
643  epsilonm == mix(epsilonl, epsilong);
644 
645  // Dissipation equation
646  tmp<fvScalarMatrix> epsEqn
647  (
648  fvm::ddt(epsilonm)
649  + fvm::div(phim, epsilonm)
650  + fvm::SuSp(-fvc::div(phim), epsilonm)
651  - fvm::laplacian(DepsilonEff(nutm), epsilonm)
652  ==
653  C1_*Gm*epsilonm/km
654  - fvm::SuSp(((2.0/3.0)*C1_)*divUm, epsilonm)
655  - fvm::Sp(C2_*epsilonm/km, epsilonm)
656  + epsilonSource()
657  + fvModels.source(epsilonm)
658  );
659 
660  epsEqn.ref().relax();
661  fvConstraints.constrain(epsEqn.ref());
662  epsEqn.ref().boundaryManipulate(epsilonm.boundaryFieldRef());
663  solve(epsEqn);
664  fvConstraints.constrain(epsilonm);
665 
666  const volScalarField Cc2(rhom/(alphal*rholEff() + alphag*rhogEff()*Ct2_()));
667  boundEpsilonm(Cc2);
668 
669  // Turbulent kinetic energy equation
670  tmp<fvScalarMatrix> kmEqn
671  (
672  fvm::ddt(km)
673  + fvm::div(phim, km)
674  + fvm::SuSp(-fvc::div(phim), km)
675  - fvm::laplacian(DkEff(nutm), km)
676  ==
677  Gm
678  - fvm::SuSp((2.0/3.0)*divUm, km)
679  - fvm::Sp(epsilonm/km, km)
680  + kSource()
681  + fvModels.source(km)
682  );
683 
684  kmEqn.ref().relax();
685  fvConstraints.constrain(kmEqn.ref());
686  solve(kmEqn);
688  bound(km, this->kMin_);
690  boundEpsilonm(Cc2);
691 
692  kl = Cc2*km;
694  epsilonl = Cc2*epsilonm;
695  epsilonl.correctBoundaryConditions();
696  correctNut();
697 
698  Ct2_() = Ct2();
699  kg = Ct2_()*kl;
701  epsilong = Ct2_()*epsilonl;
702  epsilong.correctBoundaryConditions();
703  nutg = Ct2_()*(this->nu()/gasTurbulence.nu())*nutl;
704 }
705 
706 
707 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
708 
709 } // End namespace RASModels
710 } // End namespace Foam
711 
712 // ************************************************************************* //
Bound the given scalar field where it is below the specified minimum.
void updateCoeffs()
Update the boundary condition coefficients.
Generic GeometricField class.
void correctBoundaryConditions()
Correct boundary field.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
static word groupName(Name name, const word &group)
Templated abstract base class for RAS turbulence models.
Definition: RASModel.H:56
Mixture k-epsilon turbulence model for two-phase gas-liquid systems.
virtual tmp< fvScalarMatrix > epsilonSource() const
tmp< volScalarField > Ct2() const
tmp< volScalarField > rhogEff() const
tmp< volScalarField > mixU(const volScalarField &fc, const volScalarField &fd) const
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
tmp< volScalarField > rhom() const
tmp< volScalarField > bubbleG() const
autoPtr< volScalarField > rhom_
autoPtr< volScalarField > km_
autoPtr< volScalarField > Ct2_
virtual void correctNut()
Correct the eddy-viscosity nut.
virtual tmp< fvScalarMatrix > kSource() const
tmp< surfaceScalarField > mixFlux(const surfaceScalarField &fc, const surfaceScalarField &fd) const
tmp< volScalarField > rholEff() const
autoPtr< volScalarField > epsilonm_
tmp< volScalarField > mix(const volScalarField &fc, const volScalarField &fd) const
void boundEpsilonm(const volScalarField &Cc2)
Bound epsilonm.
virtual bool read()
Re-read model coefficients if they have changed.
mixtureKEpsilon(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const viscosity &viscosity, const word &type=typeName)
Construct from components.
static autoPtr< dictionary > New(Istream &)
Construct top-level dictionary on freestore from Istream.
Definition: dictionaryIO.C:100
Generic dimensioned Type class.
Class to represent a interface between phases where one phase is considered dispersed within the othe...
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:52
volScalarField nut_
Definition: eddyViscosity.H:60
Finite volume constraints.
Definition: fvConstraints.H:67
bool constrain(fvMatrix< Type > &eqn) const
Apply constraints to an equation.
Finite volume models.
Definition: fvModels.H:65
tmp< fvMatrix< Type > > source(const VolField< Type > &field) const
Return source for an equation.
BasicMomentumTransportModel::alphaField alphaField
BasicMomentumTransportModel::rhoField rhoField
Generic thermophysical properties class for a liquid in which the functions and coefficients for each...
Definition: liquid.H:53
scalar rho(scalar p, scalar T) const
Liquid density [kg/m^3].
Definition: liquidI.H:26
tmp< volScalarField > d() const
Return the Sauter-mean diameter.
Definition: phaseModel.C:181
label index() const
Return the index of the phase.
Definition: phaseModel.C:157
const word & name() const
Return the name of this phase.
Definition: phaseModel.C:145
virtual const volScalarField & rho() const =0
Return the density field.
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:73
const phaseModelList & phases() const
Return the phase models.
Definition: phaseSystemI.H:102
const ModelType & lookupInterfacialModel(const phaseInterface &interface) const
Return a sub model for an interface.
A class for managing temporary objects.
Definition: tmp.H:55
void clear() const
If object pointer points to valid object:
Definition: tmpI.H:237
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
virtual tmp< volScalarField > Cvm() const =0
Return the virtual mass coefficient.
Abstract base class for all fluid physical properties.
Definition: viscosity.H:50
A class for handling words, derived from string.
Definition: word.H:62
Foam::fvConstraints & fvConstraints(Foam::fvConstraints::New(mesh))
Foam::fvModels & fvModels(Foam::fvModels::New(mesh))
Calculate the matrix for implicit and explicit sources.
Info<< "Reading strained laminar flame speed field Su\n"<< endl;volScalarField Su(IOobject("Su", runTime.name(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading field betav\n"<< endl;volScalarField betav(IOobject("betav", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Lobs\n"<< endl;volScalarField Lobs(IOobject("Lobs", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field CT\n"<< endl;volSymmTensorField CT(IOobject("CT", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Nv\n"<< endl;volScalarField Nv(IOobject("Nv", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field nsv\n"<< 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:192
U
Definition: pEqn.H:72
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
void correct(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiCorr, const SpType &Sp, const SuType &Su)
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
Definition: fvcGrad.C:46
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
Definition: fvcDiv.C:47
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
Definition: fvcMeshPhi.C:202
tmp< fvMatrix< Type > > Su(const DimensionedField< Type, volMesh > &, const VolField< Type > &)
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
Definition: fvmLaplacian.C:47
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const VolField< Type > &vf, const word &name)
Definition: fvmDiv.C:48
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const VolField< Type > &)
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const VolField< Type > &)
tmp< fvMatrix< Type > > ddt(const VolField< Type > &vf)
Definition: fvmDdt.C:46
Namespace for OpenFOAM.
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
dimensionedScalar pos(const dimensionedScalar &ds)
dimensionedScalar exp(const dimensionedScalar &ds)
bool read(const char *, int32_t &)
Definition: int32IO.C:85
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow3(const dimensionedScalar &ds)
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
const dimensionSet dimless
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
const dimensionSet dimDensity
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
bool bound(volScalarField &, const dimensionedScalar &min)
Bound the given scalar field where it is below the specified min value.
Definition: bound.C:31
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.