cyclicAMIPolyPatch.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 | Copyright (C) 2011-2015 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 "cyclicAMIPolyPatch.H"
27 #include "SubField.H"
28 #include "Time.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35  defineTypeNameAndDebug(cyclicAMIPolyPatch, 0);
36 
37  addToRunTimeSelectionTable(polyPatch, cyclicAMIPolyPatch, word);
38  addToRunTimeSelectionTable(polyPatch, cyclicAMIPolyPatch, dictionary);
39 }
40 
41 
42 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
43 
44 Foam::vector Foam::cyclicAMIPolyPatch::findFaceNormalMaxRadius
45 (
46  const pointField& faceCentres
47 ) const
48 {
49  // Determine a face furthest away from the axis
50 
51  const vectorField n((faceCentres - rotationCentre_) ^ rotationAxis_);
52 
53  const scalarField magRadSqr(magSqr(n));
54 
55  label faceI = findMax(magRadSqr);
56 
57  if (debug)
58  {
59  Info<< "findFaceMaxRadius(const pointField&) : patch: " << name() << nl
60  << " rotFace = " << faceI << nl
61  << " point = " << faceCentres[faceI] << nl
62  << " distance = " << Foam::sqrt(magRadSqr[faceI])
63  << endl;
64  }
65 
66  return n[faceI];
67 }
68 
69 
70 // * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * * //
71 
73 {
74  const cyclicAMIPolyPatch& half0 = *this;
75  vectorField half0Areas(half0.size());
76  forAll(half0, facei)
77  {
78  half0Areas[facei] = half0[facei].normal(half0.points());
79  }
80 
81  const cyclicAMIPolyPatch& half1 = neighbPatch();
82  vectorField half1Areas(half1.size());
83  forAll(half1, facei)
84  {
85  half1Areas[facei] = half1[facei].normal(half1.points());
86  }
87 
88  calcTransforms
89  (
90  half0,
91  half0.faceCentres(),
92  half0Areas,
93  half1.faceCentres(),
94  half1Areas
95  );
96 
97  if (debug)
98  {
99  Pout<< "calcTransforms() : patch: " << name() << nl
100  << " forwardT = " << forwardT() << nl
101  << " reverseT = " << reverseT() << nl
102  << " separation = " << separation() << nl
103  << " collocated = " << collocated() << nl << endl;
104  }
105 }
106 
107 
109 (
110  const primitivePatch& half0,
111  const pointField& half0Ctrs,
112  const vectorField& half0Areas,
113  const pointField& half1Ctrs,
114  const vectorField& half1Areas
115 )
116 {
117  if (transform() != neighbPatch().transform())
118  {
119  FatalErrorIn("cyclicAMIPolyPatch::calcTransforms()")
120  << "Patch " << name()
121  << " has transform type " << transformTypeNames[transform()]
122  << ", neighbour patch " << neighbPatchName()
123  << " has transform type "
124  << neighbPatch().transformTypeNames[neighbPatch().transform()]
125  << exit(FatalError);
126  }
127 
128 
129  // Calculate transformation tensors
130 
131  switch (transform())
132  {
133  case ROTATIONAL:
134  {
135  tensor revT = tensor::zero;
136 
137  if (rotationAngleDefined_)
138  {
139  const tensor T(rotationAxis_*rotationAxis_);
140 
141  const tensor S
142  (
143  0, -rotationAxis_.z(), rotationAxis_.y(),
144  rotationAxis_.z(), 0, -rotationAxis_.x(),
145  -rotationAxis_.y(), rotationAxis_.x(), 0
146  );
147 
148  const tensor revTPos
149  (
150  T
151  + cos(rotationAngle_)*(tensor::I - T)
152  + sin(rotationAngle_)*S
153  );
154 
155  const tensor revTNeg
156  (
157  T
158  + cos(-rotationAngle_)*(tensor::I - T)
159  + sin(-rotationAngle_)*S
160  );
161 
162  // Check - assume correct angle when difference in face areas
163  // is the smallest
164  const vector transformedAreaPos = gSum(half1Areas & revTPos);
165  const vector transformedAreaNeg = gSum(half1Areas & revTNeg);
166  const vector area0 = gSum(half0Areas);
167  const scalar magArea0 = mag(area0) + ROOTVSMALL;
168 
169  // Areas have opposite sign, so sum should be zero when correct
170  // rotation applied
171  const scalar errorPos = mag(transformedAreaPos + area0);
172  const scalar errorNeg = mag(transformedAreaNeg + area0);
173 
174  const scalar normErrorPos = errorPos/magArea0;
175  const scalar normErrorNeg = errorNeg/magArea0;
176 
177  if (errorPos > errorNeg && normErrorNeg < matchTolerance())
178  {
179  revT = revTNeg;
180  rotationAngle_ *= -1;
181  }
182  else
183  {
184  revT = revTPos;
185  }
186 
187  const scalar areaError = min(normErrorPos, normErrorNeg);
188 
189  if (areaError > matchTolerance())
190  {
191  WarningIn
192  (
193  "void Foam::cyclicAMIPolyPatch::calcTransforms"
194  "("
195  "const primitivePatch&, "
196  "const pointField&, "
197  "const vectorField&, "
198  "const pointField&, "
199  "const vectorField&"
200  ")"
201  ) << "Patch areas are not consistent within "
202  << 100*matchTolerance()
203  << " % indicating a possible error in the specified "
204  << "angle of rotation" << nl
205  << " owner patch : " << name() << nl
206  << " neighbour patch : " << neighbPatch().name()
207  << nl
208  << " angle : "
209  << radToDeg(rotationAngle_) << " deg" << nl
210  << " area error : " << 100*areaError << " %"
211  << " match tolerance : " << matchTolerance()
212  << endl;
213  }
214 
215  if (debug)
216  {
217  scalar theta = radToDeg(rotationAngle_);
218 
219  Pout<< "cyclicAMIPolyPatch::calcTransforms: patch:"
220  << name()
221  << " Specified rotation:"
222  << " swept angle: " << theta << " [deg]"
223  << " reverse transform: " << revT
224  << endl;
225  }
226  }
227  else
228  {
229  point n0 = vector::zero;
230  point n1 = vector::zero;
231  if (half0Ctrs.size())
232  {
233  n0 = findFaceNormalMaxRadius(half0Ctrs);
234  }
235  if (half1Ctrs.size())
236  {
237  n1 = -findFaceNormalMaxRadius(half1Ctrs);
238  }
239 
240  reduce(n0, maxMagSqrOp<point>());
241  reduce(n1, maxMagSqrOp<point>());
242 
243  n0 /= mag(n0) + VSMALL;
244  n1 /= mag(n1) + VSMALL;
245 
246  // Extended tensor from two local coordinate systems calculated
247  // using normal and rotation axis
248  const tensor E0
249  (
250  rotationAxis_,
251  (n0 ^ rotationAxis_),
252  n0
253  );
254  const tensor E1
255  (
256  rotationAxis_,
257  (-n1 ^ rotationAxis_),
258  -n1
259  );
260  revT = E1.T() & E0;
261 
262  if (debug)
263  {
264  scalar theta = radToDeg(acos(-(n0 & n1)));
265 
266  Pout<< "cyclicAMIPolyPatch::calcTransforms: patch:"
267  << name()
268  << " Specified rotation:"
269  << " n0:" << n0 << " n1:" << n1
270  << " swept angle: " << theta << " [deg]"
271  << " reverse transform: " << revT
272  << endl;
273  }
274  }
275 
276  const_cast<tensorField&>(forwardT()) = tensorField(1, revT.T());
277  const_cast<tensorField&>(reverseT()) = tensorField(1, revT);
278  const_cast<vectorField&>(separation()).setSize(0);
279  const_cast<boolList&>(collocated()) = boolList(1, false);
280 
281  break;
282  }
283  case TRANSLATIONAL:
284  {
285  if (debug)
286  {
287  Pout<< "cyclicAMIPolyPatch::calcTransforms : patch:" << name()
288  << " Specified translation : " << separationVector_
289  << endl;
290  }
291 
292  const_cast<tensorField&>(forwardT()).clear();
293  const_cast<tensorField&>(reverseT()).clear();
294  const_cast<vectorField&>(separation()) = vectorField
295  (
296  1,
297  separationVector_
298  );
299  const_cast<boolList&>(collocated()) = boolList(1, false);
300 
301  break;
302  }
303  default:
304  {
305  if (debug)
306  {
307  Pout<< "patch:" << name()
308  << " Assuming cyclic AMI pairs are colocated" << endl;
309  }
310 
311  const_cast<tensorField&>(forwardT()).clear();
312  const_cast<tensorField&>(reverseT()).clear();
313  const_cast<vectorField&>(separation()).setSize(0);
314  const_cast<boolList&>(collocated()) = boolList(1, true);
315 
316  break;
317  }
318  }
319 
320  if (debug)
321  {
322  Pout<< "patch: " << name() << nl
323  << " forwardT = " << forwardT() << nl
324  << " reverseT = " << reverseT() << nl
325  << " separation = " << separation() << nl
326  << " collocated = " << collocated() << nl << endl;
327  }
328 }
329 
330 
331 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
332 
334 (
336 ) const
337 {
338  if (owner())
339  {
340  AMIPtr_.clear();
341 
342  const polyPatch& nbr = neighbPatch();
343  pointField nbrPoints
344  (
345  neighbPatch().boundaryMesh().mesh().points(),
346  neighbPatch().meshPoints()
347  );
348 
349  if (debug)
350  {
351  const Time& t = boundaryMesh().mesh().time();
352  OFstream os(t.path()/name() + "_neighbourPatch-org.obj");
353  meshTools::writeOBJ(os, neighbPatch().localFaces(), nbrPoints);
354  }
355 
356  // Transform neighbour patch to local system
357  transformPosition(nbrPoints);
358  primitivePatch nbrPatch0
359  (
361  (
362  nbr.localFaces(),
363  nbr.size()
364  ),
365  nbrPoints
366  );
367 
368  if (debug)
369  {
370  const Time& t = boundaryMesh().mesh().time();
371  OFstream osN(t.path()/name() + "_neighbourPatch-trans.obj");
372  meshTools::writeOBJ(osN, nbrPatch0.localFaces(), nbrPoints);
373 
374  OFstream osO(t.path()/name() + "_ownerPatch.obj");
375  meshTools::writeOBJ(osO, this->localFaces(), localPoints());
376  }
377 
378  // Construct/apply AMI interpolation to determine addressing and weights
379  AMIPtr_.reset
380  (
382  (
383  *this,
384  nbrPatch0,
385  surfPtr(),
387  AMIRequireMatch_,
388  AMIMethod,
389  AMILowWeightCorrection_,
390  AMIReverse_
391  )
392  );
393 
394  if (debug)
395  {
396  Pout<< "cyclicAMIPolyPatch : " << name()
397  << " constructed AMI with " << nl
398  << " " << "srcAddress:" << AMIPtr_().srcAddress().size()
399  << nl
400  << " " << "tgAddress :" << AMIPtr_().tgtAddress().size()
401  << nl << endl;
402  }
403  }
404 }
405 
406 
408 {
409  // Clear the invalid AMI
410  AMIPtr_.clear();
411 
413 }
414 
415 
417 {
418  calcGeometry
419  (
420  *this,
421  faceCentres(),
422  faceAreas(),
423  faceCellCentres(),
424  neighbPatch().faceCentres(),
425  neighbPatch().faceAreas(),
426  neighbPatch().faceCellCentres()
427  );
428 }
429 
430 
432 (
433  PstreamBuffers& pBufs,
434  const pointField& p
435 )
436 {
437  // Clear the invalid AMI
438  AMIPtr_.clear();
439 
440  polyPatch::initMovePoints(pBufs, p);
441 
442  // See below. Clear out any local geometry
444 }
445 
446 
448 (
449  PstreamBuffers& pBufs,
450  const pointField& p
451 )
452 {
453  polyPatch::movePoints(pBufs, p);
454 
455  calcTransforms();
456 }
457 
458 
460 {
461  // Clear the invalid AMI
462  AMIPtr_.clear();
463 
465 }
466 
467 
469 {
470  polyPatch::updateMesh(pBufs);
471 }
472 
473 
475 {
476  AMIPtr_.clear();
478 }
479 
480 
481 // * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * * //
482 
484 (
485  const word& name,
486  const label size,
487  const label start,
488  const label index,
489  const polyBoundaryMesh& bm,
490  const word& patchType,
491  const transformType transform
492 )
493 :
494  coupledPolyPatch(name, size, start, index, bm, patchType, transform),
495  nbrPatchName_(word::null),
496  nbrPatchID_(-1),
497  rotationAxis_(vector::zero),
498  rotationCentre_(point::zero),
499  rotationAngleDefined_(false),
500  rotationAngle_(0.0),
501  separationVector_(vector::zero),
502  AMIPtr_(NULL),
503  AMIReverse_(false),
504  AMIRequireMatch_(true),
505  AMILowWeightCorrection_(-1.0),
506  surfPtr_(NULL),
507  surfDict_(fileName("surface"))
508 {
509  // Neighbour patch might not be valid yet so no transformation
510  // calculation possible
511 }
512 
513 
515 (
516  const word& name,
517  const dictionary& dict,
518  const label index,
519  const polyBoundaryMesh& bm,
520  const word& patchType
521 )
522 :
523  coupledPolyPatch(name, dict, index, bm, patchType),
524  nbrPatchName_(dict.lookupOrDefault<word>("neighbourPatch", "")),
525  coupleGroup_(dict),
526  nbrPatchID_(-1),
527  rotationAxis_(vector::zero),
528  rotationCentre_(point::zero),
529  rotationAngleDefined_(false),
530  rotationAngle_(0.0),
531  separationVector_(vector::zero),
532  AMIPtr_(NULL),
533  AMIReverse_(dict.lookupOrDefault<bool>("flipNormals", false)),
534  AMIRequireMatch_(true),
535  AMILowWeightCorrection_(dict.lookupOrDefault("lowWeightCorrection", -1.0)),
536  surfPtr_(NULL),
537  surfDict_(dict.subOrEmptyDict("surface"))
538 {
539  if (nbrPatchName_ == word::null && !coupleGroup_.valid())
540  {
542  (
543  "cyclicAMIPolyPatch::cyclicAMIPolyPatch"
544  "("
545  "const word&, "
546  "const dictionary&, "
547  "const label, "
548  "const polyBoundaryMesh&"
549  ")",
550  dict
551  ) << "No \"neighbourPatch\" or \"coupleGroup\" provided."
552  << exit(FatalIOError);
553  }
554 
555  if (nbrPatchName_ == name)
556  {
558  (
559  "cyclicAMIPolyPatch::cyclicAMIPolyPatch"
560  "("
561  "const word&, "
562  "const dictionary&, "
563  "const label, "
564  "const polyBoundaryMesh&"
565  ")",
566  dict
567  ) << "Neighbour patch name " << nbrPatchName_
568  << " cannot be the same as this patch " << name
569  << exit(FatalIOError);
570  }
571 
572  switch (transform())
573  {
574  case ROTATIONAL:
575  {
576  dict.lookup("rotationAxis") >> rotationAxis_;
577  dict.lookup("rotationCentre") >> rotationCentre_;
578  if (dict.readIfPresent("rotationAngle", rotationAngle_))
579  {
580  rotationAngleDefined_ = true;
581  rotationAngle_ = degToRad(rotationAngle_);
582 
583  if (debug)
584  {
585  Info<< "rotationAngle: " << rotationAngle_ << " [rad]"
586  << endl;
587  }
588  }
589 
590  scalar magRot = mag(rotationAxis_);
591  if (magRot < SMALL)
592  {
594  (
595  "cyclicAMIPolyPatch::cyclicAMIPolyPatch"
596  "("
597  "const word&, "
598  "const dictionary&, "
599  "const label, "
600  "const polyBoundaryMesh&"
601  ")",
602  dict
603  ) << "Illegal rotationAxis " << rotationAxis_ << endl
604  << "Please supply a non-zero vector."
605  << exit(FatalIOError);
606  }
607  rotationAxis_ /= magRot;
608 
609  break;
610  }
611  case TRANSLATIONAL:
612  {
613  dict.lookup("separationVector") >> separationVector_;
614  break;
615  }
616  default:
617  {
618  // No additional info required
619  }
620  }
621 
622  // Neighbour patch might not be valid yet so no transformation
623  // calculation possible
624 }
625 
626 
628 (
629  const cyclicAMIPolyPatch& pp,
630  const polyBoundaryMesh& bm
631 )
632 :
633  coupledPolyPatch(pp, bm),
634  nbrPatchName_(pp.nbrPatchName_),
635  coupleGroup_(pp.coupleGroup_),
636  nbrPatchID_(-1),
637  rotationAxis_(pp.rotationAxis_),
638  rotationCentre_(pp.rotationCentre_),
639  rotationAngleDefined_(pp.rotationAngleDefined_),
640  rotationAngle_(pp.rotationAngle_),
641  separationVector_(pp.separationVector_),
642  AMIPtr_(NULL),
643  AMIReverse_(pp.AMIReverse_),
644  AMIRequireMatch_(pp.AMIRequireMatch_),
645  AMILowWeightCorrection_(pp.AMILowWeightCorrection_),
646  surfPtr_(NULL),
647  surfDict_(pp.surfDict_)
648 {
649  // Neighbour patch might not be valid yet so no transformation
650  // calculation possible
651 }
652 
653 
655 (
656  const cyclicAMIPolyPatch& pp,
657  const polyBoundaryMesh& bm,
658  const label index,
659  const label newSize,
660  const label newStart,
661  const word& nbrPatchName
662 )
663 :
664  coupledPolyPatch(pp, bm, index, newSize, newStart),
665  nbrPatchName_(nbrPatchName),
666  coupleGroup_(pp.coupleGroup_),
667  nbrPatchID_(-1),
668  rotationAxis_(pp.rotationAxis_),
669  rotationCentre_(pp.rotationCentre_),
670  rotationAngleDefined_(pp.rotationAngleDefined_),
671  rotationAngle_(pp.rotationAngle_),
672  separationVector_(pp.separationVector_),
673  AMIPtr_(NULL),
674  AMIReverse_(pp.AMIReverse_),
675  AMIRequireMatch_(pp.AMIRequireMatch_),
676  AMILowWeightCorrection_(pp.AMILowWeightCorrection_),
677  surfPtr_(NULL),
678  surfDict_(pp.surfDict_)
679 {
680  if (nbrPatchName_ == name())
681  {
683  (
684  "const cyclicAMIPolyPatch& "
685  "const polyBoundaryMesh&, "
686  "const label, "
687  "const label, "
688  "const label, "
689  "const word&"
690  ) << "Neighbour patch name " << nbrPatchName_
691  << " cannot be the same as this patch " << name()
692  << exit(FatalError);
693  }
694 
695  // Neighbour patch might not be valid yet so no transformation
696  // calculation possible
697 }
698 
699 
701 (
702  const cyclicAMIPolyPatch& pp,
703  const polyBoundaryMesh& bm,
704  const label index,
705  const labelUList& mapAddressing,
706  const label newStart
707 )
708 :
709  coupledPolyPatch(pp, bm, index, mapAddressing, newStart),
710  nbrPatchName_(pp.nbrPatchName_),
711  coupleGroup_(pp.coupleGroup_),
712  nbrPatchID_(-1),
713  rotationAxis_(pp.rotationAxis_),
714  rotationCentre_(pp.rotationCentre_),
715  rotationAngleDefined_(pp.rotationAngleDefined_),
716  rotationAngle_(pp.rotationAngle_),
717  separationVector_(pp.separationVector_),
718  AMIPtr_(NULL),
719  AMIReverse_(pp.AMIReverse_),
720  AMIRequireMatch_(pp.AMIRequireMatch_),
721  AMILowWeightCorrection_(pp.AMILowWeightCorrection_),
722  surfPtr_(NULL),
723  surfDict_(pp.surfDict_)
724 {}
725 
726 
727 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
728 
730 {}
731 
732 
733 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
734 
736 {
737  if (nbrPatchID_ == -1)
738  {
739  nbrPatchID_ = this->boundaryMesh().findPatchID(neighbPatchName());
740 
741  if (nbrPatchID_ == -1)
742  {
743  FatalErrorIn("cyclicPolyAMIPatch::neighbPatchID() const")
744  << "Illegal neighbourPatch name " << neighbPatchName()
745  << nl << "Valid patch names are "
746  << this->boundaryMesh().names()
747  << exit(FatalError);
748  }
749 
750  // Check that it is a cyclic AMI patch
751  const cyclicAMIPolyPatch& nbrPatch =
752  refCast<const cyclicAMIPolyPatch>
753  (
754  this->boundaryMesh()[nbrPatchID_]
755  );
756 
757  if (nbrPatch.neighbPatchName() != name())
758  {
759  WarningIn("cyclicAMIPolyPatch::neighbPatchID() const")
760  << "Patch " << name()
761  << " specifies neighbour patch " << neighbPatchName()
762  << nl << " but that in return specifies "
763  << nbrPatch.neighbPatchName() << endl;
764  }
765  }
766 
767  return nbrPatchID_;
768 }
769 
770 
772 {
773  return index() < neighbPatchID();
774 }
775 
776 
778 {
779  const polyPatch& pp = this->boundaryMesh()[neighbPatchID()];
780  return refCast<const cyclicAMIPolyPatch>(pp);
781 }
782 
783 
786 {
787  const word surfType(surfDict_.lookupOrDefault<word>("type", "none"));
788 
789  if (!surfPtr_.valid() && owner() && surfType != "none")
790  {
791  word surfName(surfDict_.lookupOrDefault("name", name()));
792 
793  const polyMesh& mesh = boundaryMesh().mesh();
794 
795  surfPtr_ =
797  (
798  surfType,
799  IOobject
800  (
801  surfName,
802  mesh.time().constant(),
803  "triSurface",
804  mesh,
807  ),
808  surfDict_
809  );
810  }
811 
812  return surfPtr_;
813 }
814 
815 
817 {
818  if (!owner())
819  {
821  (
822  "const AMIPatchToPatchInterpolation& cyclicAMIPolyPatch::AMI()"
823  )
824  << "AMI interpolator only available to owner patch"
825  << abort(FatalError);
826  }
827 
828  if (!AMIPtr_.valid())
829  {
830  resetAMI();
831  }
832 
833  return AMIPtr_();
834 }
835 
836 
838 {
839  if (owner())
840  {
841  return AMI().applyLowWeightCorrection();
842  }
843  else
844  {
845  return neighbPatch().AMI().applyLowWeightCorrection();
846  }
847 }
848 
849 
851 {
852  if (!parallel())
853  {
854  if (transform() == ROTATIONAL)
855  {
856  l = Foam::transform(forwardT(), l - rotationCentre_)
857  + rotationCentre_;
858  }
859  else
860  {
861  l = Foam::transform(forwardT(), l);
862  }
863  }
864  else if (separated())
865  {
866  // transformPosition gets called on the receiving side,
867  // separation gets calculated on the sending side so subtract
868 
869  const vectorField& s = separation();
870  if (s.size() == 1)
871  {
872  forAll(l, i)
873  {
874  l[i] -= s[0];
875  }
876  }
877  else
878  {
879  l -= s;
880  }
881  }
882 }
883 
884 
886 (
887  point& l,
888  const label faceI
889 ) const
890 {
891  if (!parallel())
892  {
893  const tensor& T =
894  (
895  forwardT().size() == 1
896  ? forwardT()[0]
897  : forwardT()[faceI]
898  );
899 
900  if (transform() == ROTATIONAL)
901  {
902  l = Foam::transform(T, l - rotationCentre_) + rotationCentre_;
903  }
904  else
905  {
906  l = Foam::transform(T, l);
907  }
908  }
909  else if (separated())
910  {
911  const vector& s =
912  (
913  separation().size() == 1
914  ? separation()[0]
915  : separation()[faceI]
916  );
917 
918  l -= s;
919  }
920 }
921 
922 
924 (
925  point& l,
926  const label faceI
927 ) const
928 {
929  if (!parallel())
930  {
931  const tensor& T =
932  (
933  reverseT().size() == 1
934  ? reverseT()[0]
935  : reverseT()[faceI]
936  );
937 
938  if (transform() == ROTATIONAL)
939  {
940  l = Foam::transform(T, l - rotationCentre_) + rotationCentre_;
941  }
942  else
943  {
944  l = Foam::transform(T, l);
945  }
946  }
947  else if (separated())
948  {
949  const vector& s =
950  (
951  separation().size() == 1
952  ? separation()[0]
953  : separation()[faceI]
954  );
955 
956  l += s;
957  }
958 }
959 
960 
962 (
963  vector& d,
964  const label faceI
965 ) const
966 {
967  if (!parallel())
968  {
969  const tensor& T =
970  (
971  reverseT().size() == 1
972  ? reverseT()[0]
973  : reverseT()[faceI]
974  );
975 
976  d = Foam::transform(T, d);
977  }
978 }
979 
980 
982 (
983  const primitivePatch& referPatch,
984  const pointField& thisCtrs,
985  const vectorField& thisAreas,
986  const pointField& thisCc,
987  const pointField& nbrCtrs,
988  const vectorField& nbrAreas,
989  const pointField& nbrCc
990 )
991 {
992  calcTransforms
993  (
994  referPatch,
995  thisCtrs,
996  thisAreas,
997  nbrCtrs,
998  nbrAreas
999  );
1000 }
1001 
1002 
1005  PstreamBuffers& pBufs,
1006  const primitivePatch& pp
1007 ) const
1008 {}
1009 
1010 
1013  PstreamBuffers& pBufs,
1014  const primitivePatch& pp,
1015  labelList& faceMap,
1016  labelList& rotation
1017 ) const
1018 {
1019  faceMap.setSize(pp.size());
1020  faceMap = -1;
1021 
1022  rotation.setSize(pp.size());
1023  rotation = 0;
1024 
1025  return false;
1026 }
1027 
1028 
1031  const label faceI,
1032  const vector& n,
1033  point& p
1034 ) const
1035 {
1036  point prt(p);
1037  reverseTransformPosition(prt, faceI);
1038 
1039  vector nrt(n);
1040  reverseTransformDirection(nrt, faceI);
1041 
1042  label nbrFaceI = -1;
1043 
1044  if (owner())
1045  {
1046  nbrFaceI = AMI().tgtPointFace
1047  (
1048  *this,
1049  neighbPatch(),
1050  nrt,
1051  faceI,
1052  prt
1053  );
1054  }
1055  else
1056  {
1057  nbrFaceI = neighbPatch().AMI().srcPointFace
1058  (
1059  neighbPatch(),
1060  *this,
1061  nrt,
1062  faceI,
1063  prt
1064  );
1065  }
1066 
1067  if (nbrFaceI >= 0)
1068  {
1069  p = prt;
1070  }
1071 
1072  return nbrFaceI;
1073 }
1074 
1075 
1077 {
1079  if (!nbrPatchName_.empty())
1080  {
1081  os.writeKeyword("neighbourPatch") << nbrPatchName_
1082  << token::END_STATEMENT << nl;
1083  }
1084  coupleGroup_.write(os);
1085 
1086  switch (transform())
1087  {
1088  case ROTATIONAL:
1089  {
1090  os.writeKeyword("rotationAxis") << rotationAxis_
1091  << token::END_STATEMENT << nl;
1092  os.writeKeyword("rotationCentre") << rotationCentre_
1093  << token::END_STATEMENT << nl;
1094 
1095  if (rotationAngleDefined_)
1096  {
1097  os.writeKeyword("rotationAngle") << radToDeg(rotationAngle_)
1098  << token::END_STATEMENT << nl;
1099  }
1100 
1101  break;
1102  }
1103  case TRANSLATIONAL:
1104  {
1105  os.writeKeyword("separationVector") << separationVector_
1106  << token::END_STATEMENT << nl;
1107  break;
1108  }
1109  case NOORDERING:
1110  {
1111  break;
1112  }
1113  default:
1114  {
1115  // No additional info to write
1116  }
1117  }
1118 
1119  if (AMIReverse_)
1120  {
1121  os.writeKeyword("flipNormals") << AMIReverse_
1122  << token::END_STATEMENT << nl;
1123  }
1124 
1125  if (AMILowWeightCorrection_ > 0)
1126  {
1127  os.writeKeyword("lowWeightCorrection") << AMILowWeightCorrection_
1128  << token::END_STATEMENT << nl;
1129  }
1130 
1131  if (!surfDict_.empty())
1132  {
1133  os.writeKeyword(surfDict_.dictName());
1134  os << surfDict_;
1135  }
1136 }
1137 
1138 
1139 // ************************************************************************* //
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
dimensionedScalar sqrt(const dimensionedScalar &ds)
Output to file stream.
Definition: OFstream.H:81
virtual bool owner() const
Does this side own the patch?
const pointField & points
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
Buffers for inter-processor communications streams (UOPstream, UIPstream).
bool rotationAngleDefined_
Flag to show whether the rotation angle is defined.
const autoPtr< searchableSurface > & surfPtr() const
Return a reference to the projection surface.
Cyclic patch for Arbitrary Mesh Interface (AMI)
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
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){word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject( name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE ))
virtual void reverseTransformPosition(point &l, const label faceI) const
Transform a patch-based position from this side to nbr side.
fileName path() const
Return path.
Definition: Time.H:281
dimensioned< scalar > mag(const dimensioned< Type > &)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
virtual void clearGeom()
Clear geometry.
virtual void initMovePoints(PstreamBuffers &pBufs, const pointField &)
Initialise the patches for moving points.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
vector separationVector_
Translation vector.
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
static const Tensor zero
Definition: Tensor.H:80
virtual void initUpdateMesh(PstreamBuffers &)
Initialise the update of the patch topology.
cyclicAMIPolyPatch(const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm, const word &patchType, const transformType transform=UNKNOWN)
Construct from (base couped patch) components.
virtual ~cyclicAMIPolyPatch()
Destructor.
label findMax(const ListType &, const label start=0)
Find index of max element (and larger than given element).
A class for handling words, derived from string.
Definition: word.H:59
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
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface...
Definition: boundaryMesh.H:59
virtual void movePoints(PstreamBuffers &pBufs, const pointField &)
Correct patches after moving points.
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
virtual void initMovePoints(PstreamBuffers &, const pointField &)
Initialise the patches for moving points.
Definition: polyPatch.H:108
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
Definition: polyPatch.C:62
messageStream Info
virtual label neighbPatchID() const
Neighbour patch ID.
dynamicFvMesh & mesh
bool AMIRequireMatch_
Flag to indicate that patches should match/overlap.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
word nbrPatchName_
Name of other half.
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
Namespace for OpenFOAM.
Field< tensor > tensorField
Specialisation of Field<T> for tensor.
scalar rotationAngle_
Rotation angle.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const vectorField::subField faceCentres() const
Return face centres.
Definition: polyPatch.C:284
Type gSum(const FieldField< Field, Type > &f)
label n
interpolationMethod
Enumeration specifying interpolation method.
point rotationCentre_
Point on axis of rotation for rotational cyclics.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
virtual const cyclicAMIPolyPatch & neighbPatch() const
Return a reference to the neighbour patch.
static const char nl
Definition: Ostream.H:260
static autoPtr< searchableSurface > New(const word &surfaceType, const IOobject &io, const dictionary &dict)
Return a reference to the selected searchableSurface.
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
IOerror FatalIOError
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
#define WarningIn(functionName)
Report a warning using Foam::Warning.
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
Definition: polyPatch.H:100
virtual void transformPosition(pointField &) const
Transform patch-based positions from nbr side to this side.
Foam::polyBoundaryMesh.
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:465
const coupleGroupIdentifier coupleGroup_
Optional patchGroup to find neighbPatch.
static bool valid(char)
Is this character valid for a word.
Definition: wordI.H:117
const bMesh & mesh() const
Definition: boundaryMesh.H:202
virtual void calcTransforms()
Recalculate the transformation tensors.
Field< vector > vectorField
Specialisation of Field<T> for vector.
Interpolation class dealing with transfer of data between two primitive patches with an arbitrary mes...
const scalar AMILowWeightCorrection_
Low weight correction threshold for AMI.
#define forAll(list, i)
Definition: UList.H:421
label pointFace(const label faceI, const vector &n, point &p) const
Return face index on neighbour patch which shares point p.
dimensionedScalar cos(const dimensionedScalar &ds)
scalar radToDeg(const scalar rad)
Conversion from radians to degrees.
Macros for easy insertion into run-time selection tables.
dimensionedScalar acos(const dimensionedScalar &ds)
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:59
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
A list of faces which address into the list of points.
virtual void clearGeom()
Clear geometry.
Definition: polyPatch.C:69
bool applyLowWeightCorrection() const
Return true if applying the low weight correction.
points setSize(newPointi)
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
const AMIPatchToPatchInterpolation & AMI() const
Return a reference to the AMI interpolator.
virtual void resetAMI(const AMIPatchToPatchInterpolation::interpolationMethod &AMIMethod=AMIPatchToPatchInterpolation::imFaceAreaWeight) const
Reset the AMI interpolator.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
error FatalError
const bool AMIReverse_
Flag to indicate that slave patch should be reversed for AMI.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
A class for handling file names.
Definition: fileName.H:69
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
static const Vector zero
Definition: Vector.H:80
scalar degToRad(const scalar deg)
Conversion from degrees to radians.
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
virtual void initOrder(PstreamBuffers &, const primitivePatch &) const
Initialize ordering for primitivePatch. Does not.
Tensor< Cmpt > T() const
Transpose.
Definition: TensorI.H:286
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
virtual void reverseTransformDirection(vector &d, const label faceI) const
Transform a patch-based direction from this side to nbr side.
const List< Face > & localFaces() const
Return patch faces addressing into local point list.
A List obtained as a section of another List.
Definition: SubList.H:53
static const Tensor I
Definition: Tensor.H:84
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
virtual void initUpdateMesh(PstreamBuffers &)
Initialise the update of the patch topology.
Definition: polyPatch.H:115
virtual void movePoints(const Field< PointType > &)
Correct patch after moving points.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
The coupledPolyPatch is an abstract base class for patches that couple regions of the computational d...
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of point.
Definition: meshTools.C:209
UEqn clear()
#define FatalIOErrorIn(functionName, ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:325
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
vector rotationAxis_
Axis of rotation for rotational cyclics.
const Field< PointType > & points() const
Return reference to global points.
const dictionary surfDict_
Dictionary used during projection surface construction.
static const word null
An empty word.
Definition: word.H:77
defineTypeNameAndDebug(combustionModel, 0)
virtual bool order(PstreamBuffers &, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for primitivePatch.
const word & neighbPatchName() const
Neighbour patch name.
dictionary subOrEmptyDict(const word &, const bool mustRead=false) const
Find and return a sub-dictionary as a copy, or.
Definition: dictionary.C:675
prefixOSstream Pout(cout,"Pout")
Definition: IOstreams.H:53
List< bool > boolList
Bool container classes.
Definition: boolList.H:50
virtual void movePoints(PstreamBuffers &, const pointField &p)
Correct patches after moving points.
Definition: polyPatch.C:57
dimensionedScalar sin(const dimensionedScalar &ds)