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-2021 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 (
58  const fvPatch& p,
60  const dictionary& dict
61 )
62 :
64 {}
65 
66 
69 (
71  const fvPatch& p,
73  const fvPatchFieldMapper& mapper
74 )
75 :
77 {}
78 
79 
82 (
85 )
86 :
88 {}
89 
90 
91 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
92 
95 {}
96 
97 
98 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
99 
101 (
102  OFstream& os
103 ) const
104 {
105  if (log())
106  {
107  Info<< type() << ": " << this->patch().name()
108  << ": writing data to " << os.name()
109  << endl;
110  }
111 
112  const label patchi = patch().index();
113 
114  // heat flux [W/m^2]
115  scalarField qDot(this->patch().size(), 0.0);
116 
117  static word ttmName
118  (
120  (
121  thermophysicalTransportModel::typeName,
122  internalField().group()
123  )
124  );
125 
126  static word thermoName(basicThermo::dictName);
127 
128  if (db().foundObject<thermophysicalTransportModel>(ttmName))
129  {
130  const thermophysicalTransportModel& ttm =
132 
133  const basicThermo& thermo = ttm.thermo();
134 
135  const fvPatchScalarField& hep = thermo.he().boundaryField()[patchi];
136 
137  qDot = ttm.alphaEff(patchi)*hep.snGrad();
138  }
139  else if (db().foundObject<basicThermo>(thermoName))
140  {
141  const basicThermo& thermo = db().lookupObject<basicThermo>(thermoName);
142 
143  const fvPatchScalarField& hep = thermo.he().boundaryField()[patchi];
144 
145  qDot = thermo.alpha().boundaryField()[patchi]*hep.snGrad();
146  }
147  else
148  {
150  << "Condition requires either compressible turbulence and/or "
151  << "thermo model to be available" << exit(FatalError);
152  }
153 
154  // patch temperature [K]
155  const scalarField Tp(*this);
156 
157  // near wall cell temperature [K]
158  const scalarField Tc(patchInternalField());
159 
160  // heat transfer coefficient [W/m^2/K]
161  const scalarField htc(qDot/(Tp - Tc + rootVSmall));
162 
163  if (Pstream::parRun())
164  {
165  int tag = Pstream::msgType() + 1;
166 
168  magSfs[Pstream::myProcNo()].setSize(this->patch().size());
169  magSfs[Pstream::myProcNo()] = this->patch().magSf();
170  Pstream::gatherList(magSfs, tag);
171 
173  values[Pstream::myProcNo()].setSize(this->patch().size());
174  values[Pstream::myProcNo()] = Tp;
175  Pstream::gatherList(values, tag);
176 
178  qDots[Pstream::myProcNo()].setSize(this->patch().size());
179  qDots[Pstream::myProcNo()] = qDot;
180  Pstream::gatherList(qDots, tag);
181 
183  htcs[Pstream::myProcNo()].setSize(this->patch().size());
184  htcs[Pstream::myProcNo()] = htc;
185  Pstream::gatherList(htcs, tag);
186 
187  if (Pstream::master())
188  {
189  forAll(values, proci)
190  {
191  const Field<scalar>& magSf = magSfs[proci];
192  const Field<scalar>& value = values[proci];
193  const Field<scalar>& qDot = qDots[proci];
194  const Field<scalar>& htc = htcs[proci];
195 
196  forAll(magSf, facei)
197  {
198  os << magSf[facei] << token::SPACE
199  << value[facei] << token::SPACE
200  << qDot[facei] << token::SPACE
201  << htc[facei] << token::SPACE
202  << nl;
203  }
204  }
205 
206  os.flush();
207  }
208  }
209  else
210  {
211  const Field<scalar>& magSf(this->patch().magSf());
212 
213  forAll(patch(), facei)
214  {
215  os << magSf[facei] << token::SPACE
216  << Tp[facei] << token::SPACE
217  << qDot[facei] << token::SPACE
218  << htc[facei] << token::SPACE
219  << nl;
220  }
221 
222  os.flush();
223  }
224 }
225 
226 
228 (
229  const Pstream::commsTypes comms
230 )
231 {
233 }
234 
235 
237 (
238  Ostream& os
239 ) const
240 {
242 }
243 
244 
245 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
246 
247 namespace Foam
248 {
250  (
253  );
254 }
255 
256 
257 // ************************************************************************* //
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
fluidReactionThermo & thermo
Definition: createFields.H:28
Base-class for fluid and solid thermodynamic properties.
Definition: basicThermo.H:77
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Definition: fvPatchField.C:166
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:142
label index() const
Return the index of this patch in the fvBoundaryMesh.
Definition: fvPatch.H:174
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
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:62
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
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
virtual const volScalarField & alpha() const =0
Thermal diffusivity for enthalpy of mixture [kg/m/s].
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].
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:337
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 word dictName
Name of the thermophysical properties dictionary.
Definition: basicThermo.H:129
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:174
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:147
messageStream Info
virtual const word & name() const
Return name.
Definition: fvPatch.H:144
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:343