externalCoupledTemperatureMixedFvPatchScalarField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2013-2020 OpenFOAM Foundation
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 
28 #include "OFstream.H"
30 
31 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
32 
34 (
35  OFstream& os
36 ) const
37 {
38  os << "# Values: magSf value qDot htc" << endl;
39 }
40 
41 
42 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
43 
46 (
47  const fvPatch& p,
49 )
50 :
52 {}
53 
54 
57 (
59  const fvPatch& p,
61  const fvPatchFieldMapper& mapper
62 )
63 :
65 {}
66 
67 
70 (
71  const fvPatch& p,
73  const dictionary& dict
74 )
75 :
77 {}
78 
79 
82 (
84 )
85 :
87 {}
88 
89 
92 (
95 )
96 :
98 {}
99 
100 
101 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
102 
105 {}
106 
107 
108 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
109 
111 (
112  OFstream& os
113 ) const
114 {
115  if (log())
116  {
117  Info<< type() << ": " << this->patch().name()
118  << ": writing data to " << os.name()
119  << endl;
120  }
121 
122  const label patchi = patch().index();
123 
124  // heat flux [W/m^2]
125  scalarField qDot(this->patch().size(), 0.0);
126 
127  static word ttmName
128  (
130  (
131  thermophysicalTransportModel::typeName,
132  internalField().group()
133  )
134  );
135 
136  static word thermoName(basicThermo::dictName);
137 
138  if (db().foundObject<thermophysicalTransportModel>(ttmName))
139  {
140  const thermophysicalTransportModel& ttm =
142 
143  const basicThermo& thermo = ttm.thermo();
144 
145  const fvPatchScalarField& hep = thermo.he().boundaryField()[patchi];
146 
147  qDot = ttm.alphaEff(patchi)*hep.snGrad();
148  }
149  else if (db().foundObject<basicThermo>(thermoName))
150  {
151  const basicThermo& thermo = db().lookupObject<basicThermo>(thermoName);
152 
153  const fvPatchScalarField& hep = thermo.he().boundaryField()[patchi];
154 
155  qDot = thermo.alpha().boundaryField()[patchi]*hep.snGrad();
156  }
157  else
158  {
160  << "Condition requires either compressible turbulence and/or "
161  << "thermo model to be available" << exit(FatalError);
162  }
163 
164  // patch temperature [K]
165  const scalarField Tp(*this);
166 
167  // near wall cell temperature [K]
168  const scalarField Tc(patchInternalField());
169 
170  // heat transfer coefficient [W/m^2/K]
171  const scalarField htc(qDot/(Tp - Tc + rootVSmall));
172 
173  if (Pstream::parRun())
174  {
175  int tag = Pstream::msgType() + 1;
176 
178  magSfs[Pstream::myProcNo()].setSize(this->patch().size());
179  magSfs[Pstream::myProcNo()] = this->patch().magSf();
180  Pstream::gatherList(magSfs, tag);
181 
183  values[Pstream::myProcNo()].setSize(this->patch().size());
184  values[Pstream::myProcNo()] = Tp;
185  Pstream::gatherList(values, tag);
186 
188  qDots[Pstream::myProcNo()].setSize(this->patch().size());
189  qDots[Pstream::myProcNo()] = qDot;
190  Pstream::gatherList(qDots, tag);
191 
193  htcs[Pstream::myProcNo()].setSize(this->patch().size());
194  htcs[Pstream::myProcNo()] = htc;
195  Pstream::gatherList(htcs, tag);
196 
197  if (Pstream::master())
198  {
199  forAll(values, proci)
200  {
201  const Field<scalar>& magSf = magSfs[proci];
202  const Field<scalar>& value = values[proci];
203  const Field<scalar>& qDot = qDots[proci];
204  const Field<scalar>& htc = htcs[proci];
205 
206  forAll(magSf, facei)
207  {
208  os << magSf[facei] << token::SPACE
209  << value[facei] << token::SPACE
210  << qDot[facei] << token::SPACE
211  << htc[facei] << token::SPACE
212  << nl;
213  }
214  }
215 
216  os.flush();
217  }
218  }
219  else
220  {
221  const Field<scalar>& magSf(this->patch().magSf());
222 
223  forAll(patch(), facei)
224  {
225  os << magSf[facei] << token::SPACE
226  << Tp[facei] << token::SPACE
227  << qDot[facei] << token::SPACE
228  << htc[facei] << token::SPACE
229  << nl;
230  }
231 
232  os.flush();
233  }
234 }
235 
236 
238 (
239  const Pstream::commsTypes comms
240 )
241 {
243 }
244 
245 
247 (
248  Ostream& os
249 ) const
250 {
252 }
253 
254 
255 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
256 
257 namespace Foam
258 {
260  (
263  );
264 }
265 
266 
267 // ************************************************************************* //
const char *const group
Group name for atomic constants.
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
Abstract base-class for fluid and solid thermodynamic properties.
Definition: basicThermo.H:52
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Definition: fvPatchField.C:186
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const scalarField & magSf() const
Return face area magnitudes.
Definition: fvPatch.C:136
label index() const
Return the index of this patch in the fvBoundaryMesh.
Definition: fvPatch.H:173
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
virtual const fluidThermo & thermo() const =0
Access function to incompressible transport model.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
commsTypes
Types of communications.
Definition: UPstream.H:64
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
Output to file stream.
Definition: OFstream.H:82
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
externalCoupledTemperatureMixedFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
virtual const volScalarField & alpha() const
Thermal diffusivity for enthalpy of mixture [kg/m/s].
Definition: basicThermo.C:494
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::commsTypes::blocking)
Evaluate the patch field.
static int & msgType()
Message tag of standard messages.
Definition: UPstream.H:476
rhoReactionThermo & thermo
Definition: createFields.H:28
const fileName & name() const
Return the name of the stream.
Definition: OFstream.H:120
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::commsTypes::blocking)
Evaluate the patch field.
Macros for easy insertion into run-time selection tables.
virtual volScalarField & he()=0
Enthalpy/Internal energy [J/kg].
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:123
A class for handling words, derived from string.
Definition: word.H:59
static word groupName(Name name, const word &group)
Foam::fvPatchFieldMapper.
virtual void flush()
Flush stream.
Definition: OSstream.C:207
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:325
This boundary condition provides a temperature interface to an external application. Values are transferred as plain text files, where OpenFOAM data is written as:
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
static const char nl
Definition: Ostream.H:260
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
Definition: fvPatchField.C:194
virtual void writeHeader(OFstream &os) const
Write header to transfer file.
virtual void transferData(OFstream &os) const
Transfer data for external source.
Abstract base class for thermophysical transport models (RAS, LES and laminar).
void setSize(const label)
Reset size of List.
Definition: List.C:281
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:411
label patchi
virtual tmp< volScalarField > alphaEff() const =0
Effective thermal turbulent diffusivity of mixture [kg/m/s].
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const objectRegistry & db() const
Return local objectRegistry.
Definition: fvPatchField.C:167
messageStream Info
virtual const word & name() const
Return name.
Definition: fvPatch.H:143
volScalarField & p
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
label size() const
Return the number of elements in the UList.
Definition: ListI.H:171
Namespace for OpenFOAM.
const DimensionedField< Type, volMesh > & internalField() const
Return dimensioned internal field reference.
Definition: fvPatchField.H:332