35#define Log if (verbose_) Info
46 { expansionType::EXPAND_UNIFORM,
"uniform" },
47 { expansionType::EXPAND_RATIO,
"ratio" },
48 { expansionType::EXPAND_RELATIVE,
"relative"}
58inline scalar
calcGexp(
const scalar expRatio,
const label nDiv)
60 return nDiv > 1 ?
pow(expRatio, 1.0/(nDiv - 1)) : 0.0;
66 const scalar expRatio,
70 return nDiv > 1 ?
pow(expRatio, (nDiv - 1)) : 1.0;
78bool Foam::PDRblock::checkMonotonic
81 const UList<scalar>& pts
84 const label len = pts.size();
91 const scalar& minVal = pts[0];
93 for (label i=1; i < len; ++i)
99 <<
" direction do not increase monotonically" <<
nl
111 return NullObjectRef<PDRblock>();
117void Foam::PDRblock::adjustSizes()
120 sizes().x() = grid_.x().nCells();
121 sizes().y() = grid_.y().nCells();
122 sizes().z() = grid_.z().nCells();
124 if (sizes().
x() <= 0 || sizes().
y() <= 0 || sizes().z() <= 0)
134 edgeLimits_.
min() = 0;
135 edgeLimits_.
max() = 0;
140 bounds_ = bounds(grid_.x(), grid_.y(), grid_.z());
145 edgeLimits_.add(grid_.x().edgeLimits());
146 edgeLimits_.add(grid_.y().edgeLimits());
147 edgeLimits_.add(grid_.z().edgeLimits());
151void Foam::PDRblock::readGridControl
154 const dictionary&
dict,
155 const scalar scaleFactor,
156 expansionType expandType
159 gridControl& ctrl = control_[cmpt];
172 Log <<
"Reading grid control for "
180 const label nSegments = (knots.size()-1);
185 <<
"Must define at least two control points:"
187 <<
" direction" <<
nl
194 for (scalar& pt : knots)
201 checkMonotonic(cmpt, knots);
211 label nTotalDivs = 0;
212 for (
const label ndiv : divisions)
219 <<
"Negative or zero divisions:"
221 <<
" direction" <<
nl
226 if (divisions.size() != nSegments)
229 <<
"Mismatch in number of divisions and number of points:"
231 <<
" direction" <<
nl
234 else if (!nTotalDivs)
237 <<
"No grid divisions at all:"
239 <<
" direction" <<
nl
252 for (scalar& expRatio : expansion)
256 expRatio = 1.0/(-expRatio);
260 if (expansion.size() && divisions.size() != expansion.size())
263 <<
"Mismatch in number of expansion ratios and divisions:"
265 <<
" direction" <<
nl
269 if (expansion.empty())
271 expansion.resize(nSegments, scalar(1));
273 if (expandType != expansionType::EXPAND_UNIFORM)
275 expandType = expansionType::EXPAND_UNIFORM;
277 Log <<
"Warning: no 'ratios', use uniform spacing" <<
nl;
284 case expansionType::EXPAND_UNIFORM:
286 expansion = scalar(1);
290 case expansionType::EXPAND_RATIO:
296 case expansionType::EXPAND_RELATIVE:
300 auto divIter = divisions.cbegin();
302 for (scalar& expRatio : expansion)
318 gridPoint.
resize(nTotalDivs+1);
322 for (label segmenti=0; segmenti < nSegments; ++segmenti)
324 const label nDiv = divisions[segmenti];
325 const scalar expRatio = expansion[segmenti];
327 SubList<scalar> subPoint(gridPoint, nDiv+1, start);
330 subPoint[0] = knots[segmenti];
331 subPoint[nDiv] = knots[segmenti+1];
333 const scalar dist = (subPoint.last() - subPoint.first());
337 expandType == expansionType::EXPAND_UNIFORM
338 ||
equal(expRatio, 1)
341 for (label i=1; i < nDiv; ++i)
343 subPoint[i] = (subPoint[0] + (dist * i)/nDiv);
349 const scalar expFact =
calcGexp(expRatio, nDiv);
351 for (label i=1; i < nDiv; ++i)
356 + dist * (1.0 -
pow(expFact, i))/(1.0 -
pow(expFact, nDiv))
364void Foam::PDRblock::readBoundary(
const dictionary&
dict)
366 Log <<
"Reading boundary entries" <<
nl;
368 PtrList<entry> patchEntries;
372 PtrList<entry> newEntries(
dict.
lookup(
"boundary"));
373 patchEntries.transfer(newEntries);
377 patches_.resize(patchEntries.size());
380 const labelRange validFaceId(0, 6);
385 forAll(patchEntries, patchi)
387 if (!patchEntries.set(patchi))
392 const entry&
e = patchEntries[patchi];
398 <<
" in boundary section is not a valid dictionary."
402 const dictionary& patchDict =
e.dict();
404 const word& patchName =
e.keyword();
406 const word patchType(patchDict.get<word>(
"type"));
412 if (faceLabels.size() != patchFaces.size())
415 <<
"Patch: " << patchName
416 <<
": Ignoring duplicate face ids" <<
nl;
419 label nErased = patchFaces.filterKeys(validFaceId);
423 <<
"Patch: " << patchName <<
": Ignoring " << nErased
424 <<
" faces with invalid identifiers" <<
nl;
427 nErased = patchFaces.erase(handled);
431 <<
"Patch: " << patchName <<
": Ignoring " << nErased
432 <<
" faces, which were already handled" <<
nl;
436 handled += patchFaces;
441 patches_.set(patchi,
new boundaryEntry());
443 boundaryEntry& bentry = patches_[patchi];
445 bentry.name_ = patchName;
446 bentry.type_ = patchType;
448 bentry.faces_ = patchFaces.sortedToc();
451 for (
const label shapeFacei : bentry.faces_)
453 bentry.size_ += nBoundaryFaces(shapeFacei);
461 missed.erase(handled);
465 patches_.append(
new boundaryEntry());
467 boundaryEntry& bentry = patches_.last();
469 bentry.name_ =
"defaultFaces";
472 bentry.faces_ = missed.sortedToc();
475 for (
const label shapeFacei : bentry.faces_)
477 bentry.size_ += nBoundaryFaces(shapeFacei);
484 const dictionary* defaultEntry =
dict.
findDict(
"defaultPatch");
488 defaultEntry->readIfPresent(
"type", bentry.type_);
495 for (
const boundaryEntry& bentry : patches_)
497 Info<<
" patch: " << patchi
498 <<
" (size: " << bentry.size_
499 <<
" type: " << bentry.type_
500 <<
") name: " << bentry.name_
528 for (label patchi=0; patchi < 6; ++patchi)
530 patches_.
set(patchi,
new boundaryEntry());
532 boundaryEntry& bentry = patches_[patchi];
535 bentry.type_ =
"patch";
540 reset(xgrid, ygrid, zgrid);
553 verbose_(verboseOutput)
570 expansionNames_.getOrDefault
574 expansionType::EXPAND_RATIO
578 readGridControl(0,
dict.
subDict(
"x"), scaleFactor, expandType);
579 readGridControl(1,
dict.
subDict(
"y"), scaleFactor, expandType);
580 readGridControl(2,
dict.
subDict(
"z"), scaleFactor, expandType);
592 outer_.
read(*outerDictPtr);
614 checkMonotonic(cmpt, grid_[cmpt]);
621 for (boundaryEntry& bentry : patches_)
626 for (
const label shapeFacei : bentry.faces_)
628 bentry.size_ += nBoundaryFaces(shapeFacei);
643 if (!bounds_.contains(pt))
658bool Foam::PDRblock::gridIndex
665 const scalar tol = relTol * edgeLimits_.min();
670 pos[cmpt] = grid_[cmpt].findIndex(pt[cmpt], tol);
672 if (
pos[cmpt] < 0)
return false;
683 if (findCell(pt,
pos))
700 if (gridIndex(pt,
pos, relTol))
711 return grading(control_);
723 return control_[cmpt].grading();
729 <<
"Not gridControl for direction " << label(cmpt) <<
endl
Various functions to operate on Lists.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
void resize(const label len)
Adjust allocated size of list.
void clear()
Clear the list, i.e. set size to zero.
A single block x-y-z rectilinear mesh addressable as i,j,k with simplified creation....
static const PDRblock & null()
Return a PDRblock reference to a nullObject.
static const Enum< expansionType > expansionNames_
Named enumerations for the expansion type.
bool read(const dictionary &dict)
Read dictionary.
void reset(const UList< scalar > &xgrid, const UList< scalar > &ygrid, const UList< scalar > &zgrid)
Reset grid locations and mesh i-j-k sizing.
PDRblock()
Default construct, zero-size, inverted bounds etc.
Vector< gradingDescriptors > grading() const
Equivalent edge grading descriptors in (x,y,z) directions.
expansionType
The expansion type.
const T * set(const label i) const
void resize(const label newLen)
Adjust size of PtrList.
virtual bool read()
Re-read model coefficients if they have changed.
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
const point & min() const
Minimum describing the bounding box.
static const boundBox invertedBox
A large inverted boundBox: min/max == +/- ROOTVGREAT.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
bool isNullDict() const noexcept
The dictionary is actually dictionary::null (root dictionary)
dictionary * findDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX)
Find and return a sub-dictionary pointer if present.
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
ITstream & lookup(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
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
bool read(Istream &is)
Read dictionary from Istream. Discards the header.
void reset()
Reset to defaults.
List of gradingDescriptor for the sections of a block with additional IO functionality.
A simple i-j-k (row-major order) to linear addressing for a rectilinear mesh. Since the underlying me...
A class representing the concept of 1 (one) that can be used to avoid manipulating objects known to b...
static const char *const componentNames[]
static constexpr direction nComponents
Number of components in bool is 1.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
#define WarningInFunction
Report a warning using Foam::Warning.
tmp< GeometricField< Type, faPatchField, areaMesh > > ndiv(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
dimensionedScalar pos(const dimensionedScalar &ds)
scalar calcGexp(const scalar expRatio, const label nDiv)
Calculate the geometric expansion factor from the expansion ratio.
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i)
List< label > labelList
A List of labels.
Vector< label > labelVector
Vector of labels.
messageStream Info
Information stream (stdout output on master, null elsewhere)
List< scalar > scalarList
A List of scalars.
label findLower(const ListType &input, const T &val, const label start, const ComparePredicate &comp)
vector point
Point is a vector.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
bool equal(const T &s1, const T &s2)
Compare two values for equality.
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
scalar relativeToGeometricRatio(const scalar expRatio, const label nDiv)
Calculate geometric ratio from relative ratio.
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.
static const char *const typeName
The type name used in ensight case files.