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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2014-2017 OpenFOAM Foundation
9 Copyright (C) 2019-2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
31#include "slipFvPatchFields.H"
33
34// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35
36namespace Foam
37{
38 defineTypeNameAndDebug(relativeVelocityModel, 0);
39 defineRunTimeSelectionTable(relativeVelocityModel, dictionary);
40}
41
42// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
43
44Foam::wordList Foam::relativeVelocityModel::UdmPatchFieldTypes() const
45{
46 const volVectorField& U = mixture_.U();
47
48 wordList UdmTypes
49 (
50 U.boundaryField().size(),
51 calculatedFvPatchScalarField::typeName
52 );
53
54 forAll(U.boundaryField(), i)
55 {
56 if
57 (
58 isA<fixedValueFvPatchVectorField>(U.boundaryField()[i])
59 || isA<slipFvPatchVectorField>(U.boundaryField()[i])
60 || isA<partialSlipFvPatchVectorField>(U.boundaryField()[i])
61 )
62 {
63 UdmTypes[i] = fixedValueFvPatchVectorField::typeName;
64 }
65 }
66
67 return UdmTypes;
68}
69
70
71// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
72
74(
75 const dictionary& dict,
76 const incompressibleTwoPhaseInteractingMixture& mixture
77)
78:
79 mixture_(mixture),
80 alphac_(mixture.alpha2()),
81 alphad_(mixture.alpha1()),
82 rhoc_(mixture.rhoc()),
83 rhod_(mixture.rhod()),
84
85 Udm_
86 (
87 IOobject
88 (
89 "Udm",
90 alphac_.time().timeName(),
91 alphac_.mesh(),
92 IOobject::READ_IF_PRESENT,
93 IOobject::AUTO_WRITE
94 ),
95 alphac_.mesh(),
97 UdmPatchFieldTypes()
98 )
99{}
100
101
102// * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
103
105(
106 const dictionary& dict,
107 const incompressibleTwoPhaseInteractingMixture& mixture
108)
109{
110 const word modelType(dict.get<word>(typeName));
111
112 Info<< "Selecting relative velocity model " << modelType << endl;
113
114 auto* ctorPtr = dictionaryConstructorTable(modelType);
115
116 if (!ctorPtr)
117 {
119 (
120 dict,
121 "relative velocity",
122 modelType,
123 *dictionaryConstructorTablePtr_
124 ) << abort(FatalIOError);
125 }
126
127 return
128 autoPtr<relativeVelocityModel>
129 (
130 ctorPtr
131 (
132 dict.optionalSubDict(modelType + "Coeffs"),
133 mixture
134 )
135 );
136}
137
138
139// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
140
142{}
143
144
145// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
146
148{
149 return alphac_*rhoc_ + alphad_*rhod_;
150}
151
152
154{
155 volScalarField betac(alphac_*rhoc_);
156 volScalarField betad(alphad_*rhod_);
157
158 // Calculate the relative velocity of the continuous phase w.r.t the mean
159 volVectorField Ucm(betad*Udm_/betac);
160
161 return tmp<volSymmTensorField>
162 (
164 (
165 "tauDm",
166 betad*sqr(Udm_) + betac*sqr(Ucm)
167 )
168 );
169}
170
171
172// ************************************************************************* //
const volScalarField & alpha1
const volScalarField & alpha2
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
const dictionary & optionalSubDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary, otherwise return this dictionary.
Definition: dictionary.C:577
const volVectorField & U() const
Return const-access to the mixture velocity.
virtual ~relativeVelocityModel()
Destructor.
tmp< volSymmTensorField > tauDm() const
Return the stress tensor due to the phase transport.
tmp< volScalarField > rho() const
Return the mixture mean density.
const incompressibleTwoPhaseInteractingMixture & mixture_
Mixture properties.
const incompressibleTwoPhaseInteractingMixture & mixture() const
Mixture properties.
A class for managing temporary objects.
Definition: tmp.H:65
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
U
Definition: pEqn.H:72
dynamicFvMesh & mesh
#define FatalIOErrorInLookup(ios, lookupTag, lookupName, lookupTable)
Report an error message using Foam::FatalIOError.
Definition: error.H:478
word timeName
Definition: getTimeIndex.H:3
Namespace for OpenFOAM.
List< word > wordList
A List of words.
Definition: fileName.H:63
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:83
dimensionedSymmTensor sqr(const dimensionedVector &dv)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
const dimensionSet dimVelocity
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
errorManip< error > abort(error &err)
Definition: errorManip.H:144
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
IOerror FatalIOError
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
Definition: volFieldsFwd.H:86
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
#define defineRunTimeSelectionTable(baseType, argNames)
Define run-time selection table.
dictionary dict
Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n"<< endl;autoPtr< temperaturePhaseChangeTwoPhaseMixture > mixture
Definition: createFields.H:39
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333