segregated.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 | Copyright (C) 2014 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 "segregated.H"
27 #include "phasePair.H"
29 #include "fvc.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace dragModels
36 {
37  defineTypeNameAndDebug(segregated, 0);
38  addToRunTimeSelectionTable(dragModel, segregated, dictionary);
39 }
40 }
41 
42 
43 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
44 
46 (
47  const dictionary& dict,
48  const phasePair& pair,
49  const bool registerObject
50 )
51 :
52  dragModel(dict, pair, registerObject),
53  m_("m", dimless, dict.lookup("m")),
54  n_("n", dimless, dict.lookup("n"))
55 {}
56 
57 
58 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
59 
61 {}
62 
63 
64 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
65 
67 {
68  FatalErrorIn("Foam::dragModels::segregated::CdRe() const")
69  << "Not implemented."
70  << "Drag coefficient not defined for the segregated model."
71  << exit(FatalError);
72 
73  return pair_.phase1();
74 }
75 
76 
78 {
79  const fvMesh& mesh(pair_.phase1().mesh());
80 
81  const volScalarField& alpha1(pair_.phase1());
82  const volScalarField& alpha2(pair_.phase2());
83 
84  const volScalarField& rho1(pair_.phase1().rho());
85  const volScalarField& rho2(pair_.phase2().rho());
86 
87  tmp<volScalarField> tnu1(pair_.phase1().nu());
88  tmp<volScalarField> tnu2(pair_.phase2().nu());
89 
90  const volScalarField& nu1(tnu1());
91  const volScalarField& nu2(tnu2());
92 
94  (
95  IOobject
96  (
97  "L",
98  mesh.time().timeName(),
99  mesh
100  ),
101  mesh,
102  dimensionedScalar("L", dimLength, 0),
103  zeroGradientFvPatchField<scalar>::typeName
104  );
105  L.internalField() = cbrt(mesh.V());
106  L.correctBoundaryConditions();
107 
109  (
110  alpha1
111  /max
112  (
113  alpha1 + alpha2,
114  residualAlpha_
115  )
116  );
117  volScalarField magGradI
118  (
119  max
120  (
121  mag(fvc::grad(I)),
122  residualAlpha_/L
123  )
124  );
125 
126  volScalarField muI
127  (
128  rho1*nu1*rho2*nu2
129  /(rho1*nu1 + rho2*nu2)
130  );
131  volScalarField muAlphaI
132  (
133  alpha1*rho1*nu1*alpha2*rho2*nu2
134  /(alpha1*rho1*nu1 + alpha2*rho2*nu2)
135  );
136 
137  volScalarField ReI
138  (
139  pair_.rho()
140  *pair_.magUr()
141  /(magGradI*muI)
142  );
143 
144  volScalarField lambda(m_*ReI + n_*muAlphaI/muI);
145 
146  return lambda*sqr(magGradI)*muI;
147 }
148 
149 
150 // ************************************************************************* //
InternalField & internalField()
Return internal field.
virtual tmp< volScalarField > K() const
The drag function used in the momentum equation.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
Macros for easy insertion into run-time selection tables.
segregated(const dictionary &dict, const phasePair &pair, const bool registerObject)
Construct from components.
error FatalError
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual tmp< volScalarField > CdRe() const
Drag coefficient.
rho1
Definition: pEqn.H:124
rho2
Definition: pEqn.H:125
alpha2
Definition: alphaEqn.H:149
dictionary dict
dimensionedScalar lambda(laminarTransport.lookup("lambda"))
dynamicFvMesh & mesh
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual ~segregated()
Destructor.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
static const sphericalTensor I(1)
volScalarField & alpha1(mixture->alpha1())
defineTypeNameAndDebug(combustionModel, 0)
#define FatalErrorIn(fn)
Definition: error.H:317
dimensionedScalar cbrt(const dimensionedScalar &ds)
A class for managing temporary objects.
Definition: PtrList.H:89
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47