IATEsource.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-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 "IATEsource.H"
27 #include "twoPhaseSystem.H"
28 #include "fvMatrix.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 namespace diameterModels
37 {
38  defineTypeNameAndDebug(IATEsource, 0);
39  defineRunTimeSelectionTable(IATEsource, dictionary);
40 }
41 }
42 
43 
44 // * * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
45 
48 (
49  const word& type,
50  const IATE& iate,
51  const dictionary& dict
52 )
53 {
54  dictionaryConstructorTable::iterator cstrIter =
55  dictionaryConstructorTablePtr_->find(type);
56 
57  if (cstrIter == dictionaryConstructorTablePtr_->end())
58  {
60  << "Unknown IATE source type "
61  << type << nl << nl
62  << "Valid IATE source types : " << endl
63  << dictionaryConstructorTablePtr_->sortedToc()
64  << exit(FatalError);
65  }
66 
67  return autoPtr<IATEsource>(cstrIter()(iate, dict));
68 }
69 
70 
71 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
72 
74 {
76  phase().U().db().lookupObject<uniformDimensionedVectorField>("g");
77 
78  return
79  sqrt(2.0)
80  *pow025
81  (
82  fluid().sigma()*mag(g)
83  *(otherPhase().rho() - phase().rho())
84  /sqr(otherPhase().rho())
85  )
86  *pow(max(1 - phase(), scalar(0)), 1.75);
87 }
88 
90 {
91  return sqrt(2*otherPhase().turbulence().k());
92 }
93 
95 {
96  return max(Ur()*phase().d()/otherPhase().nu(), scalar(1.0e-3));
97 }
98 
100 {
101  const volScalarField Eo(this->Eo());
102  const volScalarField Re(this->Re());
103 
104  return
105  max
106  (
107  min
108  (
109  (16/Re)*(1 + 0.15*pow(Re, 0.687)),
110  48/Re
111  ),
112  8*Eo/(3*(Eo + 4))
113  );
114 }
115 
117 {
119  phase().U().db().lookupObject<uniformDimensionedVectorField>("g");
120 
121  return
122  mag(g)*pow4(otherPhase().nu())*sqr(otherPhase().rho())
123  *(otherPhase().rho() - phase().rho())
124  /pow3(fluid().sigma());
125 }
126 
128 {
130  phase().U().db().lookupObject<uniformDimensionedVectorField>("g");
131 
132  return
133  mag(g)*sqr(phase().d())
134  *(otherPhase().rho() - phase().rho())
135  /fluid().sigma();
136 }
137 
139 {
140  return otherPhase().rho()*sqr(Ur())*phase().d()/fluid().sigma();
141 }
142 
143 
144 // ************************************************************************* //
static autoPtr< IATEsource > New(const word &type, const IATE &iate, const dictionary &dict)
dictionary dict
tmp< volScalarField > Ut() const
Return the bubble turbulent velocity.
multiphaseSystem & fluid
Definition: createFields.H:11
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
UniformDimensionedField< vector > uniformDimensionedVectorField
error FatalError
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
dimensionedScalar pow025(const dimensionedScalar &ds)
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
tmp< volScalarField > We() const
Return the bubble Webber number.
tmp< volScalarField > Eo() const
Return the bubble Eotvos number.
label k
Boltzmann constant.
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
tmp< volScalarField > CD() const
Return the bubble drag coefficient.
tmp< volScalarField > Ur() const
Return the bubble relative velocity.
static const char nl
Definition: Ostream.H:265
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
tmp< volScalarField > Re() const
Return the bubble Reynolds number.
Info<< "Reading field U\"<< endl;volVectorField U(IOobject("U", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);volScalarField rho(IOobject("rho", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), thermo.rho());volVectorField rhoU(IOobject("rhoU", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *U);volScalarField rhoE(IOobject("rhoE", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *(e+0.5 *magSqr(U)));surfaceScalarField pos(IOobject("pos", runTime.timeName(), mesh), mesh, dimensionedScalar("pos", dimless, 1.0));surfaceScalarField neg(IOobject("neg", runTime.timeName(), mesh), mesh, dimensionedScalar("neg", dimless, -1.0));surfaceScalarField phi("phi", fvc::flux(rhoU));Info<< "Creating turbulence model\"<< endl;autoPtr< compressible::turbulenceModel > turbulence(compressible::turbulenceModel::New(rho, U, phi, thermo))
Definition: createFields.H:94
defineTypeNameAndDebug(combustionModel, 0)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
tmp< volScalarField > Mo() const
Return the bubble Morton number.
dimensionedScalar pow4(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:98
A class for managing temporary objects.
Definition: PtrList.H:53
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:361
scalarField Re(const UList< complex > &cf)
Definition: complexFields.C:97
volScalarField & nu
Namespace for OpenFOAM.