alphatWallBoilingWallFunctionFvPatchScalarField.H
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) 2015-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 Class
25  Foam::compressible::alphatWallBoilingWallFunctionFvPatchScalarField
26 
27 Description
28  A thermal wall function for simulation of subcooled nucleate wall boiling
29  with runtime selectable submodels for:
30  - wall heat flux partitioning model
31  - nucleation site density
32  - bubble departure frequency
33  - bubble departure diameter
34 
35  Implements a version of the well-known RPI wall boiling model
36  (Kurul & Podowski, 1991). The model implementation based on implementation
37  described in Peltola et al. (2019) and is similar to the model described by
38  Peltola & Pättikangas (2012). The present implementation includes simplified
39  support for presence of non-volatile components in addition to a single
40  volatile component.
41 
42  References:
43  \verbatim
44  Kurul, N., & Podowski, M.Z. (1991).
45  On the modeling of multidimensional effects in boiling channels.
46  ANS. Proc. National Heat Transfer Con. Minneapolis, Minnesota, USA,
47  1991.
48  ISBN: 0-89448-162-1, pp. 30-40.
49  \endverbatim
50 
51  \verbatim
52  Peltola, J., Pättikangas, T., Bainbridge, W., Lehnigk, R., Schlegel, F.
53  (2019).
54  On Development and validation of subcooled nucleate boiling models for
55  OpenFOAM Foundation Release.
56  NURETH-18 Conference Proceedings, Portland, Oregon, United States, 2019.
57  \endverbatim
58 
59  \verbatim
60  Peltola, J., & Pättikangas, T.J.H. (2012).
61  Development and validation of a boiling model for OpenFOAM multiphase
62  solver.
63  CFD4NRS-4 Conference Proceedings, Daejeon, Korea, 2012.
64  paper 59.
65  \endverbatim
66 
67 Usage
68  \table
69  Property | Description | Required | Default value
70  phaseType | 'vapour' or 'liquid' | yes |
71  useLiquidTemperatureWallFunction | \\
72  Use wall function to calculate liquid temperature? | no | yes
73  tolerance | solution tolerance | no | rootSmall
74  Prt | turbulent Prandtl number | no | 0.85
75  \endtable
76 
77  if phaseType 'vapour':
78  \table
79  partitioningModel | | yes |
80  \endtable
81 
82  if phaseType 'liquid':
83  \table
84  partitioningModel | | yes |
85  nucleationSiteModel | | yes |
86  departureDiameterModel | | yes |
87  departureFrequenctModel | | yes |
88  bubbleWaitingTimeRatio | | no | 0.8
89  \endtable
90 
91  NOTE: Runtime selectable submodels may require model specific entries
92 
93  Example usage:
94  \verbatim
95  hotWall
96  {
97  type compressible::alphatWallBoilingWallFunction;
98  phaseType liquid;
99  dmdt uniform 0;
100  Prt 0.85;
101  partitioningModel
102  {
103  type Lavieville;
104  alphaCrit 0.2;
105  }
106  nucleationSiteModel
107  {
108  type LemmertChawla;
109  }
110  departureDiameterModel
111  {
112  type TolubinskiKostanchuk;
113  }
114  departureFrequencyModel
115  {
116  type Cole;
117  }
118  value uniform 0.01;
119  \endverbatim
120 
121 SourceFiles
122  alphatWallBoilingWallFunctionFvPatchScalarField.C
123 
124 \*---------------------------------------------------------------------------*/
125 
126 #ifndef compressible_alphatWallBoilingWallFunctionFvPatchScalarField_H
127 #define compressible_alphatWallBoilingWallFunctionFvPatchScalarField_H
128 
130 #include "fixedValueFvPatchFields.H"
131 #include "partitioningModel.H"
132 #include "nucleationSiteModel.H"
133 #include "departureDiameterModel.H"
134 #include "departureFrequencyModel.H"
135 
136 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
137 
138 namespace Foam
139 {
140 namespace compressible
141 {
142 
143 /*---------------------------------------------------------------------------*\
144  Class alphatWallBoilingWallFunctionFvPatchScalarField Declaration
145 \*---------------------------------------------------------------------------*/
146 
147 class alphatWallBoilingWallFunctionFvPatchScalarField
148 :
149  public fixedValueFvPatchScalarField,
150  public alphatPhaseChangeWallFunctionBase
151 {
152 public:
153 
154  // Data types
155 
156  //- Enumeration listing the possible operational modes
157  enum phaseType
158  {
159  vapourPhase,
160  vaporPhase, // <-- backwards compatibility
162  };
163 
164  //- Heat source type names
165  static const NamedEnum<phaseType, 3> phaseTypeNames_;
166 
167 
168 private:
169 
170  // Private classes
171 
172  //- Struct to hold cached properties for one of the phases
173  struct properties;
174 
175  //- Struct to hold cached properties for the liquid when it is boiling
176  struct boilingLiquidProperties;
177 
178 
179  // Private Data
180 
181  // Controls
182 
183  //- Heat source type
184  phaseType phaseType_;
185 
186  //- Estimate liquid temperature using logarithmic wall function?
187  Switch useLiquidTemperatureWallFunction_;
188 
189  //- Solution tolerance
190  scalar tolerance_;
191 
192 
193  // Parameters
194 
195  //- Turbulent Prandtl number
196  scalar Prt_;
197 
198  //- Bubble waiting time ratio
199  scalar tau_;
200 
201 
202  // Sub-Models
203 
204  //- Heat flux partitioning model
206  partitioningModel_;
207 
208  //- Nucleation site density model
210  nucleationSiteModel_;
211 
212  //- Bubble departure diameter model
214  departureDiameterModel_;
215 
216  //- Bubble departure frequency model
218  departureFrequencyModel_;
219 
220 
221  // State
222 
223  //- Wall liquid fraction
224  scalarField wetFraction_;
225 
226  //- Departure diameter
227  scalarField dDeparture_;
228 
229  //- Departure frequency
230  scalarField fDeparture_;
231 
232  //- Nucleation site density
233  scalarField nucleationSiteDensity_;
234 
235  //- Quenching surface heat flux
236  scalarField qQuenching_;
237 
238  //- Evaporative surface heat flux
239  scalarField qEvaporative_;
240 
241  //- Mass transfer rate
242  scalarField dmdtf_;
243 
244 
245  // Private Member Functions
246 
247  //- Calculate the boiling for the given wall temperature. Return
248  // the total sum of all heat fluxes. Set the properties passed by
249  // non-const reference. Used by the functions below.
250  tmp<scalarField> calcBoiling
251  (
252  const boilingLiquidProperties& props,
253  const scalarField& Tw,
260  ) const;
261 
262  //- Calculate the boiling for the given wall temperature. Return
263  // the total sum of all heat fluxes. Use this to solve the balance
264  // between the heat fluxes specified by the boiling models and the
265  // temperature boundary condition without changing the stored
266  // boiling state.
267  tmp<scalarField> calcBoiling
268  (
269  const boilingLiquidProperties& props,
270  const scalarField& Tw
271  ) const;
272 
273  //- Calculate the boiling for the given wall temperature. Return
274  // the total sum of all heat fluxes. Also set the stored boiling
275  // state. Use this after solving with the final wall temperature
276  // to set the boiling state.
277  tmp<scalarField> evaluateBoiling
278  (
279  const boilingLiquidProperties& props,
280  const scalarField& Tw
281  );
282 
283  //- Get the temperature patch field and the parameters associated
284  // with its boundary condition that are necessary for
285  // approximately evaluating the boundary condition's heat flux at
286  // a given wall temperature.
287  const fvPatchScalarField& getTemperaturePatchField
288  (
289  const boilingLiquidProperties& props,
290  scalarField& isFixed,
291  scalarField& h,
292  scalarField& hTaPlusQa
293  ) const;
294 
295 
296 public:
297 
298  //- Runtime type information
299  TypeName("compressible::alphatWallBoilingWallFunction");
300 
301 
302  // Constructors
303 
304  //- Construct from patch, internal field and dictionary
306  (
307  const fvPatch&,
309  const dictionary&
310  );
311 
312  //- Construct by mapping given
313  // alphatWallBoilingWallFunctionFvPatchScalarField
314  // onto a new patch
316  (
318  const fvPatch&,
320  const fieldMapper&
321  );
322 
323  //- Disallow copy without setting internal field reference
325  (
327  ) = delete;
328 
329  //- Copy constructor setting internal field reference
331  (
334  );
335 
336  //- Construct and return a clone setting internal field reference
338  (
340  ) const
341  {
343  (
345  );
346  }
347 
348 
349  // Member Functions
350 
351  //- Return the wall liquid fraction field [-]
352  const scalarField& wetFraction() const
353  {
354  return wetFraction_;
355  }
356 
357  //- Return the departure diameter field [m]
358  const scalarField& dDeparture() const
359  {
360  return dDeparture_;
361  }
362 
363  //- Return the departure frequency field [Hz]
364  const scalarField& fDeparture() const
365  {
366  return fDeparture_;
367  }
368 
369  //- Return the nucleation site density field [1/m^2]
370  const scalarField& nucleationSiteDensity() const
371  {
372  return nucleationSiteDensity_;
373  }
374 
375  //- Return the quenching surface heat flux field [W/m^2]
376  const scalarField& qQuenching() const
377  {
378  return qQuenching_;
379  }
380 
381  //- Return the evaporative surface heat flux field [W/m^2]
382  const scalarField& qEvaporative() const
383  {
384  return qEvaporative_;
385  }
386 
387  //- Return the rate of phase change
388  virtual const scalarField& dmdtf() const
389  {
390  return dmdtf_;
391  }
392 
393 
394  // Mapping functions
395 
396  //- Map the given fvPatchField onto this fvPatchField
397  virtual void map(const fvPatchScalarField&, const fieldMapper&);
398 
399  //- Reset the fvPatchField to the given fvPatchField
400  // Used for mesh to mesh mapping
401  virtual void reset(const fvPatchScalarField&);
402 
403 
404  // Evaluation functions
405 
406  //- Update the coefficients associated with the patch field
407  virtual void updateCoeffs();
408 
409 
410  // I-O
411 
412  //- Write
413  virtual void write(Ostream&) const;
414 };
415 
416 
417 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
418 
419 } // End namespace compressible
420 } // End namespace Foam
421 
422 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
423 
424 #endif
425 
426 // ************************************************************************* //
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
A thermal wall function for simulation of subcooled nucleate wall boiling with runtime selectable sub...
TypeName("compressible::alphatWallBoilingWallFunction")
Runtime type information.
const scalarField & dDeparture() const
Return the departure diameter field [m].
virtual void reset(const fvPatchScalarField &)
Reset the fvPatchField to the given fvPatchField.
virtual tmp< fvPatchScalarField > clone(const DimensionedField< scalar, volMesh > &iF) const
Construct and return a clone setting internal field reference.
const scalarField & wetFraction() const
Return the wall liquid fraction field [-].
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual void map(const fvPatchScalarField &, const fieldMapper &)
Map the given fvPatchField onto this fvPatchField.
const scalarField & qQuenching() const
Return the quenching surface heat flux field [W/m^2].
alphatWallBoilingWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
Construct from patch, internal field and dictionary.
const scalarField & qEvaporative() const
Return the evaporative surface heat flux field [W/m^2].
const scalarField & fDeparture() const
Return the departure frequency field [Hz].
virtual const scalarField & dmdtf() const
Return the rate of phase change.
const scalarField & nucleationSiteDensity() const
Return the nucleation site density field [1/m^2].
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Abstract base class for field mapping.
Definition: fieldMapper.H:48
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
A class for managing temporary objects.
Definition: tmp.H:55
const dimensionedScalar h
Planck constant.
Namespace for OpenFOAM.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.