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-2023 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  const dictionary& dict
50 )
51 :
53 {}
54 
55 
58 (
60  const fvPatch& p,
62  const fieldMapper& mapper
63 )
64 :
65  externalCoupledMixedFvPatchField<scalar>(ptf, p, iF, mapper)
66 {}
67 
68 
71 (
74 )
75 :
76  externalCoupledMixedFvPatchField<scalar>(ecmpf, iF)
77 {}
78 
79 
80 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
81 
84 {}
85 
86 
87 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
88 
90 (
91  OFstream& os
92 ) const
93 {
94  if (log())
95  {
96  Info<< type() << ": " << this->patch().name()
97  << ": writing data to " << os.name()
98  << endl;
99  }
100 
101  const label patchi = patch().index();
102 
103  // heat flux [W/m^2]
104  scalarField qDot(this->patch().size(), 0.0);
105 
107  db().lookupType<fluidThermophysicalTransportModel>
108  (
109  internalField().group()
110  );
111 
112  // patch temperature [K]
113  const scalarField Tp(*this);
114 
115  qDot = ttm.kappaEff(patchi)*snGrad();
116 
117  // near wall cell temperature [K]
118  const scalarField Tc(patchInternalField());
119 
120  // heat transfer coefficient [W/m^2/K]
121  const scalarField htc(qDot/(Tp - Tc + rootVSmall));
122 
123  if (Pstream::parRun())
124  {
125  int tag = Pstream::msgType() + 1;
126 
128  magSfs[Pstream::myProcNo()].setSize(this->patch().size());
129  magSfs[Pstream::myProcNo()] = this->patch().magSf();
130  Pstream::gatherList(magSfs, tag);
131 
133  values[Pstream::myProcNo()].setSize(this->patch().size());
134  values[Pstream::myProcNo()] = Tp;
135  Pstream::gatherList(values, tag);
136 
138  qDots[Pstream::myProcNo()].setSize(this->patch().size());
139  qDots[Pstream::myProcNo()] = qDot;
140  Pstream::gatherList(qDots, tag);
141 
143  htcs[Pstream::myProcNo()].setSize(this->patch().size());
144  htcs[Pstream::myProcNo()] = htc;
145  Pstream::gatherList(htcs, tag);
146 
147  if (Pstream::master())
148  {
149  forAll(values, proci)
150  {
151  const Field<scalar>& magSf = magSfs[proci];
152  const Field<scalar>& value = values[proci];
153  const Field<scalar>& qDot = qDots[proci];
154  const Field<scalar>& htc = htcs[proci];
155 
156  forAll(magSf, facei)
157  {
158  os << magSf[facei] << token::SPACE
159  << value[facei] << token::SPACE
160  << qDot[facei] << token::SPACE
161  << htc[facei] << token::SPACE
162  << nl;
163  }
164  }
165 
166  os.flush();
167  }
168  }
169  else
170  {
171  const Field<scalar>& magSf(this->patch().magSf());
172 
173  forAll(patch(), facei)
174  {
175  os << magSf[facei] << token::SPACE
176  << Tp[facei] << token::SPACE
177  << qDot[facei] << token::SPACE
178  << htc[facei] << token::SPACE
179  << nl;
180  }
181 
182  os.flush();
183  }
184 }
185 
186 
188 (
189  const Pstream::commsTypes comms
190 )
191 {
193 }
194 
195 
197 (
198  Ostream& os
199 ) const
200 {
202 }
203 
204 
205 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
206 
207 namespace Foam
208 {
210  (
213  );
214 }
215 
216 
217 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:91
void setSize(const label)
Reset size of List.
Definition: List.C:281
Output to file stream.
Definition: OFstream.H:86
const fileName & name() const
Return the name of the stream.
Definition: OFstream.H:120
virtual void flush()
Flush stream.
Definition: OSstream.C:223
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
commsTypes
Types of communications.
Definition: UPstream.H:65
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:411
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
static int & msgType()
Message tag of standard messages.
Definition: UPstream.H:476
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
This boundary condition provides an interface to an external application. Values are transferred as p...
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::commsTypes::blocking)
Evaluate the patch field.
This boundary condition provides a temperature interface to an external application....
externalCoupledTemperatureMixedFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
Construct from patch, internal field and dictionary.
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::commsTypes::blocking)
Evaluate the patch field.
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 field mapping.
Definition: fieldMapper.H:48
Abstract base class for fluid thermophysical transport models RAS, LES and laminar.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:88
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:64
virtual tmp< volScalarField > kappaEff() const =0
Effective thermal turbulent conductivity.
label patchi
tmp< SurfaceField< Type > > snGrad(const VolField< Type > &vf, const word &name)
Definition: fvcSnGrad.C:45
Namespace for OpenFOAM.
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
messageStream Info
dimensionedScalar log(const dimensionedScalar &ds)
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
static const char nl
Definition: Ostream.H:266
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
dictionary dict
volScalarField & p