KocamustafaogullariIshii.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2016-2019 OpenFOAM Foundation
9 Copyright (C) 2018-2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
33#include "ThermalDiffusivity.H"
35#include "phaseSystem.H"
36
37// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
38
39namespace Foam
40{
41namespace wallBoilingModels
42{
43namespace departureDiameterModels
44{
47 (
51 );
52}
53}
54}
55
56
57// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
58
59Foam::wallBoilingModels::departureDiameterModels::
60KocamustafaogullariIshii::KocamustafaogullariIshii
61(
62 const dictionary& dict
63)
64:
66 phi_(dict.get<scalar>("phi"))
67{}
68
69
70// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
71
75(
76 const phaseModel& liquid,
77 const phaseModel& vapor,
78 const label patchi,
79 const scalarField& Tl,
80 const scalarField& Tsatw,
81 const scalarField& L
82) const
83{
84 // Gravitational acceleration
85 const auto& g =
86 liquid.mesh().time().lookupObject<uniformDimensionedVectorField>("g");
87
88 const scalarField rhoLiquid(liquid.thermo().rho(patchi));
89 const scalarField rhoVapor(vapor.thermo().rho(patchi));
90
91 const scalarField rhoM((rhoLiquid - rhoVapor)/rhoVapor);
92
93 const tmp<volScalarField>& tsigma =
94 liquid.fluid().sigma(phasePairKey(liquid.name(), vapor.name()));
95
96 const volScalarField& sigma = tsigma();
97 const fvPatchScalarField& sigmaw = sigma.boundaryField()[patchi];
98
99 return
100 0.0012*pow(rhoM, 0.9)*0.0208*phi_
101 *sqrt(sigmaw/(mag(g.value())*(rhoLiquid - rhoVapor)));
102}
103
104
107{
109 os.writeEntry("phi", phi_);
110}
111
112
113// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
const uniformDimensionedVectorField & g
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
virtual bool write()
Write the output fields.
Generic thermophysical properties class for a liquid in which the functions and coefficients for each...
Definition: liquid.H:57
scalar rho(scalar p, scalar T) const
Liquid density [kg/m^3].
Definition: liquidI.H:28
scalar sigma(scalar p, scalar T) const
Surface tension [N/m].
Definition: liquidI.H:94
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:61
const word & name() const
Definition: phaseModel.H:144
virtual const rhoThermo & thermo() const =0
Return the thermophysical model.
An ordered or unorder pair of phase names. Typically specified as follows.
Definition: phasePairKey.H:68
virtual tmp< volScalarField > rho() const
Density [kg/m^3].
Definition: rhoThermo.C:224
A class for managing temporary objects.
Definition: tmp.H:65
Base class for bubble departure diameter models for boiling flows.
A correlation for bubble departure diameter modelling based on Kocamustafaogullari-Ishii (1983) for b...
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
OBJstream os(runTime.globalPath()/outputName)
Namespace for OpenFOAM.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
dictionary dict
const vector L(dict.get< vector >("L"))