Go to the documentation of this file.
44 Foam::blockMesh::strategyNames_
46 { mergeStrategy::MERGE_TOPOLOGY,
"topology" },
47 { mergeStrategy::MERGE_POINTS,
"points" },
97 if ((val > 0) && !
equal(val, 1))
109 bool nonUnity =
false;
112 if (scale[cmpt] <= 0)
116 else if (!
equal(scale[cmpt], 1))
124 if (
equal(scale.x(), scale.y()) &&
equal(scale.x(), scale.z()))
146 bool Foam::blockMesh::readPointTransforms(
const dictionary&
dict)
148 transformType_ = transformTypes::NO_TRANSFORM;
151 if (
const dictionary* dictptr =
dict.findDict(
"transform"))
153 transform_ = coordSystem::cartesian(*dictptr);
155 constexpr scalar tol = VSMALL;
162 || (
mag(o.y()) > tol)
163 || (
mag(o.z()) > tol)
166 transformType_ |= transformTypes::TRANSLATION;
170 const tensor& r = transform_.
R();
173 (
mag(r.xx() - 1) > tol)
174 || (
mag(r.xy()) > tol)
175 || (
mag(r.xz()) > tol)
177 || (
mag(r.yx()) > tol)
178 || (
mag(r.yy() - 1) > tol)
179 || (
mag(r.yz()) > tol)
181 || (
mag(r.zx()) > tol)
182 || (
mag(r.zy()) > tol)
183 || (
mag(r.zz() - 1) > tol)
186 transformType_ |= transformTypes::ROTATION;
207 transformType_ |= transformTypes::PRESCALING;
209 else if (scaleType == 3)
211 transformType_ |= transformTypes::PRESCALING3;
226 {{
"convertToMeters", 1012}},
234 transformType_ |= transformTypes::SCALING;
236 else if (scaleType == 3)
238 transformType_ |= transformTypes::SCALING3;
242 return bool(transformType_);
248 Foam::blockMesh::blockMesh
258 checkFaceCorrespondence_
260 meshDict_.getOrDefault(
"checkFaceCorrespondence",
true)
262 mergeStrategy_(strategy),
263 transformType_(transformTypes::NO_TRANSFORM),
269 meshDict_.time().constant(),
275 meshDict_.found(
"geometry")
276 ? meshDict_.subDict(
"geometry")
282 meshDict_.lookup(
"vertices"),
289 topologyPtr_(createTopology(meshDict_,
regionName))
291 if (mergeStrategy_ == mergeStrategy::DEFAULT_MERGE)
293 strategyNames_.readIfPresent(
"mergeType", meshDict_, mergeStrategy_);
299 mergeStrategy_ == mergeStrategy::DEFAULT_MERGE
304 <<
"Detected collapsed blocks "
305 <<
"- using merge points instead of merge topology" <<
nl
308 mergeStrategy_ = mergeStrategy::MERGE_POINTS;
312 if (mergeStrategy_ == mergeStrategy::MERGE_POINTS)
315 calcGeometricalMerge();
320 calcTopologicalMerge();
329 return bool(topologyPtr_);
356 if (applyTransform && hasPointTransforms())
360 inplacePointTransforms(tpts.ref());
371 const polyPatchList& topoPatches = topology().boundaryMesh();
375 forAll(topoPatches, patchi)
378 topoPatches[patchi].write(
os);
410 if (patches_.empty())
421 return topology().boundaryMesh().names();
443 for (
const block& blk : blocks)
445 if (blk.zoneName().size())
457 return bool(transformType_);
468 if (transformType_ & transformTypes::PRESCALING)
472 p *= prescaling_.x();
475 else if (transformType_ & transformTypes::PRESCALING3)
483 if (transformType_ & transformTypes::ROTATION)
485 const tensor rot(transform_.R());
487 if (transformType_ & transformTypes::TRANSLATION)
489 const point origin(transform_.origin());
504 else if (transformType_ & transformTypes::TRANSLATION)
506 const point origin(transform_.origin());
514 if (transformType_ & transformTypes::SCALING)
521 else if (transformType_ & transformTypes::SCALING3)
536 if (hasPointTransforms())
540 inplacePointTransforms(tpts.ref());
A keyword and a list of tokens is an 'entry'.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
PtrList< dictionary > patchDicts() const
Patch information from the topology mesh.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Creates a single block of cells from point coordinates, numbers of cells in each direction and an exp...
bool valid() const noexcept
True if the blockMesh topology exists.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
const faceListList & patches() const
Return the patch face lists.
A class for handling words, derived from Foam::string.
static Vector< scalar > uniform(const scalar &s)
Return a VectorSpace with all elements = s.
virtual ITstream & stream() const =0
Return token stream, if entry is a primitive entry.
dimensioned< Type > cmptMultiply(const dimensioned< Type > &, const dimensioned< Type > &)
virtual void clear()
Reset origin and rotation to an identity coordinateSystem.
int verbose() const noexcept
Output verbosity level.
A class for managing temporary objects.
mergeStrategy
The block merging strategy.
const pointField & vertices() const noexcept
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
int readScaling(const entry *eptr, vector &scale)
#define forAll(list, i)
Loop across all elements in list.
PtrList< dictionary > patchDicts
bool isNumber() const noexcept
Token is LABEL, FLOAT or DOUBLE.
virtual const point & origin() const
Return origin.
An input stream of tokens.
messageStream Info
Information stream (stdout output on master, null elsewhere)
defineDebugSwitch(blockMesh, 0)
bool inplacePointTransforms(pointField &pts) const
Apply coordinate transforms and scaling.
tmp< pointField > globalPosition(const pointField &localPoints) const
Apply coordinate transforms and scaling.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
pointField vertices(const blockVertexList &bvl)
OBJstream os(runTime.globalPath()/outputName)
Input from string buffer, using a ISstream. Always UNCOMPRESSED.
Vector< scalar > vector
A scalar version of the templated Vector.
bool hasPointTransforms() const noexcept
True if scaling and/or transformations are needed.
label numZonedBlocks() const
Number of blocks with specified zones.
const pointField & points() const
unsigned int count(const UList< bool > &bools, const bool val=true)
Count number of 'true' entries.
wordList patchNames() const
Return the patch names.
Output to string buffer, using a OSstream. Always UNCOMPRESSED.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
Class used for the read-construction of.
void checkITstream(const ITstream &is) const
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
const cellShapeList & cells() const
Return cell shapes list.
static int getVerbosity(const dictionary &dict, int verbosity)
const token & peek() const
Failsafe peek at what the next read would return,.
bool equal(const T &s1, const T &s2)
Compare two values for equality.
Tensor< scalar > tensor
Tensor of scalars, i.e. Tensor<scalar>.
static bool verboseOutput
The default verbosity (true)
static constexpr direction nComponents
Number of components in this vector space.
virtual const tensor & R() const
Return const reference to the rotation tensor.