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 | Copyright (C) 2011-2013 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 "polyMesh.H"
27 
28 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
29 
30 bool Foam::globalIndexAndTransform::less::operator()
31 (
32  const labelPair& a,
33  const labelPair& b
34 ) const
35 {
36  label procA = globalIndexAndTransform::processor(a);
37  label procB = globalIndexAndTransform::processor(b);
38 
39  if (procA < procB)
40  {
41  return true;
42  }
43  else if (procA > procB)
44  {
45  return false;
46  }
47  else
48  {
49  // Equal proc.
50  label indexA = globalIndexAndTransform::index(a);
51  label indexB = globalIndexAndTransform::index(b);
52 
53  if (indexA < indexB)
54  {
55  return true;
56  }
57  else if (indexA > indexB)
58  {
59  return false;
60  }
61  else
62  {
63  // Equal index
64  label transformA = globalIndexAndTransform::transformIndex(a);
65  label transformB = globalIndexAndTransform::transformIndex(b);
66 
67  return transformA < transformB;
68  }
69  }
70 }
71 
72 
73 Foam::label Foam::globalIndexAndTransform::encodeTransformIndex
74 (
75  const List<label>& permutationIndices
76 ) const
77 {
78  if (permutationIndices.size() != transforms_.size())
79  {
81  (
82  "Foam::label encodeTransformIndex"
83  "("
84  "const List<label>& permutationIndices,"
85  ") const"
86  )
87  << "permutationIndices " << permutationIndices
88  << "are of a different size to the number of independent transforms"
89  << abort(FatalError);
90  }
91 
92  label transformIndex = 0;
93 
94  label w = 1;
95 
96  forAll(transforms_, b)
97  {
98  if (mag(permutationIndices[b]) > 1)
99  {
101  (
102  "Foam::label encodeTransformIndex"
103  "("
104  "const List<label>& permutationIndices,"
105  ") const"
106  )
107  << "permutationIndices " << permutationIndices
108  << "are illegal, they must all be only -1, 0 or +1"
109  << abort(FatalError);
110  }
111 
112  transformIndex += (permutationIndices[b] + 1)*w;
113 
114  w *= 3;
115  }
116 
117  return transformIndex;
118 }
119 
120 
121 Foam::label Foam::globalIndexAndTransform::encodeTransformIndex
122 (
123  const FixedList<Foam::label, 3>& permutation
124 ) const
125 {
126  if (nIndependentTransforms() == 0)
127  {
128  return 0;
129  }
130  if (nIndependentTransforms() == 1)
131  {
132  return permutation[0]+1;
133  }
134  else if (nIndependentTransforms() == 2)
135  {
136  return (permutation[1]+1)*3 + (permutation[0]+1);
137  }
138  else
139  {
140  return
141  (permutation[2]+1)*9
142  + (permutation[1]+1)*3
143  + (permutation[0]+1);
144  }
145 }
146 
147 
149 Foam::globalIndexAndTransform::decodeTransformIndex
150 (
151  const label transformIndex
152 ) const
153 {
154  FixedList<label, 3> permutation(label(0));
155 
156  label t = transformIndex;
157  if (nIndependentTransforms() > 0)
158  {
159  permutation[0] = (t%3)-1;
160  if (nIndependentTransforms() > 1)
161  {
162  t /= 3;
163  permutation[1] = (t%3)-1;
164  if (nIndependentTransforms() > 2)
165  {
166  t /= 3;
167  permutation[2] = (t%3)-1;
168  }
169  }
170  }
171 
172 # ifdef FULLDEBUG
173  t /= 3;
174  if (t != 0)
175  {
177  (
178  "globalIndexAndTransform::decodeTransformIndex(const label)"
179  ) << "transformIndex : " << transformIndex
180  << " has more than 3 fields."
181  << abort(FatalError);
182  }
183 # endif
184 
185  return permutation;
186 }
187 
188 
190 (
191  const label transformIndex,
192  const label patchI,
193  const bool isSendingSide,
194  const scalar tol
195 ) const
196 {
197  const Pair<label>& transSign = patchTransformSign_[patchI];
198 
199  label matchTransI = transSign.first();
200 
201  // Hardcoded for max 3 transforms only!
202 
203  if (matchTransI > -1 && matchTransI < 3)
204  {
205  FixedList<label, 3> permutation = decodeTransformIndex(transformIndex);
206 
207 
208  // Add patch transform
209  // ~~~~~~~~~~~~~~~~~~~
210 
211  label sign = transSign.second();
212  if (!isSendingSide)
213  {
214  sign = -sign;
215  }
216 
217 
218  // If this transform been found already by a patch?
219  if (permutation[matchTransI] != 0)
220  {
221  if (sign == 0)
222  {
223  // sent from patch without a transformation. Do nothing.
224  FatalErrorIn("globalIndexAndTransform::addToTransformIndex(..)")
225  << "patch:" << mesh_.boundaryMesh()[patchI].name()
226  << " transform:" << matchTransI << " sign:" << sign
227  << " current transforms:" << permutation
228  << exit(FatalError);
229  }
230  else if (sign == permutation[matchTransI])
231  {
232  // This is usually illegal. The only exception is for points
233  // on the axis of a 180 degree cyclic wedge when the
234  // transformation is going to be (-1 0 0 0 -1 0 0 0 +1)
235  // (or a different permutation but always two times -1 and
236  // once +1)
237  bool antiCyclic = false;
238 
239  const vectorTensorTransform& vt = transforms_[matchTransI];
240  if (mag(vt.t()) < SMALL && vt.hasR())
241  {
242  const tensor& R = vt.R();
243  scalar sumDiag = tr(R);
244  scalar sumMagDiag = mag(R.xx())+mag(R.yy())+mag(R.zz());
245 
246  if (mag(sumMagDiag-3) < tol && mag(sumDiag+1) < tol)
247  {
248  antiCyclic = true;
249  }
250  }
251 
252  if (antiCyclic)
253  {
254  // 180 degree rotational. Reset transformation.
255  permutation[matchTransI] = 0;
256  }
257  else
258  {
260  (
261  "Foam::label "
262  "Foam::globalIndexAndTransform::addToTransformIndex\n"
263  "(\n"
264  "const label,\n"
265  "const label,\n"
266  "const bool\n"
267  ") const\n"
268  ) << "More than one patch accessing the same transform "
269  << "but not of the same sign." << endl
270  << "patch:" << mesh_.boundaryMesh()[patchI].name()
271  << " transform:" << matchTransI << " sign:" << sign
272  << " current transforms:" << permutation
273  << exit(FatalError);
274  }
275  }
276  else
277  {
278  permutation[matchTransI] = 0;
279  }
280  }
281  else
282  {
283  permutation[matchTransI] = sign;
284  }
285 
286 
287  // Re-encode permutation
288  // ~~~~~~~~~~~~~~~~~~~~~
289 
290  return encodeTransformIndex(permutation);
291  }
292  else
293  {
294  return transformIndex;
295  }
296 }
297 
298 
300 (
301  const label transformIndex0,
302  const label transformIndex1
303 ) const
304 {
305  if (transformIndex0 == transformIndex1)
306  {
307  return transformIndex0;
308  }
309 
310 
311  // Count number of transforms
312  FixedList<label, 3> permutation0 = decodeTransformIndex(transformIndex0);
313  label n0 = 0;
314  forAll(permutation0, i)
315  {
316  if (permutation0[i] != 0)
317  {
318  n0++;
319  }
320  }
321 
322  FixedList<label, 3> permutation1 = decodeTransformIndex(transformIndex1);
323  label n1 = 0;
324  forAll(permutation1, i)
325  {
326  if (permutation1[i] != 0)
327  {
328  n1++;
329  }
330  }
331 
332  if (n0 <= n1)
333  {
334  return transformIndex0;
335  }
336  else
337  {
338  return transformIndex1;
339  }
340 }
341 
342 
344 (
345  const label transformIndex0,
346  const label transformIndex1
347 ) const
348 {
349  FixedList<label, 3> permutation0 = decodeTransformIndex(transformIndex0);
350  FixedList<label, 3> permutation1 = decodeTransformIndex(transformIndex1);
351 
352  forAll(permutation0, i)
353  {
354  permutation0[i] -= permutation1[i];
355  }
356 
357  return encodeTransformIndex(permutation0);
358 }
359 
360 
362 (
363  const label index,
364  const label transformIndex
365 )
366 {
367  return encode(Pstream::myProcNo(), index, transformIndex);
368 }
369 
370 
372 (
373  const label procI,
374  const label index,
375  const label transformIndex
376 )
377 {
378  if (transformIndex < 0 || transformIndex >= base_)
379  {
381  (
382  "Foam::labelPair Foam::globalIndexAndTransform::encode"
383  "("
384  "const label procI, "
385  "const label index, "
386  "const label transformIndex"
387  ")"
388  )
389  << "TransformIndex " << transformIndex
390  << " is outside allowed range of 0 to "
391  << base_ - 1
392  << abort(FatalError);
393  }
394 
395  if (procI > labelMax/base_)
396  {
398  (
399  "Foam::labelPair Foam::globalIndexAndTransform::encode"
400  "("
401  "const label procI, "
402  "const label index, "
403  "const label transformIndex"
404  ")"
405  )
406  << "Overflow : encoding processor " << procI << " in base " << base_
407  << " exceeds capability of label (" << labelMax
408  << "). Please recompile with larger datatype for label."
409  << exit(FatalError);
410  }
411 
412  return labelPair
413  (
414  index,
415  transformIndex + procI*base_
416  );
417 }
418 
419 
421 (
422  const labelPair& globalIAndTransform
423 )
424 {
425  return globalIAndTransform.first();
426 }
427 
428 
430 (
431  const labelPair& globalIAndTransform
432 )
433 {
434  return globalIAndTransform.second()/base_;
435 }
436 
437 
439 (
440  const labelPair& globalIAndTransform
441 )
442 {
443  return globalIAndTransform.second() % base_;
444 }
445 
446 
448 {
449  return transforms_.size();
450 }
451 
452 
455 {
456  return transforms_;
457 }
458 
459 
462 {
463  return transformPermutations_;
464 }
465 
466 
468 {
469  return nullTransformIndex_;
470 }
471 
472 
475 {
476  return patchTransformSign_;
477 }
478 
479 
481 (
482  label transformIndex
483 ) const
484 {
485  return transformPermutations_[transformIndex];
486 }
487 
488 
490 (
491  const labelHashSet& patchIs
492 ) const
493 {
494  List<label> permutation(transforms_.size(), 0);
495 
496  labelList selectedTransformIs(0);
497 
498  if (patchIs.empty() || transforms_.empty())
499  {
500  return selectedTransformIs;
501  }
502 
503  forAllConstIter(labelHashSet, patchIs, iter)
504  {
505  label patchI = iter.key();
506 
507  const Pair<label>& transSign = patchTransformSign_[patchI];
508 
509  label matchTransI = transSign.first();
510 
511  if (matchTransI > -1)
512  {
513  label sign = transSign.second();
514 
515  // If this transform been found already by a patch?
516  if (permutation[matchTransI] != 0)
517  {
518  // If so, if they have opposite signs, then this is
519  // considered an error. They are allowed to be the
520  // same sign, but this only results in a single
521  // transform.
522  if (permutation[matchTransI] != sign)
523  {
525  (
526  "const Foam::List<Foam::vectorTensorTransform>& "
527  "Foam::globalIndexAndTransform::transformsForPatches"
528  "("
529  "const labelList& patchIs"
530  ") const"
531  )
532  << "More than one patch accessing the same transform "
533  << "but not of the same sign."
534  << exit(FatalError);
535  }
536  }
537  else
538  {
539  permutation[matchTransI] = sign;
540  }
541  }
542  }
543 
544  label nUsedTrans = round(sum(mag(permutation)));
545 
546  if (nUsedTrans == 0)
547  {
548  return selectedTransformIs;
549  }
550 
551  // Number of selected transformations
552  label nSelTrans = pow(label(2), nUsedTrans) - 1;
553 
554  // Pout<< nl << permutation << nl << endl;
555 
556  selectedTransformIs.setSize(nSelTrans);
557 
558  switch (nUsedTrans)
559  {
560  case 1:
561  {
562  selectedTransformIs[0] = encodeTransformIndex(permutation);
563 
564  break;
565  }
566  case 2:
567  {
568  List<label> tempPermutation = permutation;
569 
570  label a = 0;
571  label b = 1;
572 
573  // When there are two selected transforms out of three, we
574  // need to choose which of them are being permuted
575  if (transforms_.size() > nUsedTrans)
576  {
577  if (permutation[0] == 0)
578  {
579  a = 1;
580  b = 2;
581  }
582  else if (permutation[1] == 0)
583  {
584  a = 0;
585  b = 2;
586  }
587  else if (permutation[2] == 0)
588  {
589  a = 0;
590  b = 1;
591  }
592  }
593 
594  tempPermutation[a] = a;
595  tempPermutation[b] = permutation[b];
596 
597  selectedTransformIs[0] = encodeTransformIndex(tempPermutation);
598 
599  tempPermutation[a] = permutation[a];
600  tempPermutation[b] = a;
601 
602  selectedTransformIs[1] = encodeTransformIndex(tempPermutation);
603 
604  tempPermutation[a] = permutation[a];
605  tempPermutation[b] = permutation[b];
606 
607  selectedTransformIs[2] = encodeTransformIndex(tempPermutation);
608 
609  break;
610  }
611  case 3:
612  {
613  List<label> tempPermutation = permutation;
614 
615  tempPermutation[0] = 0;
616  tempPermutation[1] = 0;
617  tempPermutation[2] = permutation[2];
618 
619  selectedTransformIs[0] = encodeTransformIndex(tempPermutation);
620 
621  tempPermutation[0] = 0;
622  tempPermutation[1] = permutation[1];
623  tempPermutation[2] = 0;
624 
625  selectedTransformIs[1] = encodeTransformIndex(tempPermutation);
626 
627  tempPermutation[0] = 0;
628  tempPermutation[1] = permutation[1];
629  tempPermutation[2] = permutation[2];
630 
631  selectedTransformIs[2] = encodeTransformIndex(tempPermutation);
632 
633  tempPermutation[0] = permutation[0];
634  tempPermutation[1] = 0;
635  tempPermutation[2] = 0;
636 
637  selectedTransformIs[3] = encodeTransformIndex(tempPermutation);
638 
639  tempPermutation[0] = permutation[0];
640  tempPermutation[1] = 0;
641  tempPermutation[2] = permutation[2];
642 
643  selectedTransformIs[4] = encodeTransformIndex(tempPermutation);
644 
645  tempPermutation[0] = permutation[0];
646  tempPermutation[1] = permutation[1];
647  tempPermutation[2] = 0;
648 
649  selectedTransformIs[5] = encodeTransformIndex(tempPermutation);
650 
651  tempPermutation[0] = permutation[0];
652  tempPermutation[1] = permutation[1];
653  tempPermutation[2] = permutation[2];
654 
655  selectedTransformIs[6] = encodeTransformIndex(tempPermutation);
656 
657  break;
658  }
659  default:
660  {
662  (
663  "const Foam::List<Foam::vectorTensorTransform>& "
664  "Foam::globalIndexAndTransform::transformsForPatches"
665  "("
666  "const labelList& patchIs"
667  ") const"
668  )
669  << "Only 1-3 transforms are possible."
670  << exit(FatalError);
671  }
672  }
673 
674  return selectedTransformIs;
675 }
676 
677 
679 (
680  const labelHashSet& patchIs,
681  const point& pt
682 ) const
683 {
684  labelList transIs = transformIndicesForPatches(patchIs);
685 
686  // Pout<< patchIs << nl << transIs << endl;
687 
688  pointField transPts(transIs.size());
689 
690  forAll(transIs, tII)
691  {
692  transPts[tII] = transformPermutations_[transIs[tII]].transformPosition
693  (
694  pt
695  );
696  }
697 
698  return transPts;
699 }
700 
701 
702 // ************************************************************************* //
labelList transformIndicesForPatches(const labelHashSet &patchIs) const
Access the all of the indices of the transform.
label nIndependentTransforms() const
Return the number of independent transforms.
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: contiguous.H:49
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:28
const Cmpt & zz() const
Definition: TensorI.H:216
dimensioned< scalar > mag(const dimensioned< Type > &)
const List< vectorTensorTransform > & transforms() const
Return access to the stored independent transforms.
const List< vectorTensorTransform > & transformPermutations() const
Return access to the permuted transforms.
#define R(A, B, C, D, E, F, K, M)
const Type & second() const
Return second.
Definition: Pair.H:99
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
const Cmpt & yy() const
Definition: TensorI.H:188
dimensionedScalar sign(const dimensionedScalar &ds)
const List< Pair< label > > & patchTransformSign() const
Return access to the per-patch transform-sign pairs.
A 1D vector of objects of type <T> with a fixed size <Size>.
Definition: FixedList.H:53
label addToTransformIndex(const label transformIndex, const label patchI, const bool isSendingSide=true, const scalar tol=SMALL) const
Add patch transformation to transformIndex. Return new.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
void setSize(const label)
Reset size of List.
Definition: List.C:318
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const vectorTensorTransform & transform(label transformIndex) const
Access the overall (permuted) transform corresponding.
#define forAll(list, i)
Definition: UList.H:421
label minimumTransformIndex(const label transformIndex0, const label transformIndex1) const
Combine two transformIndices.
label nullTransformIndex() const
Return the transformIndex (index in transformPermutations)
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const word & name() const
Return name.
Definition: IOobject.H:260
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:404
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Vector-tensor class used to perform translations and rotations in 3D space.
error FatalError
Pair< label > labelPair
Label pair.
Definition: labelPair.H:48
const Type & first() const
Return first.
Definition: Pair.H:87
label subtractTransformIndex(const label transformIndex0, const label transformIndex1) const
Subtract two transformIndices.
static label index(const labelPair &globalIAndTransform)
Index carried by the object.
bool empty() const
Return true if the hash table is empty.
Definition: HashTableI.H:72
const Cmpt & xx() const
Definition: TensorI.H:160
pointField transformPatches(const labelHashSet &patchIs, const point &pt) const
Apply all of the transform permutations.
static label transformIndex(const labelPair &globalIAndTransform)
Transform carried by the object.
static labelPair encode(const label index, const label transformIndex)
Encode index and bare index as components on own processor.
static label processor(const labelPair &globalIAndTransform)
Which processor does this come from?
static const label labelMax
Definition: label.H:62