coupledMultiphaseTemperatureFvPatchScalarField.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 
27 #include "fieldMapper.H"
28 #include "phaseSystem.H"
32 
33 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
34 
36 (
38  tmp<scalarField>& sumKappaTByDelta,
39  tmp<scalarField>& sumKappaByDelta,
40  scalarField& sumq,
41  tmp<scalarField>& qByKappa
42 ) const
43 {
44  // Lookup the fluid model
45  const phaseSystem& fluid =
46  patch().boundaryMesh().mesh()
48 
49  scalarField sumKappa(size(), scalar(0));
50  scalarField sumKappaT(size(), scalar(0));
51 
52  forAll(fluid.phases(), phasei)
53  {
54  const phaseModel& phase = fluid.phases()[phasei];
55  const rhoThermo& thermo = phase.thermo();
56 
57  const fvPatchScalarField& Tw =
58  thermo.T().boundaryField()[patch().index()];
59 
60  const fvPatchScalarField& alpha =
61  phase.boundaryField()[patch().index()];
62 
63  tmp<scalarField> kappaEff(phase.kappaEff(patch().index()));
64  tmp<scalarField> alphaKappaEff(alpha*kappaEff());
65 
66  if (&Tw == this)
67  {
68  kappa = alphaKappaEff;
69  qByKappa = sumq/kappaEff;
70  sumq -= alpha*sumq;
71  }
72  else
73  {
74  const scalarField T
75  (
76  thermo.T().boundaryField()[patch().index()]
77  .patchInternalField()
78  );
79 
80  sumKappa += alphaKappaEff();
81  sumKappaT += alphaKappaEff*T;
82  }
83  }
84 
85  sumKappaByDelta = sumKappa*patch().deltaCoeffs();
86  sumKappaTByDelta = sumKappaT*patch().deltaCoeffs();
87 }
88 
89 
91 (
92  tmp<scalarField>& sumKappaTByDeltaNbr,
93  tmp<scalarField>& sumKappaByDeltaNbr,
94  tmp<scalarField>& qNbr
95 ) const
96 {
97  // Lookup the fluid model
98  const phaseSystem& fluid =
99  patch().boundaryMesh().mesh()
101 
102  scalarField sumKappa(size(), scalar(0));
103  scalarField sumKappaT(size(), scalar(0));
104 
105  forAll(fluid.phases(), phasei)
106  {
107  const phaseModel& phase = fluid.phases()[phasei];
108  const rhoThermo& thermo = phase.thermo();
109 
110  const fvPatchScalarField& alpha =
111  phase.boundaryField()[patch().index()];
112 
113  const scalarField T
114  (
115  thermo.T().boundaryField()[patch().index()].patchInternalField()
116  );
117 
118  const scalarField alphaKappaEff(alpha*phase.kappaEff(patch().index()));
119 
120  sumKappa += alphaKappaEff;
121  sumKappaT += alphaKappaEff*T;
122  }
123 
124  sumKappaByDeltaNbr = sumKappa*patch().deltaCoeffs();
125  sumKappaTByDeltaNbr = sumKappaT*patch().deltaCoeffs();
126 }
127 
128 
130 (
131  tmp<scalarField>& TrefNbr,
132  tmp<scalarField>& qNbr
133 ) const
134 {
135  // Lookup the fluid model
136  const phaseSystem& fluid =
137  patch().boundaryMesh().mesh()
139 
140  TrefNbr = new scalarField(size(), scalar(0));
141  scalarField& Tw = TrefNbr.ref();
142 
143  forAll(fluid.phases(), phasei)
144  {
145  const phaseModel& phase = fluid.phases()[phasei];
146  const rhoThermo& thermo = phase.thermo();
147 
148  const fvPatchScalarField& alpha =
149  phase.boundaryField()[patch().index()];
150 
151  Tw += alpha*thermo.T().boundaryField()[patch().index()];
152 
153  /*
154  // Pending the addition of phase anisotropic thermal transport
155  tmp<scalarField> qCorr = ttm.qCorr(patch().index());
156 
157  if (qCorr.valid())
158  {
159  if (qCorr.valid())
160  {
161  qNbr.ref() += alpha*qCorr;
162  }
163  else
164  {
165  qNbr = alpha*qCorr;
166  }
167  }
168  */
169  }
170 }
171 
172 
173 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
174 
177 (
178  const fvPatch& p,
180  const dictionary& dict
181 )
182 :
184 {}
185 
186 
189 (
191  const fvPatch& p,
193  const fieldMapper& mapper
194 )
195 :
196  coupledTemperatureFvPatchScalarField(psf, p, iF, mapper)
197 {}
198 
199 
202 (
205 )
206 :
208 {}
209 
210 
211 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
212 
213 namespace Foam
214 {
216  (
219  );
220 }
221 
222 
223 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const Boundary & boundaryField() const
Return const-reference to the boundary field.
virtual const volScalarField & T() const =0
Temperature [K].
Mixed boundary condition for the phase temperature of a phase in an Euler-Euler multiphase simulation...
virtual void getNbr(tmp< scalarField > &sumKappaTByDeltaNbr, tmp< scalarField > &sumKappaByDeltaNbr, tmp< scalarField > &qNbr) const
Get the neighbour patch sum kappa*Tc/delta and kappa/delta.
coupledMultiphaseTemperatureFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
Construct from patch, internal field and dictionary.
virtual void getThis(tmp< scalarField > &kappa, tmp< scalarField > &sumKappaTByDelta, tmp< scalarField > &sumKappaByDelta, scalarField &sumq, tmp< scalarField > &qByKappa) const
Get the patch kappa, sum kappa*Tc/delta and kappa/delta for all phases.
Mixed boundary condition for temperature, to be used for heat-transfer with another region in a CHT c...
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Abstract base class for field mapping.
Definition: fieldMapper.H:48
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:88
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:64
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type and name.
virtual tmp< scalarField > kappaEff(const label patchi) const =0
Effective thermal turbulent conductivity.
virtual const rhoThermo & thermo() const =0
Return the thermophysical model.
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:73
static const word propertiesName
Default name of the phase properties dictionary.
Definition: phaseSystem.H:226
const phaseModelList & phases() const
Return the phase models.
Definition: phaseSystemI.H:102
const fvMesh & mesh() const
Return the mesh.
Definition: phaseSystemI.H:89
Base-class for thermodynamic properties based on density.
Definition: rhoThermo.H:54
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
Namespace for OpenFOAM.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dictionary dict
volScalarField & p
fluidMulticomponentThermo & thermo
Definition: createFields.H:31