IshiiZuber.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) 2014-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 "IshiiZuber.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace dragModels
34 {
37 }
38 }
39 
40 
41 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
42 
44 (
45  const dictionary& dict,
46  const phaseInterface& interface,
47  const bool registerObject
48 )
49 :
50  dispersedDragModel(dict, interface, registerObject)
51 {}
52 
53 
54 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
55 
57 {}
58 
59 
60 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
61 
64 {
65  const volScalarField Re(interface_.Re());
66  const volScalarField Eo(interface_.Eo());
67 
68  const volScalarField mud(interface_.dispersed().fluidThermo().mu());
69  const volScalarField muc(interface_.continuous().fluidThermo().mu());
70 
71  const volScalarField muStar((mud + 0.4*muc)/(mud + muc));
72 
73  const volScalarField muMix
74  (
75  muc*pow(max(1 - interface_.dispersed(), scalar(1e-3)), -2.5*muStar)
76  );
77 
78  const volScalarField ReM(Re*muc/muMix);
79  const volScalarField CdRe
80  (
81  pos0(1000 - ReM)*24*(1 + 0.1*pow(ReM, 0.75))
82  + neg(1000 - ReM)*0.44*ReM
83  );
84 
85  volScalarField F((muc/muMix)*sqrt(1 - interface_.dispersed()));
86  F.max(1e-3);
87 
88  const volScalarField Ealpha((1 + 17.67*pow(F, 0.8571428))/(18.67*F));
89 
90  const volScalarField CdReEllipse(Ealpha*0.6666*sqrt(Eo)*Re);
91 
92  return
93  pos0(CdReEllipse - CdRe)
94  *min(CdReEllipse, Re*sqr(1 - interface_.dispersed())*2.66667)
95  + neg(CdReEllipse - CdRe)*CdRe;
96 }
97 
98 
99 // ************************************************************************* //
Macros for easy insertion into run-time selection tables.
Generic GeometricField class.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Model for drag between phases.
Definition: dragModel.H:55
Ishii and Zuber (1979) drag model for dense dispersed bubbly flows.
Definition: IshiiZuber.H:61
IshiiZuber(const dictionary &dict, const phaseInterface &interface, const bool registerObject)
Construct from a dictionary and an interface.
Definition: IshiiZuber.C:44
virtual ~IshiiZuber()
Destructor.
Definition: IshiiZuber.C:56
virtual tmp< volScalarField > CdRe() const
Drag coefficient.
Definition: IshiiZuber.C:63
Class to represent an interface between phases. Derivations can further specify the configuration of ...
A class for managing temporary objects.
Definition: tmp.H:55
const dimensionedScalar F
Faraday constant: default SI units: [C/kmol].
defineTypeNameAndDebug(aerosolDrag, 0)
addToRunTimeSelectionTable(dragModel, aerosolDrag, dictionary)
Namespace for OpenFOAM.
const doubleScalar e
Definition: doubleScalar.H:106
dimensionedScalar pos0(const dimensionedScalar &ds)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar neg(const dimensionedScalar &ds)
scalarField Re(const UList< complex > &cf)
Definition: complexFields.C:97
dictionary dict