stabilityBlendingFactor.C
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-2021 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
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\*---------------------------------------------------------------------------*/
27
31#include "fvcGrad.H"
32#include "surfaceInterpolate.H"
33
34// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35
36namespace Foam
37{
38namespace functionObjects
39{
42 (
46 );
47}
48}
49
50
51// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
52
54Foam::functionObjects::stabilityBlendingFactor::indicator()
55{
56 const word fldName = "blendedIndicator" + fieldName_;
57
58 auto* fldPtr = mesh_.getObjectPtr<volScalarField>(fldName);
59
60 if (!fldPtr)
61 {
62 fldPtr = new volScalarField
63 (
64 IOobject
65 (
66 "blendedIndicator" + fieldName_,
68 mesh_,
71 ),
72 mesh_,
74 zeroGradientFvPatchScalarField::typeName
75 );
76
77 mesh_.objectRegistry::store(fldPtr);
78 }
79
80 return *fldPtr;
81}
82
83
84void Foam::functionObjects::stabilityBlendingFactor::calcStats
85(
86 label& nCellsScheme1,
87 label& nCellsScheme2,
88 label& nCellsBlended
89) const
90{
91 const auto* indicatorPtr =
92 mesh_.cfindObject<volScalarField>("blendedIndicator" + fieldName_);
93
94 if (!indicatorPtr)
95 {
97 << "Indicator field not set"
98 << abort(FatalError);
99 }
100
101 const auto& indicator = *indicatorPtr;
102
103 for (const scalar i : indicator)
104 {
105 if (i < tolerance_)
106 {
107 ++nCellsScheme2;
108 }
109 else if (i > (1 - tolerance_))
110 {
111 ++nCellsScheme1;
112 }
113 else
114 {
115 ++nCellsBlended;
116 }
117 }
118
119 reduce(nCellsScheme1, sumOp<label>());
120 reduce(nCellsScheme2, sumOp<label>());
121 reduce(nCellsBlended, sumOp<label>());
122}
123
124
126(
127 Ostream& os
128) const
129{
130 writeHeader(os, "Stabilization blending factor");
131 writeCommented(os, "Time");
132 writeTabbed(os, "Scheme1");
133 writeTabbed(os, "Scheme2");
134 writeTabbed(os, "Blended");
135 os << endl;
136}
137
138
139bool Foam::functionObjects::stabilityBlendingFactor::calc()
140{
141 init(false);
142 return true;
143}
144
145
147{
148 const auto* residualPtr = mesh_.findObject<IOField<scalar>>(residualName_);
149
150 auto& indicator = this->indicator();
151
152 if (residuals_)
153 {
154 if (!residualPtr)
155 {
157 << "Could not find residual field : " << residualName_
158 << ". The residual mode will not be " << nl
159 << " considered for the blended field in the stability "
160 << "blending factor. " << nl
161 << " Please add the corresponding 'solverInfo'"
162 << " function object with 'writeResidualFields true'." << nl
163 << " If the solverInfo function object is already enabled "
164 << "you might need to wait " << nl
165 << " for the first iteration."
166 << nl << endl;
167 }
168 else
169 {
170 scalar meanRes = gAverage(mag(*residualPtr)) + VSMALL;
171
172 Log << nl << name() << " : " << nl
173 << " Average(mag(residuals)) : " << meanRes << endl;
174
175 oldError_ = error_;
176 oldErrorIntegral_ = errorIntegral_;
177 error_ = mag(meanRes - mag(*residualPtr));
178 errorIntegral_ = oldErrorIntegral_ + 0.5*(error_ + oldError_);
179 const scalarField errorDifferential(error_ - oldError_);
180
181 const scalarField factorList
182 (
183 + P_*error_
184 + I_*errorIntegral_
185 + D_*errorDifferential
186 );
187
188 const scalarField indicatorResidual
189 (
190 max
191 (
192 min
193 (
194 mag(factorList - meanRes)/(maxResidual_*meanRes),
195 scalar(1)
196 ),
197 scalar(0)
198 )
199 );
200
201 forAll(indicator, i)
202 {
203 indicator[i] = indicatorResidual[i];
204 }
205 }
206 }
207
208 const volScalarField* nonOrthPtr =
209 mesh_.findObject<volScalarField>(nonOrthogonalityName_);
210
211 if (nonOrthogonality_)
212 {
213 if (maxNonOrthogonality_ >= minNonOrthogonality_)
214 {
216 << "minNonOrthogonality should be larger than "
217 << "maxNonOrthogonality."
218 << exit(FatalError);
219 }
220
221 indicator =
222 max
223 (
224 indicator,
225 min
226 (
227 max
228 (
229 scalar(0),
230 (*nonOrthPtr - maxNonOrthogonality_)
231 / (minNonOrthogonality_ - maxNonOrthogonality_)
232 ),
233 scalar(1)
234 )
235 );
236
237 if (first)
238 {
239 Log << " Max non-orthogonality : " << max(*nonOrthPtr).value()
240 << endl;
241 }
242 }
243
244 const auto* skewnessPtr = mesh_.cfindObject<volScalarField>(skewnessName_);
245
246 if (skewness_)
247 {
248 if (maxSkewness_ >= minSkewness_)
249 {
251 << "minSkewness should be larger than maxSkewness."
252 << exit(FatalError);
253 }
254
255 indicator =
256 max
257 (
258 indicator,
259 min
260 (
261 max
262 (
263 scalar(0),
264 (*skewnessPtr - maxSkewness_)
265 / (minSkewness_ - maxSkewness_)
266 ),
267 scalar(1)
268 )
269 );
270
271 if (first)
272 {
273 Log << " Max skewness : " << max(*skewnessPtr).value()
274 << endl;
275 }
276 }
277
278 const auto* faceWeightsPtr =
279 mesh_.cfindObject<volScalarField>(faceWeightName_);
280
281 if (faceWeight_)
282 {
283 if (maxFaceWeight_ >= minFaceWeight_)
284 {
286 << "minFaceWeight should be larger than maxFaceWeight."
287 << exit(FatalError);
288 }
289
290 indicator =
291 max
292 (
293 indicator,
294 min
295 (
296 max
297 (
298 scalar(0),
299 (minFaceWeight_ - *faceWeightsPtr)
300 / (minFaceWeight_ - maxFaceWeight_)
301 ),
302 scalar(1)
303 )
304 );
305
306 if (first)
307 {
308 Log << " Min face weight: " << min(*faceWeightsPtr).value()
309 << endl;
310 }
311 }
312
313
314 if (gradCc_)
315 {
316 if (maxGradCc_ >= minGradCc_)
317 {
319 << "minGradCc should be larger than maxGradCc."
320 << exit(FatalError);
321 }
322
323 auto tmagGradCC = tmp<volScalarField>::New
324 (
325 IOobject
326 (
327 "magGradCC",
328 time_.timeName(),
329 mesh_,
332 ),
333 mesh_,
335 zeroGradientFvPatchScalarField::typeName
336 );
337 auto& magGradCC = tmagGradCC.ref();
338
339 for (direction i=0; i<vector::nComponents; i++)
340 {
341 // Create field with zero grad
343 (
344 IOobject
345 (
346 "cc" + word(vector::componentNames[i]),
347 mesh_.time().timeName(),
348 mesh_,
351 ),
352 mesh_,
354 zeroGradientFvPatchScalarField::typeName
355 );
356 cci = mesh_.C().component(i);
357 cci.correctBoundaryConditions();
358 magGradCC += mag(fvc::grad(cci)).ref();
359 }
360
361 if (first)
362 {
363 Log << " Max magGradCc : " << max(magGradCC.ref()).value()
364 << endl;
365 }
366
367 indicator =
368 max
369 (
370 indicator,
371 min
372 (
373 max
374 (
375 scalar(0),
376 (magGradCC - maxGradCc_)
377 / (minGradCc_ - maxGradCc_)
378 ),
379 scalar(1)
380 )
381 );
382 }
383
384
385 const auto* UNamePtr = mesh_.findObject<volVectorField>(UName_);
386
387 if (Co_)
388 {
389 if (Co1_ >= Co2_)
390 {
392 << "Co2 should be larger than Co1."
393 << exit(FatalError);
394 }
395
396 auto CoPtr = tmp<volScalarField>::New
397 (
398 IOobject
399 (
400 "Co",
401 time_.timeName(),
402 mesh_,
405 ),
406 mesh_,
408 zeroGradientFvPatchScalarField::typeName
409 );
410
411 auto& Co = CoPtr.ref();
412
413 Co.primitiveFieldRef() =
414 mesh_.time().deltaT()*mag(*UNamePtr)/cbrt(mesh_.V());
415
416 indicator =
417 max
418 (
419 indicator,
420 min(max(scalar(0), (Co - Co1_)/(Co2_ - Co1_)), scalar(1))
421 );
422
423 if (first)
424 {
425 Log << " Max Co : " << max(Co).value()
426 << endl;
427 }
428 }
429
430 if (first)
431 {
432 Log << nl;
433 }
434
435 if (log)
436 {
437 label nCellsScheme1 = 0;
438 label nCellsScheme2 = 0;
439 label nCellsBlended = 0;
440
441 calcStats(nCellsScheme1, nCellsScheme2, nCellsBlended);
442
443 Log << nl << name() << " execute :" << nl
444 << " scheme 1 cells : " << nCellsScheme1 << nl
445 << " scheme 2 cells : " << nCellsScheme2 << nl
446 << " blended cells : " << nCellsBlended << nl
447 << endl;
448 }
449
450 indicator.correctBoundaryConditions();
451 indicator.min(1.0);
452 indicator.max(0.0);
453
454 // Update the blended surface field
455 auto& surBlended = mesh_.lookupObjectRef<surfaceScalarField>(resultName_);
456
457 surBlended = fvc::interpolate(indicator);
458
459 return true;
460}
461
462
463// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
464
466(
467 const word& name,
468 const Time& runTime,
469 const dictionary& dict
470)
471:
473 writeFile(obr_, name, typeName, dict),
474 nonOrthogonality_(dict.getOrDefault<Switch>("switchNonOrtho", false)),
475 gradCc_(dict.getOrDefault<Switch>("switchGradCc", false)),
476 residuals_(dict.getOrDefault<Switch>("switchResiduals", false)),
477 faceWeight_(dict.getOrDefault<Switch>("switchFaceWeight", false)),
478 skewness_(dict.getOrDefault<Switch>("switchSkewness", false)),
479 Co_(dict.getOrDefault<Switch>("switchCo", false)),
480
481 maxNonOrthogonality_
482 (
483 dict.getOrDefault<scalar>("maxNonOrthogonality", 20.0)
484 ),
485 minNonOrthogonality_
486 (
487 dict.getOrDefault<scalar>("minNonOrthogonality", 60.0)
488 ),
489 maxGradCc_(dict.getOrDefault<scalar>("maxGradCc", 3.0)),
490 minGradCc_(dict.getOrDefault<scalar>("minGradCc", 4.0)),
491 maxResidual_(dict.getOrDefault<scalar>("maxResidual", 10.0)),
492 minFaceWeight_(dict.getOrDefault<scalar>("minFaceWeight", 0.3)),
493 maxFaceWeight_(dict.getOrDefault<scalar>("maxFaceWeight", 0.2)),
494 maxSkewness_(dict.getOrDefault<scalar>("maxSkewness", 2.0)),
495 minSkewness_(dict.getOrDefault<scalar>("minSkewness", 3.0)),
496 Co1_(dict.getOrDefault<scalar>("Co1", 1.0)),
497 Co2_(dict.getOrDefault<scalar>("Co2", 10.0)),
498
499 nonOrthogonalityName_
500 (
501 dict.getOrDefault<word>("nonOrthogonality", "nonOrthoAngle")
502 ),
503 faceWeightName_
504 (
505 dict.getOrDefault<word>("faceWeight", "faceWeight")
506 ),
507 skewnessName_
508 (
509 dict.getOrDefault<word>("skewness", "skewness")
510 ),
511 residualName_
512 (
513 dict.getOrDefault<word>
514 (
515 "residual",
516 IOobject::scopedName("initialResidual", "p")
517 )
518 ),
519 UName_
520 (
521 dict.getOrDefault<word>("U", "U")
522 ),
523
524 tolerance_(0.001),
525 error_(mesh_.nCells(), Zero),
526 errorIntegral_(mesh_.nCells(), Zero),
527 oldError_(mesh_.nCells(), Zero),
528 oldErrorIntegral_(mesh_.nCells(), Zero),
529 P_(dict.getOrDefault<scalar>("P", 3)),
530 I_(dict.getOrDefault<scalar>("I", 0.0)),
531 D_(dict.getOrDefault<scalar>("D", 0.25))
532{
533 read(dict);
534 setResultName(typeName, "");
535
536 auto faceBlendedPtr = tmp<surfaceScalarField>::New
537 (
539 (
541 time_.timeName(),
542 mesh_,
545 ),
546 mesh_,
548 );
549 store(resultName_, faceBlendedPtr);
550
551 const auto* nonOrthPtr =
552 mesh_.findObject<volScalarField>(nonOrthogonalityName_);
553
554 if (nonOrthogonality_ && !nonOrthPtr)
555 {
556 IOobject fieldHeader
557 (
558 nonOrthogonalityName_,
559 mesh_.time().constant(),
560 mesh_,
563 );
564
565 if (fieldHeader.typeHeaderOk<volScalarField>(true, true, false))
566 {
567 volScalarField* vfPtr = new volScalarField(fieldHeader, mesh_);
568 mesh_.objectRegistry::store(vfPtr);
569 }
570 else
571 {
573 << "Field : " << nonOrthogonalityName_ << " not found."
574 << " The function object will not be used"
575 << exit(FatalError);
576 }
577 }
578
579
580 const auto* faceWeightsPtr =
581 mesh_.findObject<volScalarField>(faceWeightName_);
582
583 if (faceWeight_ && !faceWeightsPtr)
584 {
585 IOobject fieldHeader
586 (
587 faceWeightName_,
588 mesh_.time().constant(),
589 mesh_,
592 );
593
594 if (fieldHeader.typeHeaderOk<volScalarField>(true, true, false))
595 {
596 volScalarField* vfPtr = new volScalarField(fieldHeader, mesh_);
597 mesh_.objectRegistry::store(vfPtr);
598 }
599 else
600 {
602 << "Field : " << faceWeightName_ << " not found."
603 << " The function object will not be used"
604 << exit(FatalError);
605 }
606 }
607
608 const auto* skewnessPtr = mesh_.findObject<volScalarField>(skewnessName_);
609
610 if (skewness_ && !skewnessPtr)
611 {
612 IOobject fieldHeader
613 (
614 skewnessName_,
615 mesh_.time().constant(),
616 mesh_,
619 );
620
621 if (fieldHeader.typeHeaderOk<volScalarField>(true, true, false))
622 {
623 volScalarField* vfPtr(new volScalarField(fieldHeader, mesh_));
624 mesh_.objectRegistry::store(vfPtr);
625 }
626 else
627 {
629 << "Field : " << skewnessName_ << " not found."
630 << " The function object will not be used"
631 << exit(FatalError);
632 }
633 }
634
635 if (log)
636 {
637 indicator().writeOpt(IOobject::AUTO_WRITE);
638 }
639
640 if (writeToFile_)
641 {
643 }
644
645 init(true);
646}
647
648
649// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
650
652(
653 const dictionary& dict
654)
655{
657 {
658 dict.readEntry("switchNonOrtho", nonOrthogonality_);
659 dict.readEntry("switchGradCc", gradCc_);
660 dict.readEntry("switchResiduals", residuals_);
661 dict.readEntry("switchFaceWeight", faceWeight_);
662 dict.readEntry("switchSkewness", skewness_);
663 dict.readEntry("switchCo", Co_);
664
665 dict.readIfPresent("maxNonOrthogonality", maxNonOrthogonality_);
666 dict.readIfPresent("maxGradCc", maxGradCc_);
667 dict.readIfPresent("maxResidual", maxResidual_);
668 dict.readIfPresent("maxSkewness", maxSkewness_);
669 dict.readIfPresent("maxFaceWeight", maxFaceWeight_);
670 dict.readIfPresent("Co2", Co2_);
671
672 dict.readIfPresent("minFaceWeight", minFaceWeight_);
673 dict.readIfPresent("minNonOrthogonality", minNonOrthogonality_);
674 dict.readIfPresent("minGradCc", minGradCc_);
675 dict.readIfPresent("minSkewness", minSkewness_);
676 dict.readIfPresent("Co1", Co1_);
677
678
679 dict.readIfPresent("P", P_);
680 dict.readIfPresent("I", I_);
681 dict.readIfPresent("D", D_);
682
683 tolerance_ = 0.001;
684 if
685 (
686 dict.readIfPresent("tolerance", tolerance_)
687 && (tolerance_ < 0 || tolerance_ > 1)
688 )
689 {
691 << "tolerance must be in the range 0 to 1. Supplied value: "
692 << tolerance_ << exit(FatalError);
693 }
694
695 Info<< type() << " " << name() << ":" << nl;
696 if (nonOrthogonality_)
697 {
698 Info<< " Including nonOrthogonality between: "
699 << minNonOrthogonality_ << " and " << maxNonOrthogonality_
700 << endl;
701 }
702 if (gradCc_)
703 {
704 Info<< " Including gradient between: "
705 << minGradCc_ << " and " << maxGradCc_ << endl;
706 }
707 if (residuals_)
708 {
709 Info<< " Including residuals" << endl;
710 }
711 if (faceWeight_)
712 {
713 Info<< " Including faceWeight between: "
714 << minFaceWeight_ << " and " << maxFaceWeight_ << endl;
715 }
716 if (skewness_)
717 {
718 Info<< " Including skewness between: "
719 << minSkewness_ << " and " << maxSkewness_ << endl;
720 }
721 if (Co_)
722 {
723 Info<< " Including Co between: "
724 << Co2_ << " and " << Co1_ << endl;
725 }
726
727 return true;
728 }
729
730 return false;
731}
732
733
735{
736 label nCellsScheme1 = 0;
737 label nCellsScheme2 = 0;
738 label nCellsBlended = 0;
739
740 calcStats(nCellsScheme1, nCellsScheme2, nCellsBlended);
741
742 if (writeToFile_)
743 {
744 writeCurrentTime(file());
745
746 file()
747 << tab << nCellsScheme1
748 << tab << nCellsScheme2
749 << tab << nCellsBlended
750 << endl;
751 }
752
753 return true;
754}
755
756
757// ************************************************************************* //
#define Log
Definition: PDRblock.C:35
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
A primitive field of type <T> with automated input and output.
Definition: IOField.H:58
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
bool typeHeaderOk(const bool checkType=true, const bool search=true, const bool verbose=true)
Read header (uses typeFilePath to find file) and check its info.
writeOption writeOpt() const noexcept
The write option.
Definition: IOobjectI.H:179
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
virtual bool read()
Re-read model coefficients if they have changed.
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:78
const word & constant() const
Return constant name.
Definition: TimePathsI.H:96
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:80
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Abstract base-class for Time/database function objects.
Intermediate class for handling field expression function objects (e.g. blendingFactor etc....
word resultName_
Name of result field.
word fieldName_
Name of field to process.
void setResultName(const word &typeName, const word &defaultArg)
Set the name of result field.
const fvMesh & mesh_
Reference to the fvMesh.
Computes the natural logarithm of an input volScalarField.
Definition: log.H:230
bool store(word &fieldName, const tmp< ObjectType > &tfield, bool cacheable=false)
Store the field in the (sub) objectRegistry under the given name.
Computes the stabilityBlendingFactor to be used by the local blended convection scheme....
virtual void writeFileHeader(Ostream &os) const
Write the file header.
virtual bool write()
Write the stabilityBlendingFactor.
virtual bool read(const dictionary &)
Read the stabilityBlendingFactor data.
const Time & time_
Reference to the time database.
Base class for writing single files from the function objects.
Definition: writeFile.H:120
bool writeToFile_
Flag to enable/disable writing to file.
Definition: writeFile.H:141
virtual OFstream & file()
Return access to the file (if only 1)
Definition: writeFile.C:233
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:290
const Type * findObject(const word &name, const bool recursive=false) const
Return const pointer to the object of the given Type.
Type * getObjectPtr(const word &name, const bool recursive=false) const
static const char *const componentNames[]
Definition: bool.H:104
static constexpr direction nComponents
Number of components in bool is 1.
Definition: bool.H:98
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
engineTime & runTime
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
OBJstream os(runTime.globalPath()/outputName)
Calculate the gradient of the given field.
#define WarningInFunction
Report a warning using Foam::Warning.
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:54
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:83
const dimensionSet dimless
Dimensionless.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
static void writeHeader(Ostream &os, const word &fieldName)
messageStream Info
Information stream (stdout output on master, null elsewhere)
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
void reduce(const List< UPstream::commsStruct > &comms, T &value, const BinaryOp &bop, const int tag, const label comm)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
errorManip< error > abort(error &err)
Definition: errorManip.H:144
uint8_t direction
Definition: direction.H:56
Type gAverage(const FieldField< Field, Type > &f)
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
error FatalError
dimensionedScalar cbrt(const dimensionedScalar &ds)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
constexpr char tab
The tab '\t' character(0x09)
Definition: Ostream.H:52
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333