compressibleInterPhaseThermophysicalTransportModel.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) 2022-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 
28 #include "fvcSnGrad.H"
29 
30 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
31 
34 (
35  const compressibleInterPhaseTransportModel& momentumTransport
36 )
37 :
38  thermophysicalTransportModel(momentumTransport.mixture_.mesh(), word::null),
39  momentumTransport_(momentumTransport)
40 {}
41 
42 
43 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
44 
46 {
47  return true;
48 }
49 
50 
51 const Foam::dictionary&
53 {
54  return *this;
55 }
56 
57 
60 {
61  const compressibleTwoPhaseVoFMixture& mixture_ =
62  momentumTransport_.mixture_;
63 
64  if (momentumTransport_.twoPhaseTransport_)
65  {
66  return
67  mixture_.alpha1()
68  *(
69  mixture_.thermo1().kappa()
70  + mixture_.thermo1().rho()*mixture_.thermo1().Cp()
71  *momentumTransport_.momentumTransport1_->nut()
72  )
73  + mixture_.alpha2()
74  *(
75  mixture_.thermo2().kappa()
76  + mixture_.thermo2().rho()*mixture_.thermo2().Cp()
77  *momentumTransport_.momentumTransport2_->nut()
78  );
79  }
80  else
81  {
82  return
83  mixture_.alpha1()
84  *(
85  mixture_.thermo1().kappa()
86  + mixture_.thermo1().rho()*mixture_.thermo1().Cp()
87  *momentumTransport_.mixtureMomentumTransport_->nut()
88  )
89  + mixture_.alpha2()
90  *(
91  mixture_.thermo2().kappa()
92  + mixture_.thermo2().rho()*mixture_.thermo2().Cp()
93  *momentumTransport_.mixtureMomentumTransport_->nut()
94  );
95  }
96 }
97 
98 
101 (
102  const label patchi
103 ) const
104 {
105  const compressibleTwoPhaseVoFMixture& mixture_ =
106  momentumTransport_.mixture_;
107 
108  if (momentumTransport_.twoPhaseTransport_)
109  {
110  return
111  mixture_.alpha1().boundaryField()[patchi]
112  *(
113  mixture_.thermo1().kappa().boundaryField()[patchi]
114  + mixture_.thermo1().rho(patchi)
115  *mixture_.thermo1().Cp().boundaryField()[patchi]
116  *momentumTransport_.momentumTransport1_->nut(patchi)
117  )
118  + mixture_.alpha2().boundaryField()[patchi]
119  *(
120  mixture_.thermo2().kappa().boundaryField()[patchi]
121  + mixture_.thermo2().rho(patchi)
122  *mixture_.thermo2().Cp().boundaryField()[patchi]
123  *momentumTransport_.momentumTransport2_->nut(patchi)
124  );
125  }
126  else
127  {
128  return
129  mixture_.alpha1().boundaryField()[patchi]
130  *(
131  mixture_.thermo1().kappa().boundaryField()[patchi]
132  + mixture_.thermo1().rho(patchi)
133  *mixture_.thermo1().Cp().boundaryField()[patchi]
134  *momentumTransport_.mixtureMomentumTransport_->nut(patchi)
135  )
136  + mixture_.alpha2().boundaryField()[patchi]
137  *(
138  mixture_.thermo2().kappa().boundaryField()[patchi]
139  + mixture_.thermo2().rho(patchi)
140  *mixture_.thermo2().Cp().boundaryField()[patchi]
141  *momentumTransport_.mixtureMomentumTransport_->nut(patchi)
142  );
143  }
144 }
145 
146 
149 {
150  const compressibleTwoPhaseVoFMixture& mixture_ =
151  momentumTransport_.mixture_;
152 
153  if (momentumTransport_.twoPhaseTransport_)
154  {
155  return
156  mixture_.alpha1()
157  *(
158  mixture_.thermo1().kappa()
159  + mixture_.thermo1().rho()*mixture_.thermo1().Cp()
160  *momentumTransport_.momentumTransport1_->nut()
161  )/mixture_.thermo1().Cv()
162  + mixture_.alpha2()
163  *(
164  mixture_.thermo2().kappa()
165  + mixture_.thermo2().rho()*mixture_.thermo2().Cp()
166  *momentumTransport_.momentumTransport2_->nut()
167  )/mixture_.thermo2().Cv();
168  }
169  else
170  {
171  return
172  mixture_.alpha1()
173  *(
174  mixture_.thermo1().kappa()
175  + mixture_.thermo1().rho()*mixture_.thermo1().Cp()
176  *momentumTransport_.mixtureMomentumTransport_->nut()
177  )/mixture_.thermo1().Cv()
178  + mixture_.alpha2()
179  *(
180  mixture_.thermo2().kappa()
181  + mixture_.thermo2().rho()*mixture_.thermo2().Cp()
182  *momentumTransport_.mixtureMomentumTransport_->nut()
183  )/mixture_.thermo2().Cv();
184  }
185 }
186 
187 
190 {
192  (
193  "q",
194  -fvc::interpolate(kappaEff())
195  *fvc::snGrad(momentumTransport_.mixture_.T())
196  );
197 }
198 
199 
202 (
203  const label patchi
204 ) const
205 {
206  return tmp<scalarField>(nullptr);
207 }
208 
209 
212 (
214 ) const
215 {
217 
218  return tmp<fvScalarMatrix>(nullptr);
219 }
220 
221 
223 {}
224 
225 
227 {}
228 
229 
230 // ************************************************************************* //
Generic GeometricField class.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
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,.
virtual const volScalarField & Cv() const =0
Heat capacity at constant volume [J/kg/K].
virtual const volScalarField & Cp() const =0
Heat capacity at constant pressure [J/kg/K].
virtual const volScalarField & kappa() const =0
Thermal conductivity of mixture [W/m/K].
virtual void correct()
Solve the thermophysical transport model equations.
virtual tmp< scalarField > qCorr(const label patchi) const
Return the patch heat flux correction [W/m^2].
virtual void predict()
Predict the transport coefficients if possible.
virtual const dictionary & coeffDict() const
Const access to the coefficients dictionary.
tmp< volScalarField > alphaEff() const
Return the effective temperature transport coefficient.
compressibleInterPhaseThermophysicalTransportModel(const compressibleInterPhaseTransportModel &momentumTransport)
Construct from components.
virtual tmp< surfaceScalarField > q() const
Return the heat flux [W/m^2].
virtual tmp< volScalarField > kappaEff() const
Effective thermal turbulent conductivity.
virtual tmp< fvScalarMatrix > divq(volScalarField &he) const
Return the source term for the energy equation.
Transport model selection class for the compressibleInterFoam family of solvers.
Class to represent a mixture of two rhoFluidThermo-based phases.
const rhoFluidThermo & thermo1() const
Return the thermo for phase 1.
const rhoFluidThermo & thermo2() const
Return the thermo for phase 2.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
virtual tmp< volScalarField > rho() const =0
Density [kg/m^3].
Abstract base class for all fluid and solid thermophysical transport models.
A class for managing temporary objects.
Definition: tmp.H:55
const volScalarField & alpha1() const
Return the phase-fraction of phase 1.
const volScalarField & alpha2() const
Return the phase-fraction of phase 2.
A class for handling words, derived from string.
Definition: word.H:62
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:381
Calculate the snGrad of the given volField.
label patchi
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< SurfaceField< Type > > snGrad(const VolField< Type > &vf, const word &name)
Definition: fvcSnGrad.C:45
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
thermo he()