FreeFOAM The Cross-Platform CFD Toolkit
turbulentHeatFluxTemperatureFvPatchScalarField.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 | Copyright (C) 1991-2010 OpenCFD Ltd.
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 
29 #include <finiteVolume/volFields.H>
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 namespace compressible
37 {
38 
39 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40 
41 template<>
42 const char*
43 NamedEnum<turbulentHeatFluxTemperatureFvPatchScalarField::heatSourceType, 2>::
44 names[] =
45  {
46  "power",
47  "flux"
48  };
49 
50 const
51 NamedEnum<turbulentHeatFluxTemperatureFvPatchScalarField::heatSourceType, 2>
52  turbulentHeatFluxTemperatureFvPatchScalarField::heatSourceTypeNames_;
53 
54 
55 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
56 
57 turbulentHeatFluxTemperatureFvPatchScalarField::
58 turbulentHeatFluxTemperatureFvPatchScalarField
59 (
60  const fvPatch& p,
62 )
63 :
64  fixedGradientFvPatchScalarField(p, iF),
65  heatSource_(hsPower),
66  q_(p.size(), 0.0)
67 {}
68 
69 
70 turbulentHeatFluxTemperatureFvPatchScalarField::
71 turbulentHeatFluxTemperatureFvPatchScalarField
72 (
74  const fvPatch& p,
76  const fvPatchFieldMapper& mapper
77 )
78 :
79  fixedGradientFvPatchScalarField(ptf, p, iF, mapper),
80  heatSource_(ptf.heatSource_),
81  q_(ptf.q_, mapper)
82 {}
83 
84 
85 turbulentHeatFluxTemperatureFvPatchScalarField::
86 turbulentHeatFluxTemperatureFvPatchScalarField
87 (
88  const fvPatch& p,
90  const dictionary& dict
91 )
92 :
93  fixedGradientFvPatchScalarField(p, iF),
94  heatSource_(heatSourceTypeNames_.read(dict.lookup("heatSource"))),
95  q_("q", dict, p.size())
96 {
97  fvPatchField<scalar>::operator=(patchInternalField());
98  gradient() = 0.0;
99 }
100 
101 
102 turbulentHeatFluxTemperatureFvPatchScalarField::
103 turbulentHeatFluxTemperatureFvPatchScalarField
104 (
106 )
107 :
108  fixedGradientFvPatchScalarField(thftpsf),
109  heatSource_(thftpsf.heatSource_),
110  q_(thftpsf.q_)
111 {}
112 
113 
114 turbulentHeatFluxTemperatureFvPatchScalarField::
115 turbulentHeatFluxTemperatureFvPatchScalarField
116 (
119 )
120 :
121  fixedGradientFvPatchScalarField(thftpsf, iF),
122  heatSource_(thftpsf.heatSource_),
123  q_(thftpsf.q_)
124 {}
125 
126 
127 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
128 
129 void turbulentHeatFluxTemperatureFvPatchScalarField::autoMap
130 (
131  const fvPatchFieldMapper& m
132 )
133 {
134  fixedGradientFvPatchScalarField::autoMap(m);
135  q_.autoMap(m);
136 }
137 
138 
139 void turbulentHeatFluxTemperatureFvPatchScalarField::rmap
140 (
141  const fvPatchScalarField& ptf,
142  const labelList& addr
143 )
144 {
145  fixedGradientFvPatchScalarField::rmap(ptf, addr);
146 
148  refCast<const turbulentHeatFluxTemperatureFvPatchScalarField>
149  (
150  ptf
151  );
152 
153  q_.rmap(thftptf.q_, addr);
154 }
155 
156 
157 void turbulentHeatFluxTemperatureFvPatchScalarField::updateCoeffs()
158 {
159  if (updated())
160  {
161  return;
162  }
163 
164  const label patchI = patch().index();
165 
166  const RASModel& rasModel = db().lookupObject<RASModel>("RASProperties");
167 
168  const scalarField alphaEffp = rasModel.alphaEff()().boundaryField()[patchI];
169 
170 // const scalarField& Tp = thermo.T().boundaryField()[patchI];
171  const scalarField& Tp = *this;
172 
173  const scalarField Cpp = rasModel.thermo().Cp(Tp, patchI);
174 
175  switch (heatSource_)
176  {
177  case hsPower:
178  {
179  const scalar Ap = gSum(patch().magSf());
180  gradient() = q_/(Ap*Cpp*alphaEffp);
181  break;
182  }
183  case hsFlux:
184  {
185  gradient() = q_/(Cpp*alphaEffp);
186  break;
187  }
188  default:
189  {
191  (
192  "turbulentHeatFluxTemperatureFvPatchScalarField"
193  "("
194  "const fvPatch&, "
195  "const DimensionedField<scalar, volMesh>&, "
196  "const dictionary&"
197  ")"
198  ) << "Unknown heat source type. Valid types are: "
199  << heatSourceTypeNames_ << nl << exit(FatalError);
200  }
201  }
202 
203  fixedGradientFvPatchScalarField::updateCoeffs();
204 }
205 
206 
207 void turbulentHeatFluxTemperatureFvPatchScalarField::write
208 (
209  Ostream& os
210 ) const
211 {
213  os.writeKeyword("heatSource") << heatSourceTypeNames_[heatSource_]
214  << token::END_STATEMENT << nl;
215  q_.writeEntry("q", os);
216  gradient().writeEntry("gradient", os);
217  writeEntry("value", os);
218 }
219 
220 
221 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
222 
224 (
227 );
228 
229 
230 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
231 
232 } // End namespace compressible
233 } // End namespace Foam
234 
235 
236 // ************************ vim: set sw=4 sts=4 et: ************************ //
237