TomiyamaAnalytic.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) 2014-2018 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28#include "TomiyamaAnalytic.H"
29#include "phasePair.H"
31
32// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33
34namespace Foam
35{
36namespace dragModels
37{
40}
41}
42
43
44// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
45
47(
48 const dictionary& dict,
49 const phasePair& pair,
50 const bool registerObject
51)
52:
53 dragModel(dict, pair, registerObject),
54 residualRe_("residualRe", dimless, dict),
55 residualEo_("residualEo", dimless, dict),
56 residualE_("residualE", dimless, dict)
57{}
58
59
60// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
61
63{}
64
65
66// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
67
70{
71 volScalarField Eo(max(pair_.Eo(), residualEo_));
72 volScalarField E(max(pair_.E(), residualE_));
73
74 volScalarField OmEsq(max(scalar(1) - sqr(E), sqr(residualE_)));
75 volScalarField rtOmEsq(sqrt(OmEsq));
76
77 volScalarField F(max(asin(rtOmEsq) - E*rtOmEsq, residualE_)/OmEsq);
78
79 return
80 (8.0/3.0)
81 *Eo
82 /(
83 Eo*pow(E, 2.0/3.0)/OmEsq
84 + 16*pow(E, 4.0/3.0)
85 )
86 /sqr(F)
87 *max(pair_.Re(), residualRe_);
88}
89
90
91// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Analytical drag model of Tomiyama et al.
virtual ~TomiyamaAnalytic()
Destructor.
virtual tmp< volScalarField > CdRe() const
Drag coefficient.
Description for mass transfer between a pair of phases. The direction of the mass transfer is from th...
Definition: phasePair.H:56
A class for managing temporary objects.
Definition: tmp.H:65
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
volVectorField F(fluid.F())
Namespace for OpenFOAM.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
dimensionedScalar asin(const dimensionedScalar &ds)
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dictionary dict