compressibleInterPhaseTransportModel.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) 2017-2018 OpenFOAM Foundation
9 Copyright (C) 2020 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
30
31// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32
34(
35 const volScalarField& rho,
36 const volVectorField& U,
37 const surfaceScalarField& phi,
38 const surfaceScalarField& rhoPhi,
39 const surfaceScalarField& alphaPhi10,
40 const twoPhaseMixtureThermo& mixture
41)
42:
43 twoPhaseTransport_(false),
44 mixture_(mixture),
45 phi_(phi),
46 alphaPhi10_(alphaPhi10)
47{
48 {
49 IOdictionary turbulenceProperties
50 (
51 IOobject
52 (
53 turbulenceModel::propertiesName,
54 U.time().constant(),
55 U.db(),
56 IOobject::MUST_READ,
57 IOobject::NO_WRITE
58 )
59 );
60
61 const word simulationType
62 (
63 turbulenceProperties.get<word>("simulationType")
64 );
65
66 if (simulationType == "twoPhaseTransport")
67 {
68 twoPhaseTransport_ = true;
69 }
70 }
71
72 if (twoPhaseTransport_)
73 {
74 const volScalarField& alpha1(mixture_.alpha1());
75 const volScalarField& alpha2(mixture_.alpha2());
76
77 const volScalarField& rho1 = mixture_.thermo1().rho();
78 const volScalarField& rho2 = mixture_.thermo2().rho();
79
80 alphaRhoPhi1_ =
81 (
83 (
84 IOobject::groupName("alphaRhoPhi", alpha1.group()),
85 fvc::interpolate(rho1)*alphaPhi10_
86 )
87 );
88
89 alphaRhoPhi2_ =
90 (
92 (
93 IOobject::groupName("alphaRhoPhi", alpha2.group()),
94 fvc::interpolate(rho2)*(phi_ - alphaPhi10_)
95 )
96 );
97
98 turbulence1_ =
99 (
100 ThermalDiffusivity<PhaseCompressibleTurbulenceModel<fluidThermo>>
101 ::New
102 (
103 alpha1,
104 rho1,
105 U,
106 alphaRhoPhi1_(),
107 phi,
108 mixture.thermo1()
109 )
110 );
111
112 turbulence2_ =
113 (
114 ThermalDiffusivity<PhaseCompressibleTurbulenceModel<fluidThermo>>
115 ::New
116 (
117 alpha2,
118 rho2,
119 U,
120 alphaRhoPhi2_(),
121 phi,
122 mixture.thermo2()
123 )
124 );
125 }
126 else
127 {
128 turbulence_ = compressible::turbulenceModel::New
129 (
130 rho,
131 U,
132 rhoPhi,
133 mixture
134 );
135
136 turbulence_->validate();
137 }
138}
139
140
141// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
142
145{
146 if (twoPhaseTransport_)
147 {
148 return
149 mixture_.alpha1()*mixture_.thermo1().alphaEff
150 (
151 turbulence1_->alphat()
152 )
153 + mixture_.alpha2()*mixture_.thermo2().alphaEff
154 (
155 turbulence2_->alphat()
156 );
157 }
158 else
159 {
160 return mixture_.alphaEff(turbulence_->alphat());
161 }
162}
163
164
167(
169) const
170{
171 if (twoPhaseTransport_)
172 {
173 return
174 turbulence1_->divDevRhoReff(U)
175 + turbulence2_->divDevRhoReff(U);
176 }
177 else
178 {
179 return turbulence_->divDevRhoReff(U);
180 }
181}
182
183
185{
186 if (twoPhaseTransport_)
187 {
188 const volScalarField& rho1 = mixture_.thermo1().rho();
189 const volScalarField& rho2 = mixture_.thermo2().rho();
190
191 alphaRhoPhi1_.ref() = fvc::interpolate(rho1)*alphaPhi10_;
192 alphaRhoPhi2_.ref() = fvc::interpolate(rho2)*(phi_ - alphaPhi10_);
193 }
194}
195
196
198{
199 if (twoPhaseTransport_)
200 {
201 turbulence1_->correct();
202 turbulence2_->correct();
203 }
204 else
205 {
206 turbulence_->correct();
207 }
208}
209
210
211// ************************************************************************* //
rhoPhi
Definition: rhoEqn.H:10
surfaceScalarField & phi
const volScalarField & alpha1
volScalarField & rho2
const volScalarField & alpha2
volScalarField & rho1
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
virtual tmp< volScalarField > alphaEff(const volScalarField &alphat) const =0
Effective thermal diffusivity of mixture [kg/m/s].
Transport model selection class for the compressibleInterFoam family of solvers.
void correct()
Correct the phase or mixture transport models.
tmp< volScalarField > alphaEff() const
Return the effective temperature transport coefficient.
void correctPhasePhi()
Correct the phase mass-fluxes.
tmp< volVectorField > divDevRhoReff()
Return the effective viscous stress (laminar + turbulent).
Definition: momentumError.C:52
A class for managing temporary objects.
Definition: tmp.H:65
const rhoThermo & thermo1() const
const rhoThermo & thermo2() const
virtual tmp< volScalarField > alphaEff(const volScalarField &alphat) const
Effective thermal diffusivity of mixture [J/m/s/K].
const volScalarField & alpha2() const
Return the phase-fraction of phase 2.
const volScalarField & alpha1() const
Return the phase-fraction of phase 1.
U
Definition: pEqn.H:72
alphaPhi10
Definition: alphaEqn.H:7
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:83
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n"<< endl;autoPtr< temperaturePhaseChangeTwoPhaseMixture > mixture
Definition: createFields.H:39