lumpedPointState.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2016-2022 OpenCFD Ltd.
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
28 #include "lumpedPointState.H"
29 #include "unitConversion.H"
31 #include "IFstream.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 const Foam::Enum
36 <
38 >
40 ({
41  { inputFormatType::PLAIN, "plain" },
42  { inputFormatType::DICTIONARY, "dictionary" },
43 });
44 
45 
46 Foam::scalar Foam::lumpedPointState::visLength = 0.1;
47 
49 
50 
51 // * * * * * * * * * * * * * * * Local Functions * * * * * * * * * * * * * * //
52 
54 static Foam::string getLineNoComment
55 (
56  Foam::ISstream& is,
57  const char comment = '#'
58 )
59 {
60  Foam::string line;
61  do
62  {
63  is.getLine(line);
64  }
65  while ((line.empty() || line[0] == comment) && is.good());
66 
67  return line;
68 }
70 
71 
72 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
73 
74 void Foam::lumpedPointState::calcRotations() const
75 {
76  rotationPtr_.reset(new tensorField(angles_.size()));
77 
78  auto rotIter = rotationPtr_->begin();
79 
80  for (const vector& angles : angles_)
81  {
82  *rotIter =
84 
85  ++rotIter;
86  }
87 }
88 
89 
90 void Foam::lumpedPointState::readDict
91 (
92  const dictionary& dict,
93  const quaternion::eulerOrder rotOrder,
94  const bool degrees
95 )
96 {
97  dict.readEntry("points", points_);
98  dict.readEntry("angles", angles_);
99  order_ =
101  (
102  "rotationOrder",
103  dict,
104  rotOrder
105  );
106 
107  degrees_ = dict.getOrDefault("degrees", degrees);
109  rotationPtr_.reset(nullptr);
110 }
111 
112 
113 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
114 
116 :
117  points_(),
118  angles_(),
119  order_(quaternion::eulerOrder::ZXZ),
120  degrees_(false),
121  rotationPtr_(nullptr)
122 {}
123 
124 
126 :
127  points_(rhs.points_),
128  angles_(rhs.angles_),
129  order_(rhs.order_),
130  degrees_(rhs.degrees_),
131  rotationPtr_(nullptr)
132 {}
133 
134 
136 (
137  const pointField& pts,
138  const vectorField& ang,
139  const quaternion::eulerOrder rotOrder,
140  const bool degrees
141 )
142 :
143  points_(pts),
144  angles_(ang),
145  order_(rotOrder),
146  degrees_(degrees),
147  rotationPtr_(nullptr)
148 {
149  if (points_.size() != angles_.size())
150  {
151  #ifdef FULLDEBUG
153  << "Have " << points_.size() << " points but "
154  << angles_.size() << " angles" << nl
155  << exit(FatalError);
156  #else
158  << "Have " << points_.size() << " points but "
159  << angles_.size() << " angles, resizing angles to match" << nl;
160  #endif
162  angles_.resize(points_.size(), Zero);
163  }
164 }
165 
166 
168 (
169  const pointField& pts,
170  const quaternion::eulerOrder rotOrder,
171  const bool degrees
172 )
173 :
174  points_(pts),
175  angles_(points_.size(), Zero),
176  order_(rotOrder),
177  degrees_(degrees),
178  rotationPtr_(nullptr)
179 {}
180 
181 
183 (
185  const quaternion::eulerOrder rotOrder,
186  const bool degrees
187 )
188 :
189  points_(pts),
190  angles_(points_.size(), Zero),
191  order_(rotOrder),
192  degrees_(false),
193  rotationPtr_(nullptr)
194 {}
195 
196 
198 (
199  const dictionary& dict,
200  const quaternion::eulerOrder rotOrder,
201  const bool degrees
202 )
203 :
204  points_(),
205  angles_(),
206  order_(rotOrder),
207  degrees_(degrees),
208  rotationPtr_(nullptr)
209 {
210  readDict(dict, rotOrder, degrees);
211 }
212 
213 
214 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
215 
217 {
218  points_ = rhs.points_;
219  angles_ = rhs.angles_;
220  order_ = rhs.order_;
221  degrees_ = rhs.degrees_;
222 
223  rotationPtr_.reset(nullptr);
224 }
225 
226 
227 void Foam::lumpedPointState::operator+=(const point& origin)
228 {
229  for (point& p : points_)
230  {
231  p += origin;
232  }
233 }
234 
235 
236 void Foam::lumpedPointState::scalePoints(const scalar scaleFactor)
237 {
238  if (scaleFactor > 0)
239  {
240  points_ *= scaleFactor;
241  }
242 }
243 
244 
246 (
247  const scalar alpha,
248  const lumpedPointState& prev
249 )
250 {
251  points_ = prev.points_ + alpha*(points_ - prev.points_);
252 
253  scalar convert = 1.0;
254  if (degrees_ != prev.degrees_)
255  {
256  if (prev.degrees_)
257  {
258  // Was degrees, now radians
259  convert = degToRad();
260  }
261  else
262  {
263  // Was radians, now degrees
264  convert = radToDeg();
265  }
266  }
267 
268  angles_ = convert*prev.angles_ + alpha*(angles_ - convert*prev.angles_);
269 
270  rotationPtr_.reset(nullptr);
271 }
272 
273 
275 (
276  Istream& is,
277  const quaternion::eulerOrder rotOrder,
278  const bool degrees
279 )
280 {
281  // Assume generic input stream so we can do line-based input
282  ISstream& iss = dynamic_cast<ISstream&>(is);
283 
284  string line = getLineNoComment(iss);
285 
286  label count;
287  {
288  IStringStream isstr(line);
289  isstr >> count;
290  }
291 
292  points_.resize(count);
293  angles_.resize(count);
294 
295  count = 0;
296  forAll(points_, i)
297  {
298  iss.getLine(line);
299  IStringStream isstr(line);
300 
301  isstr
302  >> points_[count].x() >> points_[count].y() >> points_[count].z()
303  >> angles_[count].x() >> angles_[count].y() >> angles_[count].z();
304 
305  ++count;
306  }
307 
308  points_.resize(count);
309  angles_.resize(count);
310  order_ = quaternion::eulerOrder::ZXZ;
311  degrees_ = false;
312 
313  rotationPtr_.reset(nullptr);
314 
315  return count;
316 }
317 
318 
320 (
321  Istream& is,
322  const quaternion::eulerOrder rotOrder,
323  const bool degrees
324 )
325 {
327  readDict(dict, rotOrder, degrees);
328 
329  return points_.size();
330 }
331 
332 
334 {
335  writeDict(os);
336  return true;
337 }
338 
339 
341 {
342  os.writeEntry("points", points_);
343  os.writeEntry("angles", angles_);
344  if (order_ != quaternion::eulerOrder::ZXZ)
345  {
346  os.writeEntry("order", quaternion::eulerOrderNames[order_]);
347  }
348  if (degrees_)
349  {
350  os.writeEntry("degrees", "true");
351  }
352 }
353 
354 
355 void Foam::lumpedPointState::writePlain(Ostream& os) const
356 {
357  os <<"# input for OpenFOAM\n"
358  <<"# N, points, angles\n"
359  << points_.size() << "\n";
360 
361  forAll(points_, i)
362  {
363  const vector& p = points_[i];
364 
365  os << p.x() << ' ' << p.y() << ' ' << p.z();
366 
367  if (i < angles_.size())
368  {
369  os << ' ' << angles_[i].x()
370  << ' ' << angles_[i].y()
371  << ' ' << angles_[i].z() << '\n';
372  }
373  else
374  {
375  os << " 0 0 0\n";
376  }
377  }
378 }
379 
380 
382 (
383  const inputFormatType fmt,
384  const fileName& file,
385  const quaternion::eulerOrder rotOrder,
386  const bool degrees
387 )
388 {
389  bool ok = false;
390  if (Pstream::master())
391  {
392  IFstream is(file);
393 
394  if (fmt == inputFormatType::PLAIN)
395  {
396  ok = this->readPlain(is, rotOrder, degrees);
397  }
398  else
399  {
400  ok = this->readData(is, rotOrder, degrees);
401  }
402  }
403 
404  if (Pstream::parRun())
405  {
406  // Broadcast master data to everyone
407 
409  (
411  ok,
412  degrees_,
413  points_,
414  angles_
415  );
416  }
417  rotationPtr_.reset(nullptr);
418 
419  return ok;
420 }
421 
422 
423 // ************************************************************************* //
dictionary dict
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:116
A line primitive.
Definition: line.H:52
static scalar visLength
The length for visualization triangles.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:125
void resize(const label len)
Adjust allocated size of list.
Definition: ListI.H:160
inputFormatType
Input format types.
error FatalError
Error stream (stdout output on all processes), with additional &#39;FOAM FATAL ERROR&#39; header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
Definition: dictionary.H:129
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:598
void writeDict(Ostream &os) const
Output as dictionary content.
void writePlain(Ostream &os) const
Output as plain content.
Unit conversion functions.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
constexpr char nl
The newline &#39;\n&#39; character (0x0a)
Definition: Ostream.H:50
static bool & parRun() noexcept
Test if this a parallel run.
Definition: UPstream.H:1049
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:321
void operator=(const lumpedPointState &rhs)
Assignment operator.
static label worldComm
Communicator for all ranks. May differ from commGlobal() if local worlds are in use.
Definition: UPstream.H:409
static const Enum< eulerOrder > eulerOrderNames
The names for Euler-angle and Tait-Bryan angles, including "rollPitchYaw" and "yawPitchRoll" aliases...
Definition: quaternion.H:132
EnumType getOrDefault(const word &key, const dictionary &dict, const EnumType deflt, const bool warnOnly=false) const
Find the key in the dictionary and return the corresponding enumeration element based on its name...
Definition: Enum.C:173
bool writeData(Ostream &os) const
Output as dictionary content.
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:421
static bool visUnused
Enable/disable visualization of unused points.
bool readPlain(Istream &is, const quaternion::eulerOrder rotOrder=quaternion::eulerOrder::ZXZ, const bool degrees=false)
Read input as plain content.
static const Enum< inputFormatType > formatNames
Names for the input format types.
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...
unsigned int count(const UList< bool > &bools, const bool val=true)
Count number of &#39;true&#39; entries.
Definition: BitOps.H:73
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:38
void operator+=(const point &origin)
Shift points by specified origin.
const vectorField & angles() const
The orientation of the points (mass centres)
Quaternion class used to perform rotations in 3D space.
Definition: quaternion.H:53
Vector< scalar > vector
Definition: vector.H:57
eulerOrder
Euler-angle rotation order.
Definition: quaternion.H:115
lumpedPointState()
Default construct.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:56
ISstream & getLine(std::string &str, char delim='\n')
Raw, low-level getline (until delimiter) into a string.
Definition: ISstreamI.H:69
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
OBJstream os(runTime.globalPath()/outputName)
void scalePoints(const scalar scaleFactor)
Scale points by given factor.
static void broadcasts(const label comm, Type &arg1, Args &&... args)
Broadcast multiple items to all communicator ranks. Does nothing in non-parallel. ...
Field< tensor > tensorField
Specialisation of Field<T> for tensor.
Generic input stream using a standard (STL) stream.
Definition: ISstream.H:51
#define WarningInFunction
Report a warning using Foam::Warning.
static tensor rotation(const vector &angles, bool degrees=false)
Rotation tensor calculated for the intrinsic Euler angles in z-x-z order.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Definition: error.H:64
constexpr scalar radToDeg(const scalar rad) noexcept
Conversion from radians to degrees.
void relax(const scalar alpha, const lumpedPointState &prev)
Relax the state.
Input from string buffer, using a ISstream. Always UNCOMPRESSED.
Definition: StringStream.H:120
bool good() const noexcept
True if next operation might succeed.
Definition: IOstream.H:281
static bool master(const label communicator=worldComm)
True if process corresponds to the master rank in the communicator.
Definition: UPstream.H:1082
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...
volScalarField & p
A class for managing temporary objects.
Definition: HashPtrTable.H:50
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
bool readData(Istream &is, const quaternion::eulerOrder rotOrder=quaternion::eulerOrder::ZXZ, const bool degrees=false)
Read input as dictionary content.
constexpr scalar degToRad(const scalar deg) noexcept
Conversion from degrees to radians.
A class for handling character strings derived from std::string.
Definition: string.H:72
The state of lumped points corresponds to positions and rotations.
const pointField & pts
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:127