axesRotation.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "axesRotation.H"
27 #include "dictionary.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34  defineTypeNameAndDebug(axesRotation, 0);
35  addToRunTimeSelectionTable(coordinateRotation, axesRotation, dictionary);
37  (
38  coordinateRotation,
39  axesRotation,
40  objectRegistry
41  );
42 }
43 
44 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
45 
46 void Foam::axesRotation::calcTransform
47 (
48  const vector& axis1,
49  const vector& axis2,
50  const axisOrder& order
51 )
52 {
53  vector a = axis1/mag(axis1);
54  vector b = axis2;
55 
56  b = b - (b & a)*a;
57 
58  if (mag(b) < small)
59  {
61  << "axis1, axis2 appear co-linear: "
62  << axis1 << ", " << axis2 << endl
63  << abort(FatalError);
64  }
65 
66  b = b/mag(b);
67  vector c = a^b;
68 
69  tensor Rtr;
70  switch (order)
71  {
72  case e1e2:
73  {
74  Rtr = tensor(a, b, c);
75  break;
76  }
77  case e2e3:
78  {
79  Rtr = tensor(c, a, b);
80  break;
81  }
82  case e3e1:
83  {
84  Rtr = tensor(b, c, a);
85  break;
86  }
87  default:
88  {
90  << "Unhandled axes specification" << endl
91  << abort(FatalError);
92 
93  Rtr = Zero;
94  break;
95  }
96  }
97 
98  // Global->local transformation
99  Rtr_ = Rtr;
100 
101  // Local->global transformation
102  R_ = Rtr.T();
103 }
104 
105 
106 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
107 
109 :
110  R_(sphericalTensor::I),
111  Rtr_(R_)
112 {}
113 
114 
116 (
117  const vector& axis,
118  const vector& dir
119 )
120 :
121  R_(sphericalTensor::I),
122  Rtr_(R_)
123 {
124  calcTransform(axis, dir, e3e1);
125 }
126 
127 
129 (
130  const dictionary& dict
131 )
132 :
133  R_(sphericalTensor::I),
134  Rtr_(R_)
135 {
136  operator=(dict);
137 }
138 
139 
141 (
142  const dictionary& dict,
143  const objectRegistry& obr
144 )
145 :
146  R_(sphericalTensor::I),
147  Rtr_(R_)
148 {
149  operator=(dict);
150 }
151 
152 
154 :
155  R_(R),
156  Rtr_(R_.T())
157 {}
158 
159 
160 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
161 
163 {
165  return NullObjectRef<tensorField>();
166 }
167 
168 
170 (
171  const vectorField& st
172 ) const
173 {
174  return (R_ & st);
175 }
176 
177 
179 {
180  return (R_ & st);
181 }
182 
183 
185 (
186  const vectorField& st
187 ) const
188 {
189  return (Rtr_ & st);
190 }
191 
192 
194 {
195  return (Rtr_ & st);
196 }
197 
198 
200 (
201  const tensorField& st
202 ) const
203 {
205  return tmp<tensorField>(nullptr);
206 }
207 
208 
210 (
211  const tensor& st
212 ) const
213 {
214  return (R_ & st & Rtr_);
215 }
216 
217 
219 (
220  const tensorField& st,
221  const labelList& cellMap
222 ) const
223 {
225  return tmp<tensorField>(nullptr);
226 }
227 
228 
230 (
231  const vectorField& st
232 ) const
233 {
234  tmp<symmTensorField> tfld(new symmTensorField(st.size()));
235  symmTensorField& fld = tfld.ref();
236 
237  forAll(fld, i)
238  {
239  fld[i] = transformPrincipal(R_, st[i]);
240  }
241  return tfld;
242 }
243 
244 
246 (
247  const vector& st
248 ) const
249 {
250  return transformPrincipal(R_, st);
251 }
252 
253 
254 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
255 
257 {
258  if (debug)
259  {
260  Pout<< "axesRotation::operator=(const dictionary&) : "
261  << "assign from " << dict << endl;
262  }
263 
264  vector axis1, axis2;
265  axisOrder order(e3e1);
266 
267  if (dict.readIfPresent("e1", axis1) && dict.readIfPresent("e2", axis2))
268  {
269  order = e1e2;
270  }
271  else if (dict.readIfPresent("e2", axis1)&& dict.readIfPresent("e3", axis2))
272  {
273  order = e2e3;
274  }
275  else if (dict.readIfPresent("e3", axis1)&& dict.readIfPresent("e1", axis2))
276  {
277  order = e3e1;
278  }
279  else if (dict.found("axis") || dict.found("direction"))
280  {
281  // Both "axis" and "direction" are required
282  // If one is missing the appropriate error message will be generated
283  order = e3e1;
284  dict.lookup("axis") >> axis1;
285  dict.lookup("direction") >> axis2;
286  }
287  else
288  {
290  << "not entry of the type (e1, e2) or (e2, e3) or (e3, e1) "
291  << "found "
292  << exit(FatalError);
293  }
294 
295  calcTransform(axis1, axis2, order);
296 }
297 
298 
300 {
301  writeEntry(os, "e1", e1());
302  writeEntry(os, "e2", e2());
303  writeEntry(os, "e3", e3());
304 }
305 
306 
307 // ************************************************************************* //
axesRotation()
Construct null.
Definition: axesRotation.C:108
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:667
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
symmTensor transformPrincipal(const tensor &, const vector &) const
Transform principal.
Field< symmTensor > symmTensorField
Specialisation of Field<T> for symmTensor.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
virtual const vector e2() const
Return local Cartesian y-axis in global coordinates.
Definition: axesRotation.H:161
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
virtual tmp< vectorField > invTransform(const vectorField &st) const
Inverse transform vectorField using transformation tensor field.
Definition: axesRotation.C:185
virtual tmp< symmTensorField > transformVector(const vectorField &st) const
Transform vectorField using transformation tensorField and return.
Definition: axesRotation.C:230
virtual const tensor & R() const
Return local-to-global transformation tensor.
Definition: axesRotation.H:143
Macros for easy insertion into run-time selection tables.
Templated 3D SphericalTensor derived from VectorSpace adding construction from 1 component, element access using th ii() member function and the inner-product (dot-product) and outer-product operators.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< ' ';}gmvFile<< nl;forAll(lagrangianScalarNames, i){ const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
static const Identity< scalar > I
Definition: Identity.H:93
const dimensionedScalar & b
Wien displacement law constant: default SI units: [m K].
Definition: createFields.H:27
virtual const vector e3() const
Return local Cartesian z-axis in global coordinates.
Definition: axesRotation.H:167
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
virtual const tensorField & Tr() const
Return transformation tensor field.
Definition: axesRotation.C:162
static const zero Zero
Definition: zero.H:97
errorManip< error > abort(error &err)
Definition: errorManip.H:131
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
virtual tmp< vectorField > transform(const vectorField &st) const
Transform vectorField using transformation tensor field.
Definition: axesRotation.C:170
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
defineTypeNameAndDebug(combustionModel, 0)
virtual const vector e1() const
Return local Cartesian x-axis in global coordinates.
Definition: axesRotation.H:155
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
void operator=(const dictionary &)
Assign from dictionary.
Definition: axesRotation.C:256
virtual void write(Ostream &) const
Write.
Definition: axesRotation.C:299
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
dimensioned< scalar > mag(const dimensioned< Type > &)
A class for managing temporary objects.
Definition: PtrList.H:53
Registry of regIOobjects.
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:366
virtual tmp< tensorField > transformTensor(const tensorField &st) const
Transform tensor field using transformation tensorField.
Definition: axesRotation.C:200
Tensor< scalar > tensor
Tensor of scalars.
Definition: tensor.H:51
Namespace for OpenFOAM.
static const SphericalTensor I
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:812