Go to the documentation of this file.
41 { inputFormatType::PLAIN,
"plain" },
42 { inputFormatType::DICTIONARY,
"dictionary" },
57 const char comment =
'#'
65 while ((line.empty() || line[0] == comment) && is.
good());
74 void Foam::lumpedPointState::calcRotations()
const
76 rotationPtr_.reset(
new tensorField(angles_.size()));
78 auto rotIter = rotationPtr_->begin();
90 void Foam::lumpedPointState::readDict
92 const dictionary&
dict,
109 rotationPtr_.reset(
nullptr);
121 rotationPtr_(nullptr)
127 points_(rhs.points_),
128 angles_(rhs.angles_),
130 degrees_(rhs.degrees_),
131 rotationPtr_(nullptr)
147 rotationPtr_(
nullptr)
149 if (points_.size() != angles_.size())
153 <<
"Have " << points_.size() <<
" points but "
154 << angles_.size() <<
" angles" <<
nl
158 <<
"Have " << points_.size() <<
" points but "
159 << angles_.size() <<
" angles, resizing angles to match" <<
nl;
162 angles_.resize(points_.size(),
Zero);
175 angles_(points_.size(),
Zero),
178 rotationPtr_(
nullptr)
190 angles_(points_.size(),
Zero),
193 rotationPtr_(
nullptr)
208 rotationPtr_(
nullptr)
210 readDict(
dict, rotOrder, degrees);
218 points_ = rhs.points_;
219 angles_ = rhs.angles_;
221 degrees_ = rhs.degrees_;
223 rotationPtr_.reset(
nullptr);
240 points_ *= scaleFactor;
251 points_ = prev.points_ +
alpha*(points_ - prev.points_);
253 scalar convert = 1.0;
254 if (degrees_ != prev.degrees_)
268 angles_ = convert*prev.angles_ +
alpha*(angles_ - convert*prev.angles_);
270 rotationPtr_.reset(
nullptr);
284 string line = getLineNoComment(iss);
292 points_.resize(
count);
293 angles_.resize(
count);
308 points_.resize(
count);
309 angles_.resize(
count);
310 order_ = quaternion::eulerOrder::ZXZ;
313 rotationPtr_.reset(
nullptr);
327 readDict(
dict, rotOrder, degrees);
329 return points_.size();
344 if (order_ != quaternion::eulerOrder::ZXZ)
357 os <<
"# input for OpenFOAM\n"
358 <<
"# N, points, angles\n"
359 << points_.size() <<
"\n";
365 os <<
p.x() <<
' ' <<
p.y() <<
' ' <<
p.z();
367 if (i < angles_.size())
369 os <<
' ' << angles_[i].x()
370 <<
' ' << angles_[i].y()
371 <<
' ' << angles_[i].z() <<
'\n';
394 if (fmt == inputFormatType::PLAIN)
396 ok = this->readPlain(is, rotOrder, degrees);
400 ok = this->readData(is, rotOrder, degrees);
419 if (myComm.
above() != -1)
430 fromAbove >> points_ >> angles_ >> degrees_;
439 myComm.
below()[belowI],
445 toBelow << points_ << angles_ << degrees_;
448 rotationPtr_.reset(
nullptr);
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
void operator+=(const point &origin)
Shift points by specified origin.
void writeDict(Ostream &os) const
Output as dictionary content.
ISstream & getLine(std::string &str, char delim='\n')
Raw, low-level getline (until delimiter) into a string.
A class for handling file names.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
Input from file stream, using an ISstream.
inputFormatType
Input format types.
Output inter-processor communications stream.
void operator=(const lumpedPointState &rhs)
Assignment operator.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
static label nProcs(const label communicator=0)
Number of processes in parallel run.
static bool & parRun()
Is this a parallel run?
void relax(const scalar alpha, const lumpedPointState &prev)
Relax the state.
Field< tensor > tensorField
Specialisation of Field<T> for tensor.
static bool visUnused
Enable/disable visualization of unused points.
Generic input stream using a standard (STL) stream.
Unit conversion functions.
static void scatter(const List< commsStruct > &comms, T &Value, const int tag, const label comm)
Scatter data. Distribute without modification. Reverse of gather.
bool readPlain(Istream &is, const quaternion::eulerOrder rotOrder=quaternion::eulerOrder::ZXZ, const bool degrees=false)
Read input as plain content.
A class for handling character strings derived from std::string.
#define forAll(list, i)
Loop across all elements in list.
Quaternion class used to perform rotations in 3D space.
bool readData(Istream &is, const quaternion::eulerOrder rotOrder=quaternion::eulerOrder::ZXZ, const bool degrees=false)
Read input as dictionary content.
The state of lumped points corresponds to positions and rotations.
void scalePoints(const scalar scaleFactor)
Scale points by given factor.
static const Enum< inputFormatType > formatNames
Names for the input format types.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
const vectorField & angles() const
The orientation of the points (mass centres)
EnumType getOrDefault(const word &key, const dictionary &dict, const EnumType deflt, const bool failsafe=false) const
void writePlain(Ostream &os) const
Output as plain content.
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
constexpr scalar radToDeg(const scalar rad) noexcept
Conversion from radians to degrees.
eulerOrder
Euler-angle rotation order.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
bool writeData(Ostream &os) const
Output as dictionary content.
Input from string buffer, using a ISstream.
Vector< scalar > vector
A scalar version of the templated Vector.
constexpr scalar degToRad(const scalar deg) noexcept
Conversion from degrees to radians.
Structure for communicating between processors.
errorManipArg< error, int > exit(error &err, const int errNo=1)
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
static bool master(const label communicator=0)
Am I the master process.
static int nProcsSimpleSum
static const Enum< eulerOrder > eulerOrderNames
The names for Euler-angle rotation order.
static int & msgType()
Message tag of standard messages.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
unsigned int count(const UList< bool > &bools, const bool val=true)
Count number of 'true' entries.
static label worldComm
Default communicator (all processors)
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
#define forAllReverse(list, i)
Reverse loop across all elements in list.
Input inter-processor communications stream.
static tensor rotation(const vector &angles, bool degrees=false)
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
lumpedPointState()
Default construct.
static const List< commsStruct > & linearCommunication(const label communicator=0)
Communication schedule for linear all-to-master (proc 0)
bool good() const
Return true if next operation might succeed.
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
#define WarningInFunction
Report a warning using Foam::Warning.
static scalar visLength
The length for visualization triangles.
const labelList & below() const
static const List< commsStruct > & treeCommunication(const label communicator=0)
Communication schedule for tree all-to-master (proc 0)