relativeVelocityModel.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 "relativeVelocityModel.H"
27 
28 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
29 
30 namespace Foam
31 {
32  defineTypeNameAndDebug(relativeVelocityModel, 0);
33  defineRunTimeSelectionTable(relativeVelocityModel, dictionary);
34 }
35 
36 
37 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
38 
39 Foam::relativeVelocityModel::relativeVelocityModel
40 (
41  const dictionary& dict,
42  const incompressibleTwoPhaseInteractingMixture& mixture
43 )
44 :
45  mixture_(mixture),
46  alphac_(mixture.alpha2()),
47  alphad_(mixture.alpha1()),
48  rhoc_(mixture.rhoc()),
49  rhod_(mixture.rhod()),
50 
51  Udm_
52  (
53  IOobject
54  (
55  "Udm",
56  alphac_.time().timeName(),
57  alphac_.mesh(),
58  IOobject::NO_READ,
59  IOobject::AUTO_WRITE
60  ),
61  alphac_.mesh(),
62  dimensionedVector("Udm", dimVelocity, vector::zero),
63  mixture.U().boundaryField().types()
64  )
65 {}
66 
67 
68 // * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
69 
71 (
72  const dictionary& dict,
73  const incompressibleTwoPhaseInteractingMixture& mixture
74 )
75 {
76  word modelType(dict.lookup(typeName));
77 
78  Info<< "Selecting relative velocity model " << modelType << endl;
79 
80  dictionaryConstructorTable::iterator cstrIter =
81  dictionaryConstructorTablePtr_->find(modelType);
82 
83  if (cstrIter == dictionaryConstructorTablePtr_->end())
84  {
86  (
87  "relativeVelocityModel::New"
88  "("
89  "const dictionary&"
90  ")"
91  ) << "Unknown time scale model type " << modelType
92  << ", constructor not in hash table" << nl << nl
93  << " Valid time scale model types are:" << nl
94  << dictionaryConstructorTablePtr_->sortedToc()
95  << abort(FatalError);
96  }
97 
98  return
99  autoPtr<relativeVelocityModel>
100  (
101  cstrIter()
102  (
103  dict.subDict(modelType + "Coeffs"),
104  mixture
105  )
106  );
107 }
108 
109 
110 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
111 
113 {}
114 
115 
116 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
117 
118 tmp<volScalarField> Foam::relativeVelocityModel::rho() const
119 {
120  return alphac_*rhoc_ + alphad_*rhod_;
121 }
122 
123 
124 tmp<volSymmTensorField> Foam::relativeVelocityModel::tauDm() const
125 {
126  volScalarField betac(alphac_*rhoc_);
127  volScalarField betad(alphad_*rhod_);
128 
129  // Calculate the relative velocity of the continuous phase w.r.t the mean
130  volVectorField Ucm(betad*Udm_/betac);
131 
132  return tmp<volSymmTensorField>
133  (
135  (
136  "tauDm",
137  betad*sqr(Udm_) + betac*sqr(Ucm)
138  )
139  );
140 }
141 
142 
143 // ************************************************************************* //
tmp< volScalarField > rho() const
Return the mixture mean density.
Info<< "Reading field p_rgh\n"<< endl;volScalarField p_rgh(IOobject("p_rgh", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading field U\n"<< endl;volVectorField U(IOobject("U", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Creating phaseChangeTwoPhaseMixture\n"<< endl;autoPtr< phaseChangeTwoPhaseMixture > mixture
Definition: createFields.H:33
messageStream Info
dynamicFvMesh & mesh
Namespace for OpenFOAM.
virtual ~relativeVelocityModel()
Destructor.
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
Definition: volFieldsFwd.H:58
volScalarField & alpha1
Definition: createFields.H:15
errorManip< error > abort(error &err)
Definition: errorManip.H:131
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
alpha2
Definition: alphaEqn.H:112
error FatalError
static autoPtr< relativeVelocityModel > New(const dictionary &dict, const incompressibleTwoPhaseInteractingMixture &mixture)
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
tmp< volSymmTensorField > tauDm() const
Return the stress tensor due to the phase transport.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimVelocity
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:117
U
Definition: pEqn.H:82
defineTypeNameAndDebug(combustionModel, 0)
word timeName
Definition: getTimeIndex.H:3