relativeVelocityModel.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2014-2023 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"
28 #include "slipFvPatchFields.H"
30 #include "fvcGrad.H"
31 #include "fvcDiv.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
39 }
40 
41 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
42 
43 Foam::wordList Foam::relativeVelocityModel::UdmPatchFieldTypes() const
44 {
45  const volVectorField& U = mixture_.U();
46 
47  wordList UdmTypes
48  (
49  U.boundaryField().size(),
50  calculatedFvPatchScalarField::typeName
51  );
52 
53  forAll(U.boundaryField(), i)
54  {
55  if
56  (
57  isA<fixedValueFvPatchVectorField>(U.boundaryField()[i])
58  || isA<slipFvPatchVectorField>(U.boundaryField()[i])
59  || isA<partialSlipFvPatchVectorField>(U.boundaryField()[i])
60  )
61  {
62  UdmTypes[i] = fixedValueFvPatchVectorField::typeName;
63  }
64  }
65 
66  return UdmTypes;
67 }
68 
69 
70 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
71 
73 (
74  const dictionary& dict,
77 )
78 :
79  mixture_(mixture),
80  g_(g),
81  Udm_
82  (
83  IOobject
84  (
85  "Udm",
86  mixture_.time().name(),
87  mixture_.mesh(),
88  IOobject::READ_IF_PRESENT,
89  IOobject::AUTO_WRITE
90  ),
91  mixture_.mesh(),
93  UdmPatchFieldTypes()
94  )
95 {}
96 
97 
98 // * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
99 
101 (
102  const dictionary& dict,
105 )
106 {
107  const word modelType(dict.lookup(typeName));
108 
109  Info<< "Selecting relative velocity model " << modelType << endl;
110 
111  dictionaryConstructorTable::iterator cstrIter =
112  dictionaryConstructorTablePtr_->find(modelType);
113 
114  if (cstrIter == dictionaryConstructorTablePtr_->end())
115  {
117  << "Unknown time scale model type " << modelType
118  << ", constructor not in hash table" << nl << nl
119  << " Valid time scale model types are:" << nl
120  << dictionaryConstructorTablePtr_->sortedToc()
121  << abort(FatalError);
122  }
123 
124  return
126  (
127  cstrIter()
128  (
129  dict.optionalSubDict(modelType + "Coeffs"),
130  mixture,
131  g
132  )
133  );
134 }
135 
136 
137 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
138 
140 {}
141 
142 
143 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
144 
147 {
148  // Dispersed phase velocity
149  // const volVectorField Ud(mixture_.U() + Udm_);
150 
151  // Use the mixture rather than the dispersed-phase velocity to approximate
152  // the dispersed-phase acceleration to improve stability as only the mixture
153  // momentum equation is coupled to continuity and pressure
154  //
155  // This approximation is valid only in the limit of small drift-velocity.
156  // For large drift-velocity an Euler-Euler approach should be used in
157  // which both the continuous and dispersed-phase momentum equations are
158  // solved and coupled to the pressure.
159  const volVectorField& Ud = mixture_.U();
160 
161  return g_ - (Ud & fvc::grad(Ud));
162 }
163 
164 
166 {
167  const volScalarField betac(mixture_.alphac()*mixture_.rhoc());
168  const volScalarField betad(mixture_.alphad()*mixture_.rhod());
169 
170  // Calculate the relative velocity of the continuous phase w.r.t the mean
171  const volVectorField Ucm(betad*Udm_/betac);
172 
174  (
175  "tauDm",
176  betad*sqr(Udm_) + betac*sqr(Ucm)
177  );
178 }
179 
180 
182 {
183  return fvc::div(tauDm());
184 }
185 
186 
187 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Generic GeometricField class.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:710
const dictionary & optionalSubDict(const word &) const
Find and return a sub-dictionary if found.
Definition: dictionary.C:921
Class to represent a mixture of two constant density phases.
const volVectorField & U() const
Return the mixture velocity field.
tmp< volSymmTensorField > tauDm() const
Return the stress tensor due to the phase transport.
virtual ~relativeVelocityModel()
Destructor.
tmp< volVectorField > acceleration() const
Return the dispersed phase acceleration.
static autoPtr< relativeVelocityModel > New(const dictionary &dict, const incompressibleDriftFluxMixture &mixture, const uniformDimensionedVectorField &g)
tmp< volVectorField > divDevTau() const
Return the div stress tensor due to the phase transport.
const incompressibleDriftFluxMixture & mixture_
Mixture properties.
relativeVelocityModel(const dictionary &dict, const incompressibleDriftFluxMixture &mixture, const uniformDimensionedVectorField &g)
Construct from components.
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
Calculate the divergence of the given field.
Calculate the gradient of the given field.
U
Definition: pEqn.H:72
tmp< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
Definition: fvcGrad.C:46
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
Definition: fvcDiv.C:47
Namespace for OpenFOAM.
static const zero Zero
Definition: zero.H:97
List< word > wordList
A List of words.
Definition: fileName.H:54
VolField< vector > volVectorField
Definition: volFieldsFwd.H:65
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
errorManip< error > abort(error &err)
Definition: errorManip.H:131
messageStream Info
defineTypeNameAndDebug(combustionModel, 0)
const dimensionSet dimVelocity
error FatalError
static const char nl
Definition: Ostream.H:266
dictionary dict