FreeFOAM The Cross-Platform CFD Toolkit
turbulentTemperatureCoupledBaffleMixedFvPatchScalarField.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>
34 #include <OpenFOAM/mapDistribute.H>
35 
36 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
37 
38 namespace Foam
39 {
40 namespace compressible
41 {
42 
43 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
44 
45 
46 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
47 
48 turbulentTemperatureCoupledBaffleMixedFvPatchScalarField::
49 turbulentTemperatureCoupledBaffleMixedFvPatchScalarField
50 (
51  const fvPatch& p,
53 )
54 :
55  mixedFvPatchScalarField(p, iF),
56  neighbourFieldName_("undefined-neighbourFieldName"),
57  KName_("undefined-K")
58 {
59  this->refValue() = 0.0;
60  this->refGrad() = 0.0;
61  this->valueFraction() = 1.0;
62 }
63 
64 
65 turbulentTemperatureCoupledBaffleMixedFvPatchScalarField::
66 turbulentTemperatureCoupledBaffleMixedFvPatchScalarField
67 (
69  const fvPatch& p,
71  const fvPatchFieldMapper& mapper
72 )
73 :
74  mixedFvPatchScalarField(ptf, p, iF, mapper),
75  neighbourFieldName_(ptf.neighbourFieldName_),
76  KName_(ptf.KName_)
77 {}
78 
79 
80 turbulentTemperatureCoupledBaffleMixedFvPatchScalarField::
81 turbulentTemperatureCoupledBaffleMixedFvPatchScalarField
82 (
83  const fvPatch& p,
85  const dictionary& dict
86 )
87 :
88  mixedFvPatchScalarField(p, iF),
89  neighbourFieldName_(dict.lookup("neighbourFieldName")),
90  KName_(dict.lookup("Kcond"))
91 {
92  if (!isA<directMappedPatchBase>(this->patch().patch()))
93  {
95  (
96  "turbulentTemperatureCoupledBaffleMixedFvPatchScalarField::"
97  "turbulentTemperatureCoupledBaffleMixedFvPatchScalarField\n"
98  "(\n"
99  " const fvPatch& p,\n"
100  " const DimensionedField<scalar, volMesh>& iF,\n"
101  " const dictionary& dict\n"
102  ")\n"
103  ) << "\n patch type '" << p.type()
104  << "' not type '" << directMappedPatchBase::typeName << "'"
105  << "\n for patch " << p.name()
106  << " of field " << dimensionedInternalField().name()
107  << " in file " << dimensionedInternalField().objectPath()
108  << exit(FatalError);
109  }
110 
111  fvPatchScalarField::operator=(scalarField("value", dict, p.size()));
112 
113  if (dict.found("refValue"))
114  {
115  // Full restart
116  refValue() = scalarField("refValue", dict, p.size());
117  refGrad() = scalarField("refGradient", dict, p.size());
118  valueFraction() = scalarField("valueFraction", dict, p.size());
119  }
120  else
121  {
122  // Start from user entered data. Assume fixedValue.
123  refValue() = *this;
124  refGrad() = 0.0;
125  valueFraction() = 1.0;
126  }
127 }
128 
129 
130 turbulentTemperatureCoupledBaffleMixedFvPatchScalarField::
131 turbulentTemperatureCoupledBaffleMixedFvPatchScalarField
132 (
135 )
136 :
137  mixedFvPatchScalarField(wtcsf, iF),
138  neighbourFieldName_(wtcsf.neighbourFieldName_),
139  KName_(wtcsf.KName_)
140 {}
141 
142 
143 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
144 
147 {
148  const fvMesh& mesh = patch().boundaryMesh().mesh();
149 
150  if (KName_ == "none")
151  {
152  const compressible::RASModel& model =
153  db().lookupObject<compressible::RASModel>("RASProperties");
154 
155  const basicThermo& thermo =
156  db().lookupObject<basicThermo>("thermophysicalProperties");
157 
158  return
159  model.alphaEff()().boundaryField()[patch().index()]
160  *thermo.Cp()().boundaryField()[patch().index()];
161  }
162  else if (mesh.objectRegistry::foundObject<volScalarField>(KName_))
163  {
164  return patch().lookupPatchField<volScalarField, scalar>(KName_);
165  }
166  else if (mesh.objectRegistry::foundObject<volSymmTensorField>(KName_))
167  {
168  const symmTensorField& KWall =
169  patch().lookupPatchField<volSymmTensorField, scalar>(KName_);
170 
171  vectorField n = patch().nf();
172 
173  return n & KWall & n;
174  }
175  else
176  {
178  (
179  "turbulentTemperatureCoupledBaffleMixedFvPatchScalarField::K()"
180  " const"
181  ) << "Did not find field " << KName_
182  << " on mesh " << mesh.name() << " patch " << patch().name()
183  << endl
184  << "Please set 'K' to 'none', a valid volScalarField"
185  << " or a valid volSymmTensorField." << exit(FatalError);
186 
187  return scalarField(0);
188  }
189 }
190 
191 
192 void turbulentTemperatureCoupledBaffleMixedFvPatchScalarField::updateCoeffs()
193 {
194  if (updated())
195  {
196  return;
197  }
198 
199  // Get the coupling information from the directMappedPatchBase
200  const directMappedPatchBase& mpp = refCast<const directMappedPatchBase>
201  (
202  patch().patch()
203  );
204  const polyMesh& nbrMesh = mpp.sampleMesh();
205  const fvPatch& nbrPatch = refCast<const fvMesh>
206  (
207  nbrMesh
208  ).boundary()[mpp.samplePolyPatch().index()];
209 
210  // Force recalculation of mapping and schedule
211  const mapDistribute& distMap = mpp.map();
212 
213  tmp<scalarField> intFld = patchInternalField();
214 
215 
216  // Calculate the temperature by harmonic averaging
217  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
218 
220  refCast
221  <
223  >
224  (
225  nbrPatch.lookupPatchField<volScalarField, scalar>
226  (
227  neighbourFieldName_
228  )
229  );
230 
231  // Swap to obtain full local values of neighbour internal field
232  scalarField nbrIntFld = nbrField.patchInternalField();
234  (
236  distMap.schedule(),
237  distMap.constructSize(),
238  distMap.subMap(), // what to send
239  distMap.constructMap(), // what to receive
240  nbrIntFld
241  );
242 
243  // Swap to obtain full local values of neighbour K*delta
244  scalarField nbrKDelta = nbrField.K()*nbrPatch.deltaCoeffs();
246  (
248  distMap.schedule(),
249  distMap.constructSize(),
250  distMap.subMap(), // what to send
251  distMap.constructMap(), // what to receive
252  nbrKDelta
253  );
254 
255  tmp<scalarField> myKDelta = K()*patch().deltaCoeffs();
256 
257 
258  // Both sides agree on
259  // - temperature : (myKDelta*fld + nbrKDelta*nbrFld)/(myKDelta+nbrKDelta)
260  // - gradient : (temperature-fld)*delta
261  // We've got a degree of freedom in how to implement this in a mixed bc.
262  // (what gradient, what fixedValue and mixing coefficient)
263  // Two reasonable choices:
264  // 1. specify above temperature on one side (preferentially the high side)
265  // and above gradient on the other. So this will switch between pure
266  // fixedvalue and pure fixedgradient
267  // 2. specify gradient and temperature such that the equations are the
268  // same on both sides. This leads to the choice of
269  // - refGradient = zero gradient
270  // - refValue = neighbour value
271  // - mixFraction = nbrKDelta / (nbrKDelta + myKDelta())
272 
273 
274  this->refValue() = nbrIntFld;
275 
276  this->refGrad() = 0.0;
277 
278  this->valueFraction() = nbrKDelta / (nbrKDelta + myKDelta());
279 
280  mixedFvPatchScalarField::updateCoeffs();
281 
282 
283  if (debug)
284  {
285  scalar Q = gSum(K()*patch().magSf()*snGrad());
286 
287  Info<< patch().boundaryMesh().mesh().name() << ':'
288  << patch().name() << ':'
289  << this->dimensionedInternalField().name() << " <- "
290  << nbrMesh.name() << ':'
291  << nbrPatch.name() << ':'
292  << this->dimensionedInternalField().name() << " :"
293  << " heat[W]:" << Q
294  << " walltemperature "
295  << " min:" << gMin(*this)
296  << " max:" << gMax(*this)
297  << " avg:" << gAverage(*this)
298  << endl;
299  }
300 }
301 
302 
303 void turbulentTemperatureCoupledBaffleMixedFvPatchScalarField::write
304 (
305  Ostream& os
306 ) const
307 {
308  mixedFvPatchScalarField::write(os);
309  os.writeKeyword("neighbourFieldName")<< neighbourFieldName_
310  << token::END_STATEMENT << nl;
311  os.writeKeyword("Kcond") << KName_ << token::END_STATEMENT << nl;
312 }
313 
314 
315 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
316 
318 (
321 );
322 
323 
324 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
325 
326 } // End namespace compressible
327 } // End namespace Foam
328 
329 
330 // ************************ vim: set sw=4 sts=4 et: ************************ //