globalIndexAndTransformI.H
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) 2011-2016 OpenFOAM Foundation
9  Copyright (C) 2018 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "polyMesh.H"
30 
31 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
32 
33 bool Foam::globalIndexAndTransform::less::operator()
34 (
35  const labelPair& a,
36  const labelPair& b
37 ) const
38 {
39  label procA = gi_.processor(a);
40  label procB = gi_.processor(b);
41 
42  if (procA < procB)
43  {
44  return true;
45  }
46  else if (procA > procB)
47  {
48  return false;
49  }
50  else
51  {
52  // Equal proc.
53  label indexA = gi_.index(a);
54  label indexB = gi_.index(b);
55 
56  if (indexA < indexB)
57  {
58  return true;
59  }
60  else if (indexA > indexB)
61  {
62  return false;
63  }
64  else
65  {
66  // Equal index
67  label transformA = gi_.transformIndex(a);
68  label transformB = gi_.transformIndex(b);
69 
70  return transformA < transformB;
71  }
72  }
73 }
74 
75 
77 (
78  const labelList& permutationIndices
79 ) const
80 {
81  if (permutationIndices.size() != transforms_.size())
82  {
84  << "permutationIndices " << permutationIndices
85  << "are of a different size to the number of independent transforms"
86  << abort(FatalError);
87  }
88 
89  label transformIndex = 0;
90 
91  label w = 1;
92 
93  forAll(transforms_, b)
94  {
95  if (mag(permutationIndices[b]) > 1)
96  {
98  << "permutationIndices " << permutationIndices
99  << "are illegal, they must all be only -1, 0 or +1"
100  << abort(FatalError);
101  }
102 
103  transformIndex += (permutationIndices[b] + 1)*w;
104 
105  w *= 3;
106  }
107 
108  return transformIndex;
109 }
110 
111 
113 (
114  const label transformIndex
115 ) const
116 {
117  labelList permutation(transforms_.size(), Zero);
118 
119  label t = transformIndex;
120  forAll(permutation, i)
121  {
122  permutation[i] = (t%3)-1;
123  t /= 3;
124  }
125 
126  return permutation;
127 }
128 
129 
131 (
132  const label transformIndex,
133  const label patchi,
134  const bool isSendingSide,
135  const scalar tol
136 ) const
137 {
138  const labelPair& transSign = patchTransformSign_[patchi];
139 
140  label matchTransI = transSign.first();
141 
142  if (matchTransI >= transforms_.size())
143  {
145  << "patch:" << mesh_.boundaryMesh()[patchi].name()
146  << " transform:" << matchTransI
147  << " out of possible transforms:" << transforms_
148  << exit(FatalError);
149  return labelMin;
150  }
151  else if (matchTransI == -1)
152  {
153  // No additional transformation for this patch
154  return transformIndex;
155  }
156  else
157  {
158  // Decode current set of transforms
159  labelList permutation(decodeTransformIndex(transformIndex));
160 
161 
162  // Add patch transform
163  // ~~~~~~~~~~~~~~~~~~~
164 
165  label sign = transSign.second();
166  if (!isSendingSide)
167  {
168  sign = -sign;
169  }
170 
171 
172  // If this transform been found already by a patch?
173  if (permutation[matchTransI] != 0)
174  {
175  if (sign == 0)
176  {
177  // sent from patch without a transformation. Do nothing.
179  << "patch:" << mesh_.boundaryMesh()[patchi].name()
180  << " transform:" << matchTransI << " sign:" << sign
181  << " current transforms:" << permutation
182  << exit(FatalError);
183  }
184  else if (sign == permutation[matchTransI])
185  {
186  // This is usually illegal. The only exception is for points
187  // on the axis of a 180 degree cyclic wedge when the
188  // transformation is going to be (-1 0 0 0 -1 0 0 0 +1)
189  // (or a different permutation but always two times -1 and
190  // once +1)
191  bool antiCyclic = false;
192 
193  const vectorTensorTransform& vt = transforms_[matchTransI];
194  if (mag(vt.t()) < SMALL && vt.hasR())
195  {
196  const tensor& R = vt.R();
197  scalar sumDiag = tr(R);
198  scalar sumMagDiag = mag(R.xx())+mag(R.yy())+mag(R.zz());
199 
200  if (mag(sumMagDiag-3) < tol && mag(sumDiag+1) < tol)
201  {
202  antiCyclic = true;
203  }
204  }
205 
206  if (antiCyclic)
207  {
208  // 180 degree rotational. Reset transformation.
209  permutation[matchTransI] = 0;
210  }
211  else
212  {
214  << "More than one patch accessing the same transform "
215  << "but not of the same sign." << endl
216  << "patch:" << mesh_.boundaryMesh()[patchi].name()
217  << " transform:" << matchTransI << " sign:" << sign
218  << " current transforms:" << permutation
219  << exit(FatalError);
220  }
221  }
222  else
223  {
224  permutation[matchTransI] = 0;
225  }
226  }
227  else
228  {
229  permutation[matchTransI] = sign;
230  }
231 
232 
233  // Re-encode permutation
234  // ~~~~~~~~~~~~~~~~~~~~~
235 
236  return encodeTransformIndex(permutation);
237  }
238 }
239 
240 
242 (
243  const label transformIndex0,
244  const label transformIndex1
245 ) const
246 {
247  if (transformIndex0 == transformIndex1)
248  {
249  return transformIndex0;
250  }
251 
252 
253  // Count number of transforms
254  labelList permutation0(decodeTransformIndex(transformIndex0));
255  label n0 = 0;
256  forAll(permutation0, i)
257  {
258  if (permutation0[i] != 0)
259  {
260  n0++;
261  }
262  }
263 
264  labelList permutation1(decodeTransformIndex(transformIndex1));
265  label n1 = 0;
266  forAll(permutation1, i)
267  {
268  if (permutation1[i] != 0)
269  {
270  n1++;
271  }
272  }
273 
274  if (n0 <= n1)
275  {
276  return transformIndex0;
277  }
278  else
279  {
280  return transformIndex1;
281  }
282 }
283 
284 
286 (
287  const label transformIndex0,
288  const label transformIndex1
289 ) const
290 {
291  labelList permutation0(decodeTransformIndex(transformIndex0));
292  labelList permutation1(decodeTransformIndex(transformIndex1));
293 
294  forAll(permutation0, i)
295  {
296  permutation0[i] -= permutation1[i];
297  }
298 
299  return encodeTransformIndex(permutation0);
300 }
301 
302 
304 (
305  const label index,
306  const label transformIndex
307 ) const
308 {
309  return encode(Pstream::myProcNo(), index, transformIndex);
310 }
311 
312 
314 (
315  const label proci,
316  const label index,
317  const label transformIndex
318 ) const
319 {
320  if (transformIndex < 0 || transformIndex >= transformPermutations_.size())
321  {
323  << "TransformIndex " << transformIndex
324  << " is outside allowed range of 0 to "
325  << transformPermutations_.size() - 1
326  << abort(FatalError);
327  }
328 
329  if (proci > labelMax/transformPermutations_.size())
330  {
332  << "Overflow : encoding processor " << proci
333  << " in base " << transformPermutations_.size()
334  << " exceeds capability of label (" << labelMax
335  << "). Please recompile with larger datatype for label."
336  << exit(FatalError);
337  }
338 
339  return labelPair
340  (
341  index,
342  transformIndex + proci*transformPermutations_.size()
343  );
344 }
345 
346 
348 (
349  const labelPair& globalIAndTransform
350 ) const
351 {
352  return globalIAndTransform.first();
353 }
354 
355 
357 (
358  const labelPair& globalIAndTransform
359 ) const
360 {
361  return globalIAndTransform.second()/transformPermutations_.size();
362 }
363 
364 
366 (
367  const labelPair& globalIAndTransform
368 ) const
369 {
370  return globalIAndTransform.second()%transformPermutations_.size();
371 }
372 
373 
375 {
376  return transforms_.size();
377 }
378 
379 
382 {
383  return transforms_;
384 }
385 
386 
389 {
390  return transformPermutations_;
391 }
392 
393 
395 {
396  return nullTransformIndex_;
397 }
398 
399 
400 const Foam::labelPairList&
402 {
403  return patchTransformSign_;
404 }
405 
406 
408 (
409  label transformIndex
410 ) const
411 {
412  return transformPermutations_[transformIndex];
413 }
414 
415 
417 (
418  const labelHashSet& patchis
419 ) const
420 {
421  labelList permutation(transforms_.size(), Zero);
422 
423  labelList selectedTransformIs(0);
424 
425  if (patchis.empty() || transforms_.empty())
426  {
427  return selectedTransformIs;
428  }
429 
430  for (const label patchi : patchis)
431  {
432  const labelPair& transSign = patchTransformSign_[patchi];
433 
434  label matchTransI = transSign.first();
435 
436  if (matchTransI > -1)
437  {
438  label sign = transSign.second();
439 
440  // If this transform been found already by a patch?
441  if (permutation[matchTransI] != 0)
442  {
443  // If so, if they have opposite signs, then this is
444  // considered an error. They are allowed to be the
445  // same sign, but this only results in a single
446  // transform.
447  if (permutation[matchTransI] != sign)
448  {
450  << "More than one patch accessing the same transform "
451  << "but not of the same sign."
452  << exit(FatalError);
453  }
454  }
455  else
456  {
457  permutation[matchTransI] = sign;
458  }
459  }
460  }
461 
462  label nUsedTrans = round(sum(mag(permutation)));
463 
464  if (nUsedTrans == 0)
465  {
466  return selectedTransformIs;
467  }
468 
469  // Number of selected transformations
470  label nSelTrans = pow(label(2), nUsedTrans) - 1;
471 
472  // Pout<< nl << permutation << nl << endl;
473 
474  selectedTransformIs.setSize(nSelTrans);
475 
476  switch (nUsedTrans)
477  {
478  case 1:
479  {
480  selectedTransformIs[0] = encodeTransformIndex(permutation);
481 
482  break;
483  }
484  case 2:
485  {
486  labelList tempPermutation = permutation;
487 
488  label a = 0;
489  label b = 1;
490 
491  // When there are two selected transforms out of three, we
492  // need to choose which of them are being permuted
493  if (transforms_.size() > nUsedTrans)
494  {
495  if (permutation[0] == 0)
496  {
497  a = 1;
498  b = 2;
499  }
500  else if (permutation[1] == 0)
501  {
502  a = 0;
503  b = 2;
504  }
505  else if (permutation[2] == 0)
506  {
507  a = 0;
508  b = 1;
509  }
510  }
511 
512  tempPermutation[a] = a;
513  tempPermutation[b] = permutation[b];
514 
515  selectedTransformIs[0] = encodeTransformIndex(tempPermutation);
516 
517  tempPermutation[a] = permutation[a];
518  tempPermutation[b] = a;
519 
520  selectedTransformIs[1] = encodeTransformIndex(tempPermutation);
521 
522  tempPermutation[a] = permutation[a];
523  tempPermutation[b] = permutation[b];
524 
525  selectedTransformIs[2] = encodeTransformIndex(tempPermutation);
526 
527  break;
528  }
529  case 3:
530  {
531  labelList tempPermutation = permutation;
532 
533  tempPermutation[0] = 0;
534  tempPermutation[1] = 0;
535  tempPermutation[2] = permutation[2];
536 
537  selectedTransformIs[0] = encodeTransformIndex(tempPermutation);
538 
539  tempPermutation[0] = 0;
540  tempPermutation[1] = permutation[1];
541  tempPermutation[2] = 0;
542 
543  selectedTransformIs[1] = encodeTransformIndex(tempPermutation);
544 
545  tempPermutation[0] = 0;
546  tempPermutation[1] = permutation[1];
547  tempPermutation[2] = permutation[2];
548 
549  selectedTransformIs[2] = encodeTransformIndex(tempPermutation);
550 
551  tempPermutation[0] = permutation[0];
552  tempPermutation[1] = 0;
553  tempPermutation[2] = 0;
554 
555  selectedTransformIs[3] = encodeTransformIndex(tempPermutation);
556 
557  tempPermutation[0] = permutation[0];
558  tempPermutation[1] = 0;
559  tempPermutation[2] = permutation[2];
560 
561  selectedTransformIs[4] = encodeTransformIndex(tempPermutation);
562 
563  tempPermutation[0] = permutation[0];
564  tempPermutation[1] = permutation[1];
565  tempPermutation[2] = 0;
566 
567  selectedTransformIs[5] = encodeTransformIndex(tempPermutation);
568 
569  tempPermutation[0] = permutation[0];
570  tempPermutation[1] = permutation[1];
571  tempPermutation[2] = permutation[2];
572 
573  selectedTransformIs[6] = encodeTransformIndex(tempPermutation);
574 
575  break;
576  }
577  default:
578  {
580  << "Only 1-3 transforms are possible."
581  << exit(FatalError);
582  }
583  }
584 
585  return selectedTransformIs;
586 }
587 
588 
590 (
591  const labelHashSet& patchis,
592  const point& pt
593 ) const
594 {
595  labelList transIs = transformIndicesForPatches(patchis);
596 
597  // Pout<< patchis << nl << transIs << endl;
598 
599  pointField transPts(transIs.size());
600 
601  forAll(transIs, tII)
602  {
603  transPts[tII] = transformPermutations_[transIs[tII]].transformPosition
604  (
605  pt
606  );
607  }
608 
609  return transPts;
610 }
611 
612 
613 // ************************************************************************* //
Foam::globalIndexAndTransform::processor
label processor(const labelPair &globalIAndTransform) const
Which processor does this come from?
Definition: globalIndexAndTransformI.H:357
Foam::Pair::second
const T & second() const noexcept
Return second element, which is also the last element.
Definition: PairI.H:122
Foam::Tensor< scalar >
Foam::globalIndexAndTransform::decodeTransformIndex
labelList decodeTransformIndex(const label transformIndex) const
Decode transform index.
Definition: globalIndexAndTransformI.H:113
Foam::globalIndexAndTransform::patchTransformSign
const List< labelPair > & patchTransformSign() const
Return access to the per-patch transform-sign pairs.
Definition: globalIndexAndTransformI.H:401
Foam::globalIndexAndTransform::addToTransformIndex
label addToTransformIndex(const label transformIndex, const label patchi, const bool isSendingSide=true, const scalar tol=SMALL) const
Add patch transformation to transformIndex. Return new.
Definition: globalIndexAndTransformI.H:131
Foam::labelMax
constexpr label labelMax
Definition: label.H:61
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::globalIndexAndTransform::index
label index(const labelPair &globalIAndTransform) const
Index carried by the object.
Definition: globalIndexAndTransformI.H:348
Foam::globalIndexAndTransform::transformPermutations
const List< vectorTensorTransform > & transformPermutations() const
Return access to the permuted transforms.
Definition: globalIndexAndTransformI.H:388
Foam::globalIndexAndTransform::nIndependentTransforms
label nIndependentTransforms() const
Return the number of independent transforms.
Definition: globalIndexAndTransformI.H:374
Foam::globalIndexAndTransform::encode
labelPair encode(const label index, const label transformIndex) const
Encode index and bare index as components on own processor.
Definition: globalIndexAndTransformI.H:304
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::sign
dimensionedScalar sign(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:166
polyMesh.H
Foam::HashSet< label, Hash< label > >
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::labelPair
Pair< label > labelPair
A pair of labels.
Definition: Pair.H:54
R
#define R(A, B, C, D, E, F, K, M)
Foam::globalIndexAndTransform::nullTransformIndex
label nullTransformIndex() const
Return the transformIndex (index in transformPermutations)
Definition: globalIndexAndTransformI.H:394
Foam::constant::physicoChemical::b
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:27
Foam::Field< vector >
Foam::globalIndexAndTransform::transformIndex
label transformIndex(const labelPair &globalIAndTransform) const
Transform carried by the object.
Definition: globalIndexAndTransformI.H:366
Foam::globalIndexAndTransform::encodeTransformIndex
label encodeTransformIndex(const labelList &permutationIndices) const
Generate a transform index from the permutation indices of.
Definition: globalIndexAndTransformI.H:77
Foam::IOstream::name
virtual const fileName & name() const
Return the name of the stream.
Definition: IOstream.C:40
Foam::pow
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Definition: dimensionedScalar.C:75
Foam::FatalError
error FatalError
Foam::globalIndexAndTransform::transformIndicesForPatches
labelList transformIndicesForPatches(const labelHashSet &patchIs) const
Access the all of the indices of the transform.
Definition: globalIndexAndTransformI.H:417
Foam::vectorTensorTransform::t
const vector & t() const
Definition: vectorTensorTransformI.H:69
Foam::globalIndexAndTransform::transformPatches
pointField transformPatches(const labelHashSet &patchIs, const point &pt) const
Apply all of the transform permutations.
Definition: globalIndexAndTransformI.H:590
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
Foam::vectorTensorTransform::hasR
bool hasR() const
Definition: vectorTensorTransformI.H:81
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::globalIndexAndTransform::subtractTransformIndex
label subtractTransformIndex(const label transformIndex0, const label transformIndex1) const
Subtract two transformIndices.
Definition: globalIndexAndTransformI.H:286
Foam::globalIndexAndTransform::transform
const vectorTensorTransform & transform(label transformIndex) const
Access the overall (permuted) transform corresponding.
Definition: globalIndexAndTransformI.H:408
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::Pair< label >
Foam::Vector< scalar >
Foam::List< label >
Foam::vectorTensorTransform
Vector-tensor class used to perform translations and rotations in 3D space.
Definition: vectorTensorTransform.H:63
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::vectorTensorTransform::R
const tensor & R() const
Definition: vectorTensorTransformI.H:75
Foam::sum
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
Definition: DimensionedFieldFunctions.C:327
Foam::labelMin
constexpr label labelMin
Definition: label.H:60
Foam::tr
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
Definition: dimensionedSphericalTensor.C:51
Foam::globalIndexAndTransform::transforms
const List< vectorTensorTransform > & transforms() const
Return access to the stored independent transforms.
Definition: globalIndexAndTransformI.H:381
Foam::globalIndexAndTransform::minimumTransformIndex
label minimumTransformIndex(const label transformIndex0, const label transformIndex1) const
Combine two transformIndices.
Definition: globalIndexAndTransformI.H:242