35 #define Log if (verbose_) Info 46 { expansionType::EXPAND_UNIFORM,
"uniform" },
47 { expansionType::EXPAND_RATIO,
"ratio" },
48 { expansionType::EXPAND_RELATIVE,
"relative"}
58 inline 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;
78 bool 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)
98 <<
"Points in " << vector::componentNames[cmpt]
99 <<
" direction do not increase monotonically" <<
nl 111 void Foam::PDRblock::addDefaultPatches()
115 for (label patchi=0; patchi < 6; ++patchi)
117 boundaryEntry& bentry = patches_.emplace_set(patchi);
120 bentry.type_ =
"patch";
122 bentry.faces_ =
labelList(one{}, patchi);
127 void Foam::PDRblock::adjustSizes()
130 sizes().x() = grid_.x().nCells();
131 sizes().y() = grid_.y().nCells();
132 sizes().z() = grid_.z().nCells();
134 if (sizes().
x() <= 0 || sizes().
y() <= 0 || sizes().z() <= 0)
144 edgeLimits_.min() = 0;
145 edgeLimits_.max() = 0;
150 bounds_ = bounds(grid_.x(), grid_.y(), grid_.z());
155 edgeLimits_.add(grid_.x().edgeLimits());
156 edgeLimits_.add(grid_.y().edgeLimits());
157 edgeLimits_.add(grid_.z().edgeLimits());
161 void Foam::PDRblock::readGridControl
164 const dictionary&
dict,
165 const scalar scaleFactor,
166 expansionType expandType
169 gridControl& ctrl = control_[cmpt];
182 Log <<
"Reading grid control for " 183 << vector::componentNames[cmpt] <<
" direction" <<
nl;
190 const label nSegments = (knots.size()-1);
195 <<
"Must define at least two control points:" 196 <<
" in " << vector::componentNames[cmpt]
197 <<
" direction" <<
nl 204 for (scalar& pt : knots)
211 checkMonotonic(cmpt, knots);
221 label nTotalDivs = 0;
222 for (
const label
ndiv : divisions)
229 <<
"Negative or zero divisions:" 230 <<
" in " << vector::componentNames[cmpt]
231 <<
" direction" <<
nl 236 if (divisions.size() != nSegments)
239 <<
"Mismatch in number of divisions and number of points:" 240 <<
" in " << vector::componentNames[cmpt]
241 <<
" direction" <<
nl 244 else if (!nTotalDivs)
247 <<
"No grid divisions at all:" 248 <<
" in " << vector::componentNames[cmpt]
249 <<
" direction" <<
nl 262 for (scalar& expRatio : expansion)
266 expRatio = 1.0/(-expRatio);
270 if (expansion.size() && divisions.size() != expansion.size())
273 <<
"Mismatch in number of expansion ratios and divisions:" 274 <<
" in " << vector::componentNames[cmpt]
275 <<
" direction" <<
nl 279 if (expansion.empty())
281 expansion.resize(nSegments, scalar(1));
283 if (expandType != expansionType::EXPAND_UNIFORM)
285 expandType = expansionType::EXPAND_UNIFORM;
287 Log <<
"Warning: no 'ratios', use uniform spacing" <<
nl;
294 case expansionType::EXPAND_UNIFORM:
296 expansion = scalar(1);
300 case expansionType::EXPAND_RATIO:
306 case expansionType::EXPAND_RELATIVE:
310 auto divIter = divisions.cbegin();
312 for (scalar& expRatio : expansion)
328 gridPoint.
resize(nTotalDivs+1);
332 for (label segmenti=0; segmenti < nSegments; ++segmenti)
334 const label nDiv = divisions[segmenti];
335 const scalar expRatio = expansion[segmenti];
337 SubList<scalar> subPoint(gridPoint, nDiv+1, start);
340 subPoint[0] = knots[segmenti];
341 subPoint[nDiv] = knots[segmenti+1];
343 const scalar dist = (subPoint.last() - subPoint.first());
347 expandType == expansionType::EXPAND_UNIFORM
348 ||
equal(expRatio, 1)
351 for (label i=1; i < nDiv; ++i)
353 subPoint[i] = (subPoint[0] + (dist * i)/nDiv);
359 const scalar expFact =
calcGexp(expRatio, nDiv);
361 for (label i=1; i < nDiv; ++i)
366 + dist * (1.0 -
pow(expFact, i))/(1.0 -
pow(expFact, nDiv))
374 void Foam::PDRblock::readBoundary(
const dictionary&
dict)
376 Log <<
"Reading boundary entries" <<
nl;
378 PtrList<entry> patchEntries;
382 PtrList<entry> newEntries(
dict.
lookup(
"boundary"));
383 patchEntries.transfer(newEntries);
387 patches_.resize(patchEntries.size());
390 const labelRange validFaceId(0, 6);
395 forAll(patchEntries, patchi)
397 if (!patchEntries.set(patchi))
402 const entry&
e = patchEntries[patchi];
408 <<
" in boundary section is not a valid dictionary." 414 const word& patchName =
e.keyword();
416 const word patchType(patchDict.get<word>(
"type"));
425 <<
"Patch: " << patchName
426 <<
": Ignoring duplicate face ids" <<
nl;
429 label nErased = patchFaces.filterKeys(validFaceId);
433 <<
"Patch: " << patchName <<
": Ignoring " << nErased
434 <<
" faces with invalid identifiers" <<
nl;
437 nErased = patchFaces.erase(handled);
441 <<
"Patch: " << patchName <<
": Ignoring " << nErased
442 <<
" faces, which were already handled" <<
nl;
446 handled += patchFaces;
451 boundaryEntry& bentry = patches_.emplace_set(patchi);
453 bentry.name_ = patchName;
454 bentry.type_ = patchType;
456 bentry.faces_ = patchFaces.sortedToc();
459 for (
const label shapeFacei : bentry.faces_)
461 bentry.size_ += nBoundaryFaces(shapeFacei);
469 missed.erase(handled);
473 boundaryEntry& bentry = patches_.emplace_back();
475 bentry.name_ =
"defaultFaces";
476 bentry.type_ = emptyPolyPatch::typeName;
478 bentry.faces_ = missed.sortedToc();
481 for (
const label shapeFacei : bentry.faces_)
483 bentry.size_ += nBoundaryFaces(shapeFacei);
494 defaultEntry->readIfPresent(
"type", bentry.type_);
501 for (
const boundaryEntry& bentry : patches_)
503 Info<<
" patch: " << patchi
504 <<
" (size: " << bentry.size_
505 <<
" type: " << bentry.type_
506 <<
") name: " << bentry.name_
533 reset(xgrid, ygrid, zgrid);
555 verbose_(verboseOutput)
557 if (!
dict.isNullDict())
572 expansionNames_.getOrDefault
576 expansionType::EXPAND_RATIO
580 readGridControl(0,
dict.
subDict(
"x"), scaleFactor, expandType);
581 readGridControl(1,
dict.
subDict(
"y"), scaleFactor, expandType);
582 readGridControl(2,
dict.
subDict(
"z"), scaleFactor, expandType);
594 outer_.
read(*outerDictPtr);
616 checkMonotonic(cmpt, grid_[cmpt]);
623 for (boundaryEntry& bentry : patches_)
628 for (
const label shapeFacei : bentry.faces_)
630 bentry.size_ += nBoundaryFaces(shapeFacei);
638 grid_.x().reset(box.
min().x(), box.
max().x(), nCells.
x());
639 grid_.y().reset(box.
min().y(), box.
max().y(), nCells.
y());
640 grid_.z().reset(box.
min().z(), box.
max().z(), nCells.
z());
645 for (boundaryEntry& bentry : patches_)
650 for (
const label shapeFacei : bentry.faces_)
652 bentry.size_ += nBoundaryFaces(shapeFacei);
667 if (!bounds_.contains(pt))
682 bool Foam::PDRblock::gridIndex
689 const scalar tol = relTol * edgeLimits_.min();
694 pos[cmpt] = grid_[cmpt].findIndex(pt[cmpt], tol);
696 if (
pos[cmpt] < 0)
return false;
707 if (findCell(pt,
pos))
724 if (gridIndex(pt,
pos, relTol))
735 return grading(control_);
747 return control_[cmpt].grading();
753 <<
"Not gridControl for direction " << label(cmpt) <<
endl 758 return gradingDescriptors();
List< scalar > scalarList
List of scalar.
ITstream & lookup(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return an entry data stream. FatalIOError if not found, or not a stream. ...
void size(const label n)
Older name for setAddressableSize.
scalar relativeToGeometricRatio(const scalar expRatio, const label nDiv)
Calculate geometric ratio from relative ratio.
label findLower(const ListType &input, const T &val, const label start, const ComparePredicate &comp)
Binary search to find the index of the last element in a sorted list that is less than value...
errorManipArg< error, int > exit(error &err, const int errNo=1)
bool read(Istream &is)
Read dictionary from Istream (discards the header). Reads entries until EOF or when the first token i...
void resize(const label len)
Adjust allocated size of list.
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
bool equal(const T &a, const T &b)
Compare two values for equality.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
constexpr char nl
The newline '\n' character (0x0a)
label nCells() const
The number of mesh cells (nx*ny*nz) in the i-j-k mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
A bounding box defined in terms of min/max extrema points.
const Cmpt & y() const noexcept
Access to the vector y component.
::Foam::direction nComponents(const expressions::valueTypeCode) noexcept
The number of components associated with given valueTypeCode.
const point & min() const noexcept
Minimum describing the bounding box.
expansionType
The expansion type.
labelList faceLabels(nFaceLabels)
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
Various functions to operate on Lists.
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find an entry (const access) with the given keyword.
#define forAll(list, i)
Loop across all elements in list.
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, IOobjectOption::readOption readOpt=IOobjectOption::MUST_READ) const
Find entry and assign to T val. FatalIOError if it is found and the number of tokens is incorrect...
bool read(const dictionary &dict)
Read dictionary.
dimensionedScalar pos(const dimensionedScalar &ds)
const point & max() const noexcept
Maximum describing the bounding box.
PDRblock()
Default construct, zero-size, inverted bounds etc.
const dimensionedScalar e
Elementary charge.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i), works like std::iota() but returning a...
void clear()
Clear the list, i.e. set size to zero.
static const Enum< expansionType > expansionNames_
Named enumerations for the expansion type.
A single block x-y-z rectilinear mesh addressable as i,j,k with simplified creation. Some of the input is similar to blockMeshDict, but since this specialization is for a single-block that is aligned with the x-y-z directions, it provides a different means of specifying the mesh.
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
static constexpr direction nComponents
Number of components in this vector space.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
Vector< gradingDescriptors > grading() const
Equivalent edge grading descriptors in (x,y,z) directions.
const Cmpt & x() const noexcept
Access to the vector x component.
scalar calcGexp(const scalar expRatio, const label nDiv)
Calculate the geometric expansion factor from the expansion ratio.
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
Find an entry if present, and assign to T val. FatalIOError if it is found and the number of tokens i...
A simple i-j-k (row-major order) to linear addressing for a rectilinear mesh. Since the underlying me...
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void clear()
Reset to (0,0,0) sizing.
const Cmpt & z() const noexcept
Access to the vector z component.
vector point
Point is a vector.
#define WarningInFunction
Report a warning using Foam::Warning.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
PtrList< volScalarField > & Y
void reset(const UList< scalar > &xgrid, const UList< scalar > &ygrid, const UList< scalar > &zgrid)
Reset grid locations and mesh i-j-k sizing.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Vector< label > labelVector
Vector of labels.
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a T, or return the given default value. FatalIOError if it is found and the number of...
List< label > labelList
A List of labels.
tmp< GeometricField< Type, faPatchField, areaMesh > > ndiv(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
List of gradingDescriptor for the sections of a block with additional IO functionality.
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
const dictionary * findDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary pointer if present (and it is a dictionary) otherwise return nullptr...
IOerror FatalIOError
Error stream (stdout output on all processes), with additional 'FOAM FATAL IO ERROR' header text and ...