FreeFOAM The Cross-Platform CFD Toolkit
solidWallMixedTemperatureCoupledFvPatchScalarField.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 <OpenFOAM/mapDistribute.H>
33 
34 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
35 
38 (
39  const fvPatch& p,
40  const DimensionedField<scalar, volMesh>& iF
41 )
42 :
43  mixedFvPatchScalarField(p, iF),
44  neighbourFieldName_("undefined-neighbourFieldName"),
45  KName_("undefined-K")
46 {
47  this->refValue() = 0.0;
48  this->refGrad() = 0.0;
49  this->valueFraction() = 1.0;
50 }
51 
52 
55 (
57  const fvPatch& p,
58  const DimensionedField<scalar, volMesh>& iF,
59  const fvPatchFieldMapper& mapper
60 )
61 :
62  mixedFvPatchScalarField(ptf, p, iF, mapper),
63  neighbourFieldName_(ptf.neighbourFieldName_),
64  KName_(ptf.KName_)
65 {}
66 
67 
70 (
71  const fvPatch& p,
72  const DimensionedField<scalar, volMesh>& iF,
73  const dictionary& dict
74 )
75 :
76  mixedFvPatchScalarField(p, iF),
77  neighbourFieldName_(dict.lookup("neighbourFieldName")),
78  KName_(dict.lookup("Kcond"))
79 {
80  if (!isA<directMappedPatchBase>(this->patch().patch()))
81  {
83  (
84  "solidWallMixedTemperatureCoupledFvPatchScalarField::"
85  "solidWallMixedTemperatureCoupledFvPatchScalarField\n"
86  "(\n"
87  " const fvPatch& p,\n"
88  " const DimensionedField<scalar, volMesh>& iF,\n"
89  " const dictionary& dict\n"
90  ")\n"
91  ) << "\n patch type '" << p.type()
92  << "' not type '" << directMappedPatchBase::typeName << "'"
93  << "\n for patch " << p.name()
94  << " of field " << dimensionedInternalField().name()
95  << " in file " << dimensionedInternalField().objectPath()
96  << exit(FatalError);
97  }
98 
99  fvPatchScalarField::operator=(scalarField("value", dict, p.size()));
100 
101  if (dict.found("refValue"))
102  {
103  // Full restart
104  refValue() = scalarField("refValue", dict, p.size());
105  refGrad() = scalarField("refGradient", dict, p.size());
106  valueFraction() = scalarField("valueFraction", dict, p.size());
107  }
108  else
109  {
110  // Start from user entered data. Assume fixedValue.
111  refValue() = *this;
112  refGrad() = 0.0;
113  valueFraction() = 1.0;
114  }
115 }
116 
117 
120 (
122  const DimensionedField<scalar, volMesh>& iF
123 )
124 :
125  mixedFvPatchScalarField(wtcsf, iF),
126  neighbourFieldName_(wtcsf.neighbourFieldName_),
127  KName_(wtcsf.KName_)
128 {}
129 
130 
131 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
132 
135 {
136  return this->patch().lookupPatchField<volScalarField, scalar>(KName_);
137 }
138 
139 
141 {
142  if (updated())
143  {
144  return;
145  }
146 
147  // Get the coupling information from the directMappedPatchBase
148  const directMappedPatchBase& mpp = refCast<const directMappedPatchBase>
149  (
150  patch().patch()
151  );
152  const polyMesh& nbrMesh = mpp.sampleMesh();
153  const fvPatch& nbrPatch = refCast<const fvMesh>
154  (
155  nbrMesh
156  ).boundary()[mpp.samplePolyPatch().index()];
157 
158  // Force recalculation of mapping and schedule
159  const mapDistribute& distMap = mpp.map();
160 
161  tmp<scalarField> intFld = patchInternalField();
162 
163 
165  refCast<const solidWallMixedTemperatureCoupledFvPatchScalarField>
166  (
167  nbrPatch.lookupPatchField<volScalarField, scalar>
168  (
169  neighbourFieldName_
170  )
171  );
172 
173  // Swap to obtain full local values of neighbour internal field
174  scalarField nbrIntFld = nbrField.patchInternalField();
176  (
178  distMap.schedule(),
179  distMap.constructSize(),
180  distMap.subMap(), // what to send
181  distMap.constructMap(), // what to receive
182  nbrIntFld
183  );
184 
185  // Swap to obtain full local values of neighbour K*delta
186  scalarField nbrKDelta = nbrField.K()*nbrPatch.deltaCoeffs();
188  (
190  distMap.schedule(),
191  distMap.constructSize(),
192  distMap.subMap(), // what to send
193  distMap.constructMap(), // what to receive
194  nbrKDelta
195  );
196 
197  tmp<scalarField> myKDelta = K()*patch().deltaCoeffs();
198 
199 
200  // Both sides agree on
201  // - temperature : (myKDelta*fld + nbrKDelta*nbrFld)/(myKDelta+nbrKDelta)
202  // - gradient : (temperature-fld)*delta
203  // We've got a degree of freedom in how to implement this in a mixed bc.
204  // (what gradient, what fixedValue and mixing coefficient)
205  // Two reasonable choices:
206  // 1. specify above temperature on one side (preferentially the high side)
207  // and above gradient on the other. So this will switch between pure
208  // fixedvalue and pure fixedgradient
209  // 2. specify gradient and temperature such that the equations are the
210  // same on both sides. This leads to the choice of
211  // - refGradient = zero gradient
212  // - refValue = neighbour value
213  // - mixFraction = nbrKDelta / (nbrKDelta + myKDelta())
214 
215 
216  this->refValue() = nbrIntFld;
217 
218  this->refGrad() = 0.0;
219 
220  this->valueFraction() = nbrKDelta / (nbrKDelta + myKDelta());
221 
222  mixedFvPatchScalarField::updateCoeffs();
223 
224 
225  if (debug)
226  {
227  scalar Q = gSum(K()*patch().magSf()*snGrad());
228 
229  Info<< patch().boundaryMesh().mesh().name() << ':'
230  << patch().name() << ':'
231  << this->dimensionedInternalField().name() << " <- "
232  << nbrMesh.name() << ':'
233  << nbrPatch.name() << ':'
234  << this->dimensionedInternalField().name() << " :"
235  << " heat[W]:" << Q
236  << " walltemperature "
237  << " min:" << gMin(*this)
238  << " max:" << gMax(*this)
239  << " avg:" << gAverage(*this)
240  << endl;
241  }
242 }
243 
244 
246 (
247  Ostream& os
248 ) const
249 {
250  mixedFvPatchScalarField::write(os);
251  os.writeKeyword("neighbourFieldName")<< neighbourFieldName_
252  << token::END_STATEMENT << nl;
253  os.writeKeyword("Kcond") << KName_ << token::END_STATEMENT << nl;
254 }
255 
256 
257 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
258 
259 namespace Foam
260 {
261 
263 (
266 );
267 
268 } // End namespace Foam
269 
270 // ************************ vim: set sw=4 sts=4 et: ************************ //