FreeFOAM The Cross-Platform CFD Toolkit
turbulentTemperatureCoupledBaffleFvPatchScalarField.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 #include "regionProperties.H"
34 #include <OpenFOAM/mapDistribute.H>
35 
36 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
37 
38 bool Foam::turbulentTemperatureCoupledBaffleFvPatchScalarField::interfaceOwner
39 (
40  const polyMesh& nbrRegion,
41  const polyPatch& nbrPatch
42 ) const
43 {
44  const fvMesh& myRegion = patch().boundaryMesh().mesh();
45 
46  if (nbrRegion.name() == myRegion.name())
47  {
48  return patch().index() < nbrPatch.index();
49  }
50  else
51  {
52  const regionProperties& props =
53  myRegion.objectRegistry::parent().lookupObject<regionProperties>
54  (
55  "regionProperties"
56  );
57 
58  label myIndex = findIndex(props.fluidRegionNames(), myRegion.name());
59  if (myIndex == -1)
60  {
61  label i = findIndex(props.solidRegionNames(), myRegion.name());
62 
63  if (i == -1)
64  {
66  (
67  "turbulentTemperatureCoupledBaffleFvPatchScalarField"
68  "::interfaceOwner(const polyMesh&"
69  ", const polyPatch&)const"
70  ) << "Cannot find region " << myRegion.name()
71  << " neither in fluids " << props.fluidRegionNames()
72  << " nor in solids " << props.solidRegionNames()
73  << exit(FatalError);
74  }
75  myIndex = props.fluidRegionNames().size() + i;
76  }
77  label nbrIndex = findIndex
78  (
79  props.fluidRegionNames(),
80  nbrRegion.name()
81  );
82  if (nbrIndex == -1)
83  {
84  label i = findIndex(props.solidRegionNames(), nbrRegion.name());
85 
86  if (i == -1)
87  {
89  (
90  "coupleManager::interfaceOwner"
91  "(const polyMesh&, const polyPatch&) const"
92  ) << "Cannot find region " << nbrRegion.name()
93  << " neither in fluids " << props.fluidRegionNames()
94  << " nor in solids " << props.solidRegionNames()
95  << exit(FatalError);
96  }
97  nbrIndex = props.fluidRegionNames().size() + i;
98  }
99 
100  return myIndex < nbrIndex;
101  }
102 }
103 
104 
105 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
106 
109 (
110  const fvPatch& p,
112 )
113 :
114  fixedValueFvPatchScalarField(p, iF),
115  neighbourFieldName_("undefined-neighbourFieldName"),
116  KName_("undefined-K")
117 {}
118 
119 
122 (
124  const fvPatch& p,
126  const fvPatchFieldMapper& mapper
127 )
128 :
129  fixedValueFvPatchScalarField(ptf, p, iF, mapper),
130  neighbourFieldName_(ptf.neighbourFieldName_),
131  KName_(ptf.KName_)
132 {}
133 
134 
137 (
138  const fvPatch& p,
140  const dictionary& dict
141 )
142 :
143  fixedValueFvPatchScalarField(p, iF, dict),
144  neighbourFieldName_(dict.lookup("neighbourFieldName")),
145  KName_(dict.lookup("Kcond"))
146 {
147  if (!isA<directMappedPatchBase>(this->patch().patch()))
148  {
150  (
151  "turbulentTemperatureCoupledBaffleFvPatchScalarField::"
152  "turbulentTemperatureCoupledBaffleFvPatchScalarField\n"
153  "(\n"
154  " const fvPatch& p,\n"
155  " const DimensionedField<scalar, volMesh>& iF,\n"
156  " const dictionary& dict\n"
157  ")\n"
158  ) << "\n patch type '" << p.type()
159  << "' not type '" << directMappedPatchBase::typeName << "'"
160  << "\n for patch " << p.name()
161  << " of field " << dimensionedInternalField().name()
162  << " in file " << dimensionedInternalField().objectPath()
163  << exit(FatalError);
164  }
165 }
166 
167 
170 (
173 )
174 :
175  fixedValueFvPatchScalarField(wtcsf, iF),
176  neighbourFieldName_(wtcsf.neighbourFieldName_),
177  KName_(wtcsf.KName_)
178 {}
179 
180 
181 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
182 
185 {
186  const fvMesh& mesh = patch().boundaryMesh().mesh();
187 
188  if (KName_ == "none")
189  {
190  const compressible::RASModel& model =
191  db().lookupObject<compressible::RASModel>("RASProperties");
192 
193  tmp<volScalarField> talpha = model.alphaEff();
194 
195  const basicThermo& thermo =
196  db().lookupObject<basicThermo>("thermophysicalProperties");
197 
198  return
199  talpha().boundaryField()[patch().index()]
200  *thermo.Cp()().boundaryField()[patch().index()];
201  }
202  else if (mesh.objectRegistry::foundObject<volScalarField>(KName_))
203  {
204  return patch().lookupPatchField<volScalarField, scalar>(KName_);
205  }
206  else if (mesh.objectRegistry::foundObject<volSymmTensorField>(KName_))
207  {
208  const symmTensorField& KWall =
209  patch().lookupPatchField<volSymmTensorField, scalar>(KName_);
210 
211  vectorField n = patch().nf();
212 
213  return n & KWall & n;
214  }
215  else
216  {
218  (
219  "turbulentTemperatureCoupledBaffleFvPatchScalarField::K() const"
220  ) << "Did not find field " << KName_
221  << " on mesh " << mesh.name() << " patch " << patch().name()
222  << endl
223  << "Please set 'K' to 'none', a valid volScalarField"
224  << " or a valid volSymmTensorField." << exit(FatalError);
225 
226  return scalarField(0);
227  }
228 }
229 
230 
232 {
233  if (updated())
234  {
235  return;
236  }
237 
238  // Get the coupling information from the directMappedPatchBase
239  const directMappedPatchBase& mpp = refCast<const directMappedPatchBase>
240  (
241  patch().patch()
242  );
243  const polyMesh& nbrMesh = mpp.sampleMesh();
244  const fvPatch& nbrPatch = refCast<const fvMesh>
245  (
246  nbrMesh
247  ).boundary()[mpp.samplePolyPatch().index()];
248 
249  // Force recalculation of mapping and schedule
250  const mapDistribute& distMap = mpp.map();
251  (void)distMap.schedule();
252 
253  tmp<scalarField> intFld = patchInternalField();
254 
255  if (interfaceOwner(nbrMesh, nbrPatch.patch()))
256  {
257  // Note: other side information could be cached - it only needs
258  // to be updated the first time round the iteration (i.e. when
259  // switching regions) but unfortunately we don't have this information.
260 
261 
262  // Calculate the temperature by harmonic averaging
263  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
264 
266  refCast<const turbulentTemperatureCoupledBaffleFvPatchScalarField>
267  (
268  nbrPatch.lookupPatchField<volScalarField, scalar>
269  (
270  neighbourFieldName_
271  )
272  );
273 
274  // Swap to obtain full local values of neighbour internal field
275  scalarField nbrIntFld = nbrField.patchInternalField();
277  (
279  distMap.schedule(),
280  distMap.constructSize(),
281  distMap.subMap(), // what to send
282  distMap.constructMap(), // what to receive
283  nbrIntFld
284  );
285 
286  // Swap to obtain full local values of neighbour K*delta
287  scalarField nbrKDelta = nbrField.K()*nbrPatch.deltaCoeffs();
289  (
291  distMap.schedule(),
292  distMap.constructSize(),
293  distMap.subMap(), // what to send
294  distMap.constructMap(), // what to receive
295  nbrKDelta
296  );
297 
298  tmp<scalarField> myKDelta = K()*patch().deltaCoeffs();
299 
300  // Calculate common wall temperature. Reuse *this to store common value.
301  scalarField Twall
302  (
303  (myKDelta()*intFld() + nbrKDelta*nbrIntFld)
304  / (myKDelta() + nbrKDelta)
305  );
306  // Assign to me
308  // Distribute back and assign to neighbour
310  (
312  distMap.schedule(),
313  nbrField.size(),
314  distMap.constructMap(), // reverse : what to send
315  distMap.subMap(),
316  Twall
317  );
319  (
320  nbrField
322  }
323 
324  if (debug)
325  {
326  //tmp<scalarField> normalGradient =
327  // (*this-intFld())
328  // * patch().deltaCoeffs();
329 
330  scalar Q = gSum(K()*patch().magSf()*snGrad());
331 
332  Info<< patch().boundaryMesh().mesh().name() << ':'
333  << patch().name() << ':'
334  << this->dimensionedInternalField().name() << " <- "
335  << nbrMesh.name() << ':'
336  << nbrPatch.name() << ':'
337  << this->dimensionedInternalField().name() << " :"
338  << " heat[W]:" << Q
339  << " walltemperature "
340  << " min:" << gMin(*this)
341  << " max:" << gMax(*this)
342  << " avg:" << gAverage(*this)
343  << endl;
344  }
345 
346  fixedValueFvPatchScalarField::updateCoeffs();
347 }
348 
349 
351 (
352  Ostream& os
353 ) const
354 {
355  fixedValueFvPatchScalarField::write(os);
356  os.writeKeyword("neighbourFieldName")<< neighbourFieldName_
357  << token::END_STATEMENT << nl;
358  os.writeKeyword("Kcond") << KName_ << token::END_STATEMENT << nl;
359 }
360 
361 
362 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
363 
364 namespace Foam
365 {
366 
368 (
371 );
372 
373 } // End namespace Foam
374 
375 // ************************************************************************* //