stabilityBlendingFactor.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2018-2020 OpenCFD Ltd.
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 Class
27  Foam::functionObjects::stabilityBlendingFactor
28 
29 Group
30  grpFieldFunctionObjects
31 
32 Description
33  Computes the \c stabilityBlendingFactor to be used by the
34  local blended convection scheme. The output is a surface field weight
35  between 0-1.
36 
37  The weight of a blended scheme, i.e. \c w, is given by a function of
38  the blending factor, \c f:
39 
40  \f[
41  w = f_{scheme_1} + (1 - f_{scheme_2})
42  \f]
43 
44  The factor is calculated based on six criteria:
45  \verbatim
46  1. mesh non-orthogonality field
47  2. magnitude of cell centres gradient
48  3. convergence rate of residuals
49  4. faceWeight
50  5. skewness
51  6. Courant number
52  \endverbatim
53 
54  The user can enable them individually.
55 
56  For option 1, the following relation is used, where \f$\phi_1\f$ is
57  the non-orthogonality:
58  \f[
59  fNon =
60  min
61  (
62  max
63  (
64  0.0,
65  (\phi_1 - max(\phi_1))
66  /(min(\phi_1) - max(\phi_1))
67  ),
68  1.0
69  )
70  \f]
71 
72  For option 2, the following relation is used, where \f$\phi_2\f$ is
73  the magnitude of cell centres gradient (Note that \f$\phi_2 = 3\f$
74  for orthogonal meshes):
75 
76  \f[
77  fMagGradCc =
78  min
79  (
80  max
81  (
82  0.0,
83  (\phi_2 - max(\phi_2))
84  / (min(\phi_2) - max(\phi_2))
85  ),
86  1.0
87  )
88  \f]
89 
90  For option 3, a PID control is used in order to control residual
91  unbounded fluctuations for individual cells.
92 
93  \f[
94  factor =
95  P*residual
96  + I*residualIntegral
97  + D*residualDifferential
98  \f]
99 
100  where \c P, \c I and \c D are user inputs.
101 
102  The following relation is used:
103  \f[
104  fRes = (factor - meanRes)/(maxRes*meanRes);
105  \f]
106 
107  where
108  \vartable
109  meanRes | Average(residual)
110  maxRes | User input
111  \endvartable
112 
113  Note that \f$f_{Res}\f$ will blend more towards one as
114  the cell residual is larger then the domain mean residuals.
115 
116 
117  For option 4, the following relation is used, where \f$\phi_4\f$ is
118  the face weight (Note that \f$\phi_4 = 0.5\f$ for orthogonal meshes):
119 
120  \f[
121  ffaceWeight = min
122  (
123  max
124  (
125  0.0,
126  (min(\phi_4) - \phi_4)
127  / (min(\phi_4) - max(\phi_4))
128  ),
129  1.0
130  )
131  \f]
132 
133 
134  For option 5, the following relation is used, where \f$\phi_5\f$ is
135  the cell skewness:
136 
137  \f[
138  fskewness =
139  min
140  (
141  max
142  (
143  0.0,
144  (\phi_5 - max(\phi_5))
145  / (min(\phi_5) - max(\phi_5))
146  ),
147  1.0
148  )
149  \f]
150 
151 
152  For option 6, the following relation is used:
153 
154  \f[
155  fCoWeight = min(max((Co - Co1)/(Co2 - Co1), 0), 1)
156  \f]
157 
158  where
159  \vartable
160  Co1 | Courant number below which scheme2 is used
161  Co2 | Courant number above which scheme1 is used
162  \endvartable
163 
164  The final factor is determined by:
165 
166  \f[
167  f = max(fNon, fMagGradCc, fRes, ffaceWeight, fskewness, fCoWeight)
168  \f]
169 
170  An indicator (volume) field, named \c blendedIndicator
171  is generated if the log flag is on:
172  - 1 represent scheme1 as active,
173  - 0 represent scheme2 as active.
174 
175  Additional reporting is written to the standard output, providing
176  statistics as to the number of cells used by each scheme.
177 
178  Operands:
179  \table
180  Operand | Type | Location
181  input | - | -
182  output file | dat | $FOAM_CASE/postProcessing/<FO>/<time>/<file>
183  output field | volScalarField | $FOAM_CASE/<time>/<outField>
184  \endtable
185 
186 Usage
187  Minimal example by using \c system/controlDict.functions:
188  \verbatim
189  stabilityBlendingFactor1
190  {
191  // Mandatory entries (unmodifiable)
192  type stabilityBlendingFactor;
193  libs (fieldFunctionObjects);
194 
195  // Mandatory entries (unmodifiable)
196  field <field>; // U;
197  result <outField>; // UBlendingFactor;
198 
199  // Optional entries (runtime modifiable)
200  tolerance 0.001;
201 
202  // Any of the options can be chosen in combinations
203 
204  // Option-1
205  switchNonOrtho true;
206  nonOrthogonality nonOrthoAngle;
207  maxNonOrthogonality 20;
208  minNonOrthogonality 60;
209 
210  // Option-2
211  switchGradCc true;
212  maxGradCc 3;
213  minGradCc 4;
214 
215  // Option-3
216  switchResiduals true;
217  maxResidual 10;
218  residual initialResidual:p;
219  P 1.5;
220  I 0;
221  D 0.5;
222 
223  // Option-4
224  switchFaceWeight true;
225  maxFaceWeight 0.3;
226  minFaceWeight 0.2;
227 
228  // Option-5
229  switchSkewness true;
230  maxSkewness 2;
231  minSkewness 3;
232 
233  // Option-6
234  switchCo true;
235  U U;
236  Co1 1;
237  Co2 10;
238 
239  // Optional (inherited) entries
240  ...
241  }
242  \endverbatim
243 
244  Example of function object specification to calculate the \c residuals used
245  by \c stabilityBlendingFactor. The following writes 'initialResidual:p'
246  field
247  \verbatim
248  residuals
249  {
250  type residuals;
251  libs (utilityFunctionObjects);
252  writeFields true;
253  writeControl writeTime;
254  fields (p);
255  }
256  \endverbatim
257 
258  where the entries mean:
259  \table
260  Property | Description | Type | Req'd | Dflt
261  type | Type name: stabilityBlendingFactor | word | yes | -
262  libs | Library name: fieldFunctionObjects | word | yes | -
263  field | Name of operand field | word | yes | -
264  result | Name of surface field to be used in the localBlended scheme <!--
265  --> | word | yes
266  switchNonOrtho | Select non-orthogonal method | bool | no | false
267  nonOrthogonality | Name of the non-orthogonal field <!--
268  --> | word | no | nonOrthoAngle
269  maxNonOrthogonality| Maximum non-orthogonal for scheme2 | scalar | no | 20
270  minNonOrthogonality| Minimum non-orthogonal for scheme1 | scalar | no | 60
271  switchGradCc | Select cell centre gradient method | bool | no | false
272  maxGradCc| Maximum gradient for scheme2 | scalar | no | 2
273  minGradCc| Minimum gradient for scheme1 | scalar | no | 4
274  switchResiduals | Select residual evolution method | bool | no | false
275  residual | Name of the residual field | word | no | initialResidual:p
276  maxResidual| Maximum residual-mean ratio for scheme1 | scalar | no | 10
277  P | Proportional factor for PID | scalar | no | 3
278  I | Integral factor for PID | scalar | no | 0
279  D | Differential factor for PID | scalar | no | 0.25
280  switchFaceWeight | Select face weight method | bool | no | false
281  faceWeight | Name of the faceWeight field | word | no | faceWeight
282  maxFaceWeight | Maximum face weight for scheme1 | scalar | no | 0.2
283  minFaceWeight | Minimum face weight for scheme2 | scalar | no | 0.3
284  switchSkewness | Select skewness method | bool | no | false
285  skewness | Name of the skewness field | word | no | skewness
286  maxSkewness | Maximum skewness for scheme2 | scalar | no | 2
287  minSkewness | Minimum skewness for scheme1 | scalar | no | 3
288  switchCo | Select Co blended method | bool | no | false
289  U | Name of the flux field for Co blended | word | no | U
290  Co1 | Courant number below which scheme2 is used | scalar | no | 1
291  Co2 | Courant number above which scheme1 is used | scalar | no | 10
292  tolerance | Tolerance for number of blended cells | scalar | no | 0.001
293  \endtable
294 
295  The \c result entry is the field which is read by the \c localBlended scheme
296  specified in \c fvSchemes. This name is determined by the \c localBlended
297  class.
298 
299  The inherited entries are elaborated in:
300  - \link functionObject.H \endlink
301  - \link fieldExpression.H \endlink
302  - \link writeFile.H \endlink
303 
304  Usage by the \c postProcess utility is not available.
305 
306 See also
307  - Foam::functionObject
308  - Foam::functionObjects::fieldExpression
309  - Foam::functionObjects::fvMeshFunctionObject
310  - Foam::functionObjects::writeFile
311  - ExtendedCodeGuide::functionObjects::field::stabilityBlendingFactor
312 
313 SourceFiles
314  stabilityBlendingFactor.C
315 
316 \*---------------------------------------------------------------------------*/
317 
318 #ifndef functionObjects_stabilityBlendingFactor_H
319 #define functionObjects_stabilityBlendingFactor_H
320 
321 #include "fieldExpression.H"
322 #include "writeFile.H"
323 #include "volFields.H"
324 
325 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
326 
327 namespace Foam
328 {
329 namespace functionObjects
330 {
331 
332 /*---------------------------------------------------------------------------*\
333  Class stabilityBlendingFactor Declaration
334 \*---------------------------------------------------------------------------*/
335 
336 class stabilityBlendingFactor
337 :
338  public fieldExpression,
339  public writeFile
340 {
341  // Private Member Data
342 
343  //- Cell based blended indicator
344  volScalarField indicator_;
345 
346  // Switches
347 
348  //- Switch for non-orthogonality
349  Switch nonOrthogonality_;
350 
351  //- Switch for grad of cell centres
352  Switch gradCc_;
353 
354  //- Switch for residuals
355  Switch residuals_;
356 
357  //- Switch for face weight
358  Switch faceWeight_;
359 
360  //- Switch for skewness
361  Switch skewness_;
362 
363  //- Switch for Co
364  Switch Co_;
365 
366 
367  // Lower and upper limits
368 
369  //- Maximum non-orthogonality for fully scheme 2
370  scalar maxNonOrthogonality_;
371 
372  //- Minimum non-orthogonality for fully scheme 1
373  scalar minNonOrthogonality_;
374 
375  //- Maximum gradcc for fully scheme 2
376  scalar maxGradCc_;
377 
378  //- Minimum gradcc for fully scheme 1
379  scalar minGradCc_;
380 
381  //- Maximum ratio to average residual for scheme 2
382  scalar maxResidual_;
383 
384  //- Minimum face weight for fully scheme 2
385  scalar minFaceWeight_;
386 
387  //- Maximum face weight for fully scheme 1
388  scalar maxFaceWeight_;
389 
390  //- Maximum skewness for fully scheme 2
391  scalar maxSkewness_;
392 
393  //- Minimum skewness for fully scheme 1
394  scalar minSkewness_;
395 
396  //- Maximum Co for fully scheme 2
397  scalar Co1_;
398 
399  //- Minimum Co for fully scheme 1
400  scalar Co2_;
401 
402 
403  // File names
404 
405  //- Name of the non-orthogonality field
406  word nonOrthogonalityName_;
407 
408  //- Name of the face weight field
409  word faceWeightName_;
410 
411  //- Name of the skewnes field
412  word skewnessName_;
413 
414  //- Name of the residual field
415  word residualName_;
416 
417  //- Name of the U used for Co based blended
418  word UName_;
419 
420 
421  //- Tolerance used when calculating the number of blended cells
422  scalar tolerance_;
423 
424 
425  //- Error fields
426  scalarField error_;
427  scalarField errorIntegral_;
428  scalarField oldError_;
429  scalarField oldErrorIntegral_;
430 
431  //- Proportional gain
432  scalar P_;
433 
434  //- Integral gain
435  scalar I_;
436 
437  //- Derivative gain
438  scalar D_;
439 
440 
441  // Private Member Functions
442 
443  //- Init fields
444  bool init(bool first);
445 
446  //- Calculate statistics
447  void calcStats(label&, label&, label&) const ;
448 
449  //- Calculate the blending factor field and return true if successful
450  virtual bool calc();
451 
452 
453 protected:
454 
455  // Protected Member Functions
456 
457  //- Write the file header
458  virtual void writeFileHeader(Ostream& os) const;
459 
460 
461 public:
462 
463  //- Runtime type information
464  TypeName("stabilityBlendingFactor");
465 
466 
467  // Constructors
468 
469  //- Construct from Time and dictionary
471  (
472  const word& name,
473  const Time& runTime,
474  const dictionary& dict
475  );
476 
477  //- No copy construct
479 
480  //- No copy assignment
481  void operator=(const stabilityBlendingFactor&) = delete;
482 
483 
484  //- Destructor
485  virtual ~stabilityBlendingFactor() = default;
486 
487 
488  // Member Functions
489 
490  //- Read the stabilityBlendingFactor data
491  virtual bool read(const dictionary&);
492 
493  //- Write the stabilityBlendingFactor
494  virtual bool write();
495 };
496 
497 
498 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
499 
500 } // End namespace functionObjects
501 } // End namespace Foam
502 
503 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
504 
505 #endif
506 
507 // ************************************************************************* //
volFields.H
runTime
engineTime & runTime
Definition: createEngineTime.H:13
Foam::scalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Definition: primitiveFieldsFwd.H:52
writeFile.H
Foam::Switch
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:77
Foam::functionObjects::stabilityBlendingFactor::operator=
void operator=(const stabilityBlendingFactor &)=delete
No copy assignment.
Foam::Time
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:73
Foam::functionObjects::stabilityBlendingFactor::TypeName
TypeName("stabilityBlendingFactor")
Runtime type information.
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:62
Foam::functionObjects::stabilityBlendingFactor::write
virtual bool write()
Write the stabilityBlendingFactor.
Definition: stabilityBlendingFactor.C:715
Foam::functionObjects::stabilityBlendingFactor::~stabilityBlendingFactor
virtual ~stabilityBlendingFactor()=default
Destructor.
Foam::functionObjects::stabilityBlendingFactor
Computes the stabilityBlendingFactor to be used by the local blended convection scheme....
Definition: stabilityBlendingFactor.H:556
Foam::Field< scalar >
Foam::volScalarField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:121
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::functionObjects::fieldExpression
Intermediate class for handling field expression function objects (e.g. blendingFactor etc....
Definition: fieldExpression.H:103
Foam::functionObject::name
const word & name() const
Return the name of this functionObject.
Definition: functionObject.C:131
fieldExpression.H
Foam::functionObjects::writeFile
Base class for writing single files from the function objects.
Definition: writeFile.H:119
Foam::functionObjects::stabilityBlendingFactor::stabilityBlendingFactor
stabilityBlendingFactor(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
Definition: stabilityBlendingFactor.C:426
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::functionObjects::stabilityBlendingFactor::read
virtual bool read(const dictionary &)
Read the stabilityBlendingFactor data.
Definition: stabilityBlendingFactor.C:633
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::functionObjects::stabilityBlendingFactor::writeFileHeader
virtual void writeFileHeader(Ostream &os) const
Write the file header.
Definition: stabilityBlendingFactor.C:85