Moraga.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 -------------------------------------------------------------------------------
10 License
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 "Moraga.H"
29 #include "phasePair.H"
30 #include "fvcGrad.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 namespace liftModels
38 {
39  defineTypeNameAndDebug(Moraga, 0);
40  addToRunTimeSelectionTable(liftModel, Moraga, dictionary);
41 }
42 }
43 
44 
45 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
46 
48 (
49  const dictionary& dict,
50  const phasePair& pair
51 )
52 :
53  liftModel(dict, pair)
54 {}
55 
56 
57 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
58 
60 {}
61 
62 
63 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
64 
66 {
67  volScalarField Re(pair_.Re());
68 
69  volScalarField sqrSr
70  (
71  sqr(pair_.dispersed().d())
72  /pair_.continuous().nu()
73  *mag(fvc::grad(pair_.continuous().U()))
74  );
75 
76  if
77  (
78  min(Re).value() < 1200.0
79  || max(Re).value() > 18800.0
80  || min(sqrSr).value() < 0.0016
81  || max(sqrSr).value() > 0.04
82  )
83  {
85  << "Re and/or Sr are out of the range of applicability of the "
86  << "Moraga model. Clamping to range bounds"
87  << endl;
88  }
89 
90  Re.min(1200.0);
91  Re.max(18800.0);
92 
93  sqrSr.min(0.0016);
94  sqrSr.max(0.04);
95 
96  return 0.2*exp(- Re*sqrSr/3.6e5 - 0.12)*exp(Re*sqrSr/3.0e7);
97 }
98 
99 
100 // ************************************************************************* //
Foam::phasePair
Description for mass transfer between a pair of phases. The direction of the mass transfer is from th...
Definition: phasePair.H:53
Foam::fvc::grad
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:54
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::exp
dimensionedScalar exp(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:261
Foam::min
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
Moraga.H
Foam::liftModels::Moraga::Cl
virtual tmp< volScalarField > Cl() const
Lift coefficient.
Definition: Moraga.C:65
Foam::liftModels::Moraga::Moraga
Moraga(const dictionary &dict, const phasePair &pair)
Construct from a dictionary and a phase pair.
Definition: Moraga.C:48
Foam::liftModels::defineTypeNameAndDebug
defineTypeNameAndDebug(constantLiftCoefficient, 0)
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::sqr
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Definition: dimensionedSymmTensor.C:51
Foam::liftModel
Definition: liftModel.H:55
fvcGrad.H
Calculate the gradient of the given field.
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::Re
scalarField Re(const UList< complex > &cf)
Extract real component.
Definition: complexField.C:159
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::liftModels::addToRunTimeSelectionTable
addToRunTimeSelectionTable(liftModel, constantLiftCoefficient, dictionary)
WarningInFunction
#define WarningInFunction
Report a warning using Foam::Warning.
Definition: messageStream.H:328
Foam::liftModels::Moraga::~Moraga
virtual ~Moraga()
Destructor.
Definition: Moraga.C:59