42Foam::layerParameters::thicknessModelTypeNames_
44 { thicknessModelType::FIRST_AND_TOTAL,
"firstAndOverall" },
45 { thicknessModelType::FIRST_AND_EXPANSION,
"firstAndExpansion" },
46 { thicknessModelType::FINAL_AND_TOTAL,
"finalAndOverall" },
47 { thicknessModelType::FINAL_AND_EXPANSION,
"finalAndExpansion" },
48 { thicknessModelType::TOTAL_AND_EXPANSION,
"overallAndExpansion" },
49 { thicknessModelType::FIRST_AND_RELATIVE_FINAL,
"firstAndRelativeFinal" },
52const Foam::scalar Foam::layerParameters::defaultConcaveAngle = 90;
57Foam::scalar Foam::layerParameters::layerExpansionRatio
60 const scalar totalOverFirst
69 const scalar tol = 1
e-8;
71 if (
mag(
n-totalOverFirst) < tol)
80 if (totalOverFirst <
n)
83 maxR =
pow(totalOverFirst/
n, scalar(1)/(
n-1));
87 minR =
pow(totalOverFirst/
n, scalar(1)/(
n-1));
88 maxR = totalOverFirst/(
n - 1);
92 scalar r = 0.5*(minR + maxR);
96 const scalar prevr = r;
98 const scalar fx =
pow(r,
n) - totalOverFirst*r - (1 - totalOverFirst);
99 const scalar dfx =
n*
pow(r,
n - 1) - totalOverFirst;
102 if (
mag(r - prevr) < tol)
111void Foam::layerParameters::readLayerParameters
114 const dictionary&
dict,
115 const thicknessModelType& spec,
116 scalar& firstLayerThickness,
117 scalar& finalLayerThickness,
119 scalar& expansionRatio
125 case FIRST_AND_TOTAL:
128 Info<<
"Layer specification as" <<
nl
129 <<
"- first layer thickness ('firstLayerThickness')" <<
nl
130 <<
"- overall thickness ('thickness')" <<
endl;
132 firstLayerThickness =
dict.
get<scalar>(
"firstLayerThickness");
133 thickness =
dict.
get<scalar>(
"thickness");
136 case FIRST_AND_EXPANSION:
139 Info<<
"Layer specification as" <<
nl
140 <<
"- first layer thickness ('firstLayerThickness')" <<
nl
141 <<
"- expansion ratio ('expansionRatio')" <<
endl;
143 firstLayerThickness =
dict.
get<scalar>(
"firstLayerThickness");
144 expansionRatio =
dict.
get<scalar>(
"expansionRatio");
147 case FINAL_AND_TOTAL:
150 Info<<
"Layer specification as" <<
nl
151 <<
"- final layer thickness ('finalLayerThickness')" <<
nl
152 <<
"- overall thickness ('thickness')" <<
endl;
154 finalLayerThickness =
dict.
get<scalar>(
"finalLayerThickness");
155 thickness =
dict.
get<scalar>(
"thickness");
158 case FINAL_AND_EXPANSION:
161 Info<<
"Layer specification as" <<
nl
162 <<
"- final layer thickness ('finalLayerThickness')" <<
nl
163 <<
"- expansion ratio ('expansionRatio')" <<
endl;
165 finalLayerThickness =
dict.
get<scalar>(
"finalLayerThickness");
166 expansionRatio =
dict.
get<scalar>(
"expansionRatio");
169 case TOTAL_AND_EXPANSION:
172 Info<<
"Layer specification as" <<
nl
173 <<
"- overall thickness ('thickness')" <<
nl
174 <<
"- expansion ratio ('expansionRatio')" <<
endl;
176 thickness =
dict.
get<scalar>(
"thickness");
177 expansionRatio =
dict.
get<scalar>(
"expansionRatio");
180 case FIRST_AND_RELATIVE_FINAL:
183 Info<<
"Layer specification as" <<
nl
184 <<
"- absolute first layer thickness"
185 <<
" ('firstLayerThickness')"
187 <<
"- and final layer thickness"
188 <<
" ('finalLayerThickness')" <<
nl
191 firstLayerThickness =
dict.
get<scalar>(
"firstLayerThickness");
192 finalLayerThickness =
dict.
get<scalar>(
"finalLayerThickness");
203void Foam::layerParameters::calculateLayerParameters
205 const thicknessModelType& spec,
207 scalar& firstThickness,
208 scalar& finalThickness,
210 scalar& expansionRatio
216 case FIRST_AND_TOTAL:
217 expansionRatio = layerExpansionRatio
228 *finalLayerThicknessRatio(nLayers, expansionRatio);
232 case FIRST_AND_EXPANSION:
233 thickness = layerThickness
244 *finalLayerThicknessRatio(nLayers, expansionRatio);
248 case FINAL_AND_TOTAL:
249 firstThickness = firstLayerThickness
258 expansionRatio = layerExpansionRatio
270 case FINAL_AND_EXPANSION:
271 firstThickness = firstLayerThickness
280 thickness = layerThickness
291 case TOTAL_AND_EXPANSION:
292 firstThickness = firstLayerThickness
303 *finalLayerThicknessRatio(nLayers, expansionRatio);
307 case FIRST_AND_RELATIVE_FINAL:
308 thickness = layerThickness
317 expansionRatio = layerExpansionRatio
365 mergePatchFacesAngle_
367 dict.getOrDefault<scalar>
369 "mergePatchFacesAngle",
375 dict.getOrDefault<scalar>(
"concaveAngle", defaultConcaveAngle)
378 maxFaceThicknessRatio_
382 nBufferCellsNoExtrude_
388 additionalReporting_(
dict.getOrDefault(
"additionalReporting", false)),
397 nOuterIter_(
dict.getOrDefault<scalar>(
"nOuterIter", 1))
404 layerModels_ = thicknessModelTypeNames_[spec];
411 bool haveFirst =
dict.
found(
"firstLayerThickness");
416 bool haveFinal =
dict.
found(
"finalLayerThickness");
421 bool haveTotal =
dict.
found(
"thickness");
426 bool haveExp =
dict.
found(
"expansionRatio");
432 if (nSpec == 2 && haveFirst && haveTotal)
438 else if (nSpec == 2 && haveFirst && haveExp)
444 else if (nSpec == 2 && haveFinal && haveTotal)
450 else if (nSpec == 2 && haveFinal && haveExp)
456 else if (nSpec == 2 && haveTotal && haveExp)
462 else if (nSpec == 2 && haveFirst && haveFinal)
471 <<
"Over- or underspecified layer thickness."
472 <<
" Please specify" <<
nl
473 <<
" first layer thickness ('firstLayerThickness')"
474 <<
" and overall thickness ('thickness') or" <<
nl
475 <<
" first layer thickness ('firstLayerThickness')"
476 <<
" and expansion ratio ('expansionRatio') or" <<
nl
477 <<
" final layer thickness ('finalLayerThickness')"
478 <<
" and expansion ratio ('expansionRatio') or" <<
nl
479 <<
" final layer thickness ('finalLayerThickness')"
480 <<
" and overall thickness ('thickness') or" <<
nl
481 <<
" overall thickness ('thickness')"
482 <<
" and expansion ratio ('expansionRatio'"
489 scalar firstThickness;
490 scalar finalThickness;
503 firstLayerThickness_ = firstThickness;
504 finalLayerThickness_ = finalThickness;
510 if (nLayerIter_ < 0 || nRelaxedIter_ < 0)
513 <<
"Layer iterations should be >= 0." <<
nl
514 <<
"nLayerIter:" << nLayerIter_
515 <<
" nRelaxedIter:" << nRelaxedIter_
527 for (
const entry& dEntry : layersDict)
531 const keyType& key = dEntry.keyword();
539 if (patchIDs.size() == 0)
542 <<
"Layer specification for " << key
543 <<
" does not match any patch." <<
endl
548 for (
const label patchi : patchIDs)
551 layerDict.
get<label>(
"nSurfaceLayers");
558 layerModels_[patchi] = thicknessModelTypeNames_[spec];
563 layerModels_[patchi],
564 firstLayerThickness_[patchi],
565 finalLayerThickness_[patchi],
567 expansionRatio_[patchi]
569 minThickness_[patchi] =
570 layerDict.
get<scalar>(
"minThickness");
575 switch (layerModels_[patchi])
580 "firstLayerThickness",
581 firstLayerThickness_[patchi]
593 "firstLayerThickness",
594 firstLayerThickness_[patchi]
599 expansionRatio_[patchi]
606 "finalLayerThickness",
607 finalLayerThickness_[patchi]
619 "finalLayerThickness",
620 finalLayerThickness_[patchi]
625 expansionRatio_[patchi]
638 expansionRatio_[patchi]
645 "firstLayerThickness",
646 firstLayerThickness_[patchi]
650 "finalLayerThickness",
651 finalLayerThickness_[patchi]
664 minThickness_[patchi]
671 relativeSizes_[patchi]
679 forAll(numLayers_, patchi)
682 calculateLayerParameters
684 layerModels_[patchi],
686 firstLayerThickness_[patchi],
687 finalLayerThickness_[patchi],
689 expansionRatio_[patchi]
700 const scalar layerThickness,
701 const scalar expansionRatio,
703 const label layerStart,
704 const label layerSize
707 if (layerSize == 0 || nLayers == 0)
711 else if (layerSize > nLayers || layerStart >= nLayers)
714 <<
" overall nLayers:" << nLayers
715 <<
" slice nLayers:" << layerSize
716 <<
" slice start:" << layerStart
720 else if (
mag(expansionRatio-1) < SMALL)
722 return layerThickness*layerSize/nLayers;
726 const scalar firstLayerThickness =
727 finalLayerThicknessRatio(nLayers, expansionRatio)
729 /
pow(expansionRatio, nLayers-1);
732 const scalar startThickness =
734 *
pow(expansionRatio, layerStart);
737 const scalar thickness =
739 *(1.0 -
pow(expansionRatio, layerSize))
740 /(1.0 - expansionRatio);
751 const scalar firstLayerThickness,
752 const scalar finalLayerThickness,
753 const scalar totalThickness,
754 const scalar expansionRatio
759 case FIRST_AND_TOTAL:
760 case FINAL_AND_TOTAL:
761 case TOTAL_AND_EXPANSION:
763 return totalThickness;
767 case FIRST_AND_EXPANSION:
769 if (
mag(expansionRatio-1) < SMALL)
771 return firstLayerThickness * nLayers;
775 return firstLayerThickness
776 *(1.0 -
pow(expansionRatio, nLayers))
777 /(1.0 - expansionRatio);
782 case FINAL_AND_EXPANSION:
784 if (
mag(expansionRatio-1) < SMALL)
786 return finalLayerThickness * nLayers;
790 scalar invExpansion = 1.0 / expansionRatio;
791 return finalLayerThickness
792 *(1.0 -
pow(invExpansion, nLayers))
793 /(1.0 - invExpansion);
798 case FIRST_AND_RELATIVE_FINAL:
800 if (
mag(expansionRatio-1) < SMALL)
802 return firstLayerThickness * nLayers;
806 scalar ratio = layerExpansionRatio
816 if (
mag(ratio-1) < SMALL)
818 return firstLayerThickness * nLayers;
822 return firstLayerThickness *
823 (1.0 -
pow(ratio, nLayers))
840Foam::scalar Foam::layerParameters::layerExpansionRatio
844 const scalar firstLayerThickness,
845 const scalar finalLayerThickness,
846 const scalar totalThickness,
847 const scalar expansionRatio
852 case FIRST_AND_EXPANSION:
853 case FINAL_AND_EXPANSION:
854 case TOTAL_AND_EXPANSION:
856 return expansionRatio;
860 case FIRST_AND_TOTAL:
862 if (firstLayerThickness < SMALL)
869 return layerExpansionRatio
872 totalThickness/firstLayerThickness
878 case FINAL_AND_TOTAL:
880 if (finalLayerThickness < SMALL)
889 / layerExpansionRatio
892 totalThickness/finalLayerThickness
898 case FIRST_AND_RELATIVE_FINAL:
900 if (firstLayerThickness < SMALL || nLayers <= 1)
910 finalLayerThickness/firstLayerThickness,
931 const scalar firstLayerThickness,
932 const scalar finalLayerThickness,
933 const scalar totalThickness,
934 const scalar expansionRatio
939 case FIRST_AND_EXPANSION:
940 case FIRST_AND_TOTAL:
941 case FIRST_AND_RELATIVE_FINAL:
943 return firstLayerThickness;
946 case FINAL_AND_EXPANSION:
948 if (expansionRatio < SMALL)
955 return finalLayerThickness*
pow(1.0/expansionRatio, nLayers-1);
960 case FINAL_AND_TOTAL:
962 scalar r = layerExpansionRatio
971 return finalLayerThickness/
pow(r, nLayers-1);
975 case TOTAL_AND_EXPANSION:
977 scalar r = finalLayerThicknessRatio
982 scalar finalThickness = r*totalThickness;
983 return finalThickness/
pow(expansionRatio, nLayers-1);
1000 const scalar expansionRatio
1005 if (
mag(expansionRatio-1) < SMALL)
1012 pow(expansionRatio, nLayers - 1)
1013 *(1.0 - expansionRatio)
1014 /(1.0 -
pow(expansionRatio, nLayers));
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface....
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
A keyword and a list of tokens is an 'entry'.
A class for handling keywords in dictionaries.
Simple container to keep together layer specific information.
const scalarField & expansionRatio() const
const scalarField & thickness() const
Wanted overall thickness of all layers.
thicknessModelType
Enumeration defining the layer specification:
@ FIRST_AND_RELATIVE_FINAL
static scalar layerThickness(const thicknessModelType, const label nLayers, const scalar firstLayerThickness, const scalar finalLayerThickness, const scalar totalThickness, const scalar expansionRatio)
Determine overall thickness. Uses two of the four parameters.
const scalarField & firstLayerThickness() const
Wanted thickness of the layer nearest to the wall.
static scalar finalLayerThicknessRatio(const label nLayers, const scalar expansionRatio)
Determine ratio of final layer thickness to.
const dictionary & dict() const
Helper class which maintains intersections of (changing) mesh with (static) surfaces.
static const dictionary & subDict(const dictionary &dict, const word &keyword, const bool noExit, enum keyType::option matchOpt=keyType::REGEX)
Wrapper around dictionary::subDict which does not exit.
A class representing the concept of 1 (one) that can be used to avoid manipulating objects known to b...
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
A wordRe is a Foam::word, but can contain a regular expression for matching words or strings.
A List of wordRe with additional matching capabilities.
A class for handling words, derived from Foam::string.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
#define IOWarningInFunction(ios)
Report an IO warning using Foam::Warning.
messageStream Info
Information stream (stdout output on master, null elsewhere)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
static constexpr int maxIters
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
errorManipArg< error, int > exit(error &err, const int errNo=1)
constexpr char nl
The newline '\n' character (0x0a)
#define forAll(list, i)
Loop across all elements in list.
Unit conversion functions.