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 return NullObjectRef<PDRblock>();
117 void Foam::PDRblock::addDefaultPatches()
121 for (label patchi=0; patchi < 6; ++patchi)
123 boundaryEntry& bentry = patches_.emplace_set(patchi);
126 bentry.type_ =
"patch";
128 bentry.faces_ =
labelList(one{}, patchi);
133 void Foam::PDRblock::adjustSizes()
136 sizes().x() = grid_.x().nCells();
137 sizes().y() = grid_.y().nCells();
138 sizes().z() = grid_.z().nCells();
140 if (sizes().
x() <= 0 || sizes().
y() <= 0 || sizes().z() <= 0)
150 edgeLimits_.min() = 0;
151 edgeLimits_.max() = 0;
156 bounds_ = bounds(grid_.x(), grid_.y(), grid_.z());
161 edgeLimits_.add(grid_.x().edgeLimits());
162 edgeLimits_.add(grid_.y().edgeLimits());
163 edgeLimits_.add(grid_.z().edgeLimits());
167 void Foam::PDRblock::readGridControl
170 const dictionary&
dict,
171 const scalar scaleFactor,
172 expansionType expandType
175 gridControl& ctrl = control_[cmpt];
188 Log <<
"Reading grid control for " 189 << vector::componentNames[cmpt] <<
" direction" <<
nl;
196 const label nSegments = (knots.size()-1);
201 <<
"Must define at least two control points:" 202 <<
" in " << vector::componentNames[cmpt]
203 <<
" direction" <<
nl 210 for (scalar& pt : knots)
217 checkMonotonic(cmpt, knots);
227 label nTotalDivs = 0;
228 for (
const label
ndiv : divisions)
235 <<
"Negative or zero divisions:" 236 <<
" in " << vector::componentNames[cmpt]
237 <<
" direction" <<
nl 242 if (divisions.size() != nSegments)
245 <<
"Mismatch in number of divisions and number of points:" 246 <<
" in " << vector::componentNames[cmpt]
247 <<
" direction" <<
nl 250 else if (!nTotalDivs)
253 <<
"No grid divisions at all:" 254 <<
" in " << vector::componentNames[cmpt]
255 <<
" direction" <<
nl 268 for (scalar& expRatio : expansion)
272 expRatio = 1.0/(-expRatio);
276 if (expansion.size() && divisions.size() != expansion.size())
279 <<
"Mismatch in number of expansion ratios and divisions:" 280 <<
" in " << vector::componentNames[cmpt]
281 <<
" direction" <<
nl 285 if (expansion.empty())
287 expansion.resize(nSegments, scalar(1));
289 if (expandType != expansionType::EXPAND_UNIFORM)
291 expandType = expansionType::EXPAND_UNIFORM;
293 Log <<
"Warning: no 'ratios', use uniform spacing" <<
nl;
300 case expansionType::EXPAND_UNIFORM:
302 expansion = scalar(1);
306 case expansionType::EXPAND_RATIO:
312 case expansionType::EXPAND_RELATIVE:
316 auto divIter = divisions.cbegin();
318 for (scalar& expRatio : expansion)
334 gridPoint.
resize(nTotalDivs+1);
338 for (label segmenti=0; segmenti < nSegments; ++segmenti)
340 const label nDiv = divisions[segmenti];
341 const scalar expRatio = expansion[segmenti];
343 SubList<scalar> subPoint(gridPoint, nDiv+1, start);
346 subPoint[0] = knots[segmenti];
347 subPoint[nDiv] = knots[segmenti+1];
349 const scalar dist = (subPoint.last() - subPoint.first());
353 expandType == expansionType::EXPAND_UNIFORM
354 ||
equal(expRatio, 1)
357 for (label i=1; i < nDiv; ++i)
359 subPoint[i] = (subPoint[0] + (dist * i)/nDiv);
365 const scalar expFact =
calcGexp(expRatio, nDiv);
367 for (label i=1; i < nDiv; ++i)
372 + dist * (1.0 -
pow(expFact, i))/(1.0 -
pow(expFact, nDiv))
380 void Foam::PDRblock::readBoundary(
const dictionary&
dict)
382 Log <<
"Reading boundary entries" <<
nl;
384 PtrList<entry> patchEntries;
388 PtrList<entry> newEntries(
dict.
lookup(
"boundary"));
389 patchEntries.transfer(newEntries);
393 patches_.resize(patchEntries.size());
396 const labelRange validFaceId(0, 6);
401 forAll(patchEntries, patchi)
403 if (!patchEntries.set(patchi))
408 const entry&
e = patchEntries[patchi];
414 <<
" in boundary section is not a valid dictionary." 418 const dictionary& patchDict =
e.dict();
420 const word& patchName =
e.keyword();
422 const word patchType(patchDict.get<word>(
"type"));
431 <<
"Patch: " << patchName
432 <<
": Ignoring duplicate face ids" <<
nl;
435 label nErased = patchFaces.filterKeys(validFaceId);
439 <<
"Patch: " << patchName <<
": Ignoring " << nErased
440 <<
" faces with invalid identifiers" <<
nl;
443 nErased = patchFaces.erase(handled);
447 <<
"Patch: " << patchName <<
": Ignoring " << nErased
448 <<
" faces, which were already handled" <<
nl;
452 handled += patchFaces;
457 boundaryEntry& bentry = patches_.emplace_set(patchi);
459 bentry.name_ = patchName;
460 bentry.type_ = patchType;
462 bentry.faces_ = patchFaces.sortedToc();
465 for (
const label shapeFacei : bentry.faces_)
467 bentry.size_ += nBoundaryFaces(shapeFacei);
475 missed.erase(handled);
479 boundaryEntry& bentry = patches_.emplace_back();
481 bentry.name_ =
"defaultFaces";
482 bentry.type_ = emptyPolyPatch::typeName;
484 bentry.faces_ = missed.sortedToc();
487 for (
const label shapeFacei : bentry.faces_)
489 bentry.size_ += nBoundaryFaces(shapeFacei);
496 const dictionary* defaultEntry =
dict.
findDict(
"defaultPatch");
500 defaultEntry->readIfPresent(
"type", bentry.type_);
507 for (
const boundaryEntry& bentry : patches_)
509 Info<<
" patch: " << patchi
510 <<
" (size: " << bentry.size_
511 <<
" type: " << bentry.type_
512 <<
") name: " << bentry.name_
539 reset(xgrid, ygrid, zgrid);
561 verbose_(verboseOutput)
563 if (!
dict.isNullDict())
578 expansionNames_.getOrDefault
582 expansionType::EXPAND_RATIO
586 readGridControl(0,
dict.
subDict(
"x"), scaleFactor, expandType);
587 readGridControl(1,
dict.
subDict(
"y"), scaleFactor, expandType);
588 readGridControl(2,
dict.
subDict(
"z"), scaleFactor, expandType);
600 outer_.
read(*outerDictPtr);
620 for (
direction cmpt=0; cmpt < vector::nComponents; ++cmpt)
622 checkMonotonic(cmpt, grid_[cmpt]);
629 for (boundaryEntry& bentry : patches_)
634 for (
const label shapeFacei : bentry.faces_)
636 bentry.size_ += nBoundaryFaces(shapeFacei);
644 grid_.x().reset(box.
min().x(), box.
max().x(), nCells.
x());
645 grid_.y().reset(box.
min().y(), box.
max().y(), nCells.
y());
646 grid_.z().reset(box.
min().z(), box.
max().z(), nCells.
z());
651 for (boundaryEntry& bentry : patches_)
656 for (
const label shapeFacei : bentry.faces_)
658 bentry.size_ += nBoundaryFaces(shapeFacei);
673 if (!bounds_.contains(pt))
688 bool Foam::PDRblock::gridIndex
695 const scalar tol = relTol * edgeLimits_.min();
700 pos[cmpt] = grid_[cmpt].findIndex(pt[cmpt], tol);
702 if (
pos[cmpt] < 0)
return false;
713 if (findCell(pt,
pos))
730 if (gridIndex(pt,
pos, relTol))
741 return grading(control_);
753 return control_[cmpt].grading();
759 <<
"Not gridControl for direction " << label(cmpt) <<
endl 764 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.
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.
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 INVALID.
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i)
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.
static const PDRblock & null()
Return a PDRblock reference to a nullObject.
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 a sub-dictionary) otherwise return nullptr...
IOerror FatalIOError
Error stream (stdout output on all processes), with additional 'FOAM FATAL IO ERROR' header text and ...