convectiveHeatTransferFvPatchScalarField.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) 2011-2016 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
30#include "fvPatchFieldMapper.H"
32
33// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34
35namespace Foam
36{
37namespace compressible
38{
39
40// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
41
42convectiveHeatTransferFvPatchScalarField::
43convectiveHeatTransferFvPatchScalarField
44(
45 const fvPatch& p,
47)
48:
49 fixedValueFvPatchScalarField(p, iF),
50 L_(1.0)
51{}
52
53
54convectiveHeatTransferFvPatchScalarField::
55convectiveHeatTransferFvPatchScalarField
56(
58 const fvPatch& p,
60 const fvPatchFieldMapper& mapper
61)
62:
63 fixedValueFvPatchScalarField(ptf, p, iF, mapper),
64 L_(ptf.L_)
65{}
66
67
68convectiveHeatTransferFvPatchScalarField::
69convectiveHeatTransferFvPatchScalarField
70(
71 const fvPatch& p,
73 const dictionary& dict
74)
75:
76 fixedValueFvPatchScalarField(p, iF, dict),
77 L_(dict.get<scalar>("L"))
78{}
79
80
81convectiveHeatTransferFvPatchScalarField::
82convectiveHeatTransferFvPatchScalarField
83(
85)
86:
87 fixedValueFvPatchScalarField(htcpsf),
88 L_(htcpsf.L_)
89{}
90
91
92convectiveHeatTransferFvPatchScalarField::
93convectiveHeatTransferFvPatchScalarField
94(
97)
98:
99 fixedValueFvPatchScalarField(htcpsf, iF),
100 L_(htcpsf.L_)
101{}
102
103
104// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
105
107{
108 if (updated())
109 {
110 return;
111 }
112
113 const label patchi = patch().index();
114
115 const compressible::turbulenceModel& turbModel =
116 db().lookupObject<compressible::turbulenceModel>
117 (
119 (
121 internalField().group()
122 )
123 );
124
125 const scalarField alphaEffw(turbModel.alphaEff(patchi));
126
127 const tmp<scalarField> tmuw = turbModel.mu(patchi);
128 const scalarField& muw = tmuw();
129
130 const scalarField& rhow = turbModel.rho().boundaryField()[patchi];
131 const vectorField& Uc = turbModel.U();
132 const vectorField& Uw = turbModel.U().boundaryField()[patchi];
133 const scalarField& Tw = turbModel.transport().T().boundaryField()[patchi];
134 const scalarField& pw = turbModel.transport().p().boundaryField()[patchi];
135 const scalarField Cpw(turbModel.transport().Cp(pw, Tw, patchi));
136
137 const scalarField kappaw(Cpw*alphaEffw);
138 const scalarField Pr(muw*Cpw/kappaw);
139
140 scalarField& htc = *this;
141 forAll(htc, facei)
142 {
143 label celli = patch().faceCells()[facei];
144
145 scalar Re = rhow[facei]*mag(Uc[celli] - Uw[facei])*L_/muw[facei];
146
147 if (Re < 5.0E+05)
148 {
149 htc[facei] = 0.664*sqrt(Re)*cbrt(Pr[facei])*kappaw[facei]/L_;
150 }
151 else
152 {
153 htc[facei] = 0.037*pow(Re, 0.8)*cbrt(Pr[facei])*kappaw[facei]/L_;
154 }
155 }
156
157 fixedValueFvPatchScalarField::updateCoeffs();
158}
159
160
161// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
162
164{
166 os.writeEntry("L", L_);
167 writeEntry("value", os);
168}
169
170
171// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
172
174(
177);
178
179// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
180
181} // End namespace compressible
182} // End namespace Foam
183
184// ************************************************************************* //
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...
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:239
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
virtual tmp< volScalarField > alphaEff() const
Return the effective turbulent thermal diffusivity for enthalpy.
This boundary condition provides a convective heat transfer coefficient condition.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
virtual bool write()
Write the output fields.
A FieldMapper for finite-volume patch fields.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:71
static const word propertiesName
Default name of the phase properties dictionary.
Definition: phaseSystem.H:290
A class for managing temporary objects.
Definition: tmp.H:65
volScalarField & p
OBJstream os(runTime.globalPath()/outputName)
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
Definition: fvPatchField.H:676
bool compressible
Definition: pEqn.H:2
Namespace for OpenFOAM.
scalarField Re(const UList< complex > &cf)
Extract real component.
Definition: complexField.C:159
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
dimensionedScalar cbrt(const dimensionedScalar &ds)
dictionary dict
dimensionedScalar Pr("Pr", dimless, laminarTransport)
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333