fixedMultiPhaseHeatFluxFvPatchScalarField.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) 2015-2020 OpenFOAM Foundation
9  Copyright (C) 2020-2021 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
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 #include "fvPatchFieldMapper.H"
32 
33 #include "phaseSystem.H"
35 #include "ThermalDiffusivity.H"
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
42 (
43  const fvPatch& p,
45 )
46 :
47  fixedValueFvPatchScalarField(p, iF),
48  q_(p.size(), Zero),
49  relax_(1),
50  Tmin_(273)
51 {}
52 
53 
56 (
57  const fvPatch& p,
59  const dictionary& dict
60 )
61 :
62  fixedValueFvPatchScalarField(p, iF, dict),
63  q_("q", dict, p.size()),
64  relax_(dict.getOrDefault<scalar>("relax", 1)),
65  Tmin_(dict.getOrDefault<scalar>("Tmin", 273))
66 {}
67 
68 
71 (
73  const fvPatch& p,
75  const fvPatchFieldMapper& mapper
76 )
77 :
78  fixedValueFvPatchScalarField(psf, p, iF, mapper),
79  q_(psf.q_, mapper),
80  relax_(psf.relax_),
81  Tmin_(psf.Tmin_)
82 {}
83 
84 
87 (
89 )
90 :
91  fixedValueFvPatchScalarField(psf),
92  q_(psf.q_),
93  relax_(psf.relax_),
94  Tmin_(psf.Tmin_)
95 {}
96 
97 
100 (
103 )
104 :
105  fixedValueFvPatchScalarField(psf, iF),
106  q_(psf.q_),
107  relax_(psf.relax_),
108  Tmin_(psf.Tmin_)
109 {}
110 
111 
112 
113 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
114 
116 {
117  if (updated())
118  {
119  return;
120  }
121 
122  // Lookup the fluid model
123  const phaseSystem& fluid =
124  db().lookupObject<phaseSystem>("phaseProperties");
125 
126  const scalarField& Tp = *this;
127 
128  scalarField A(Tp.size(), Zero);
129  scalarField B(Tp.size(), Zero);
130  scalarField Q(Tp.size(), Zero);
131 
132  forAll(fluid.phases(), phasei)
133  {
134  const phaseModel& phase = fluid.phases()[phasei];
135  const fluidThermo& thermo = phase.thermo();
136 
137  const fvPatchScalarField& alpha =
138  phase.boundaryField()[patch().index()];
139 
140  const fvPatchScalarField& T =
141  thermo.T().boundaryField()[patch().index()];
142 
143  const scalarField kappaEff(phase.kappaEff(patch().index()));
144 
145  if (debug)
146  {
147  const scalarField q0(T.snGrad()*alpha*kappaEff);
148  Q += q0;
149 
150  Info<< patch().name() << " " << phase.name()
151  << ": Heat flux " << gMin(q0) << " - " << gMax(q0) << endl;
152  }
153 
154  A += T.patchInternalField()*alpha*kappaEff*patch().deltaCoeffs();
155  B += alpha*kappaEff*patch().deltaCoeffs();
156  }
157 
158  if (debug)
159  {
160  Info<< patch().name() << " " << ": overall heat flux "
161  << gMin(Q) << " - " << gMax(Q) << " W/m2, power: "
162  << gSum(patch().magSf()*Q) << " W" << endl;
163  }
164 
165  operator==((scalar(1) - relax_)*Tp + relax_*max(Tmin_,(q_ + A)/(B)));
166 
167  fixedValueFvPatchScalarField::updateCoeffs();
168 }
169 
170 
172 (
173  const fvPatchFieldMapper& m
174 )
175 {
176  fixedValueFvPatchScalarField::autoMap(m);
177  m(q_);
178 }
179 
180 
182 (
183  const fvPatchScalarField& ptf,
184  const labelList& addr
185 )
186 {
187  fixedValueFvPatchScalarField::rmap(ptf, addr);
188 
190  refCast<const fixedMultiPhaseHeatFluxFvPatchScalarField>(ptf);
191 
192  q_.rmap(mptf.q_, addr);
193 }
194 
195 
197 {
199  os.writeEntry("relax", relax_);
200  q_.writeEntry("q", os);
201  writeEntry("value", os);
202 }
203 
204 
205 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
206 
207 namespace Foam
208 {
210  (
213  );
214 }
215 
216 
217 // ************************************************************************* //
Foam::expressions::patchExpr::debug
int debug
Static debugging option.
Foam::fvPatchField< scalar >
Foam::fvPatchField::write
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:384
Foam::phaseModel
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:54
p
volScalarField & p
Definition: createFieldRefs.H:8
compressibleTurbulenceModel.H
Foam::phase
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phase.H:54
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
ThermalDiffusivity.H
Foam::constant::atomic::alpha
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Definition: readThermalProperties.H:212
kappaEff
kappaEff
Definition: TEqn.H:10
B
static const Foam::dimensionedScalar B("", Foam::dimless, 18.678)
Foam::fluidThermo
Fundamental fluid thermodynamic properties.
Definition: fluidThermo.H:52
thermo
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
A
static const Foam::dimensionedScalar A("", Foam::dimPressure, 611.21)
phasei
label phasei
Definition: pEqn.H:27
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
fvPatchFieldMapper.H
Foam::gSum
Type gSum(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:594
fluid
twoPhaseSystem & fluid
Definition: setRegionFluidFields.H:3
PhaseCompressibleTurbulenceModel.H
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::fixedMultiPhaseHeatFluxFvPatchScalarField::updateCoeffs
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: fixedMultiPhaseHeatFluxFvPatchScalarField.C:115
Foam::fixedMultiPhaseHeatFluxFvPatchScalarField::rmap
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Definition: fixedMultiPhaseHeatFluxFvPatchScalarField.C:182
Foam::Field< scalar >
Foam::operator==
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::fvPatch
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:65
Foam::T
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
Definition: FieldFieldFunctions.C:58
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::fixedMultiPhaseHeatFluxFvPatchScalarField
Calculates a wall temperature that produces the specified overall wall heat flux across all the phase...
Definition: fixedMultiPhaseHeatFluxFvPatchScalarField.H:121
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
os
OBJstream os(runTime.globalPath()/outputName)
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::phase::name
const word & name() const
Definition: phase.H:111
Foam::foamVersion::patch
const std::string patch
OpenFOAM patch number as a std::string.
Foam::List< label >
fixedMultiPhaseHeatFluxFvPatchScalarField.H
Foam::Ostream::writeEntry
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:236
Foam::fvPatchFieldMapper
Foam::fvPatchFieldMapper.
Definition: fvPatchFieldMapper.H:47
Foam::phaseSystem
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:66
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::makePatchTypeField
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Foam::gMin
Type gMin(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:593
Foam::fixedMultiPhaseHeatFluxFvPatchScalarField::autoMap
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Definition: fixedMultiPhaseHeatFluxFvPatchScalarField.C:172
Foam::fixedMultiPhaseHeatFluxFvPatchScalarField::fixedMultiPhaseHeatFluxFvPatchScalarField
fixedMultiPhaseHeatFluxFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
Definition: fixedMultiPhaseHeatFluxFvPatchScalarField.C:42
Foam::gMax
Type gMax(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:592
Foam::fixedMultiPhaseHeatFluxFvPatchScalarField::write
virtual void write(Ostream &) const
Write.
Definition: fixedMultiPhaseHeatFluxFvPatchScalarField.C:196
Foam::GeometricField::boundaryField
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Definition: GeometricFieldI.H:62
Foam::DimensionedField
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Definition: DimensionedField.H:54