cyclicPolyPatch.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-2018 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 "cyclicPolyPatch.H"
28 #include "polyBoundaryMesh.H"
29 #include "polyMesh.H"
30 #include "demandDrivenData.H"
31 #include "OFstream.H"
32 #include "matchPoints.H"
33 #include "EdgeMap.H"
34 #include "Time.H"
35 #include "transformField.H"
36 #include "SubField.H"
37 #include "unitConversion.H"
38 
39 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40 
41 namespace Foam
42 {
43  defineTypeNameAndDebug(cyclicPolyPatch, 0);
44 
45  addToRunTimeSelectionTable(polyPatch, cyclicPolyPatch, word);
46  addToRunTimeSelectionTable(polyPatch, cyclicPolyPatch, dictionary);
47 }
48 
49 
50 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
51 
52 Foam::label Foam::cyclicPolyPatch::findMaxArea
53 (
54  const pointField& points,
55  const faceList& faces
56 )
57 {
58  label maxI = -1;
59  scalar maxAreaSqr = -great;
60 
61  forAll(faces, facei)
62  {
63  scalar areaSqr = magSqr(faces[facei].area(points));
64 
65  if (areaSqr > maxAreaSqr)
66  {
67  maxAreaSqr = areaSqr;
68  maxI = facei;
69  }
70  }
71  return maxI;
72 }
73 
74 
76 {
77  if (size())
78  {
79  // Half0
80  const cyclicPolyPatch& half0 = *this;
81  vectorField half0Areas(half0.size());
82  forAll(half0, facei)
83  {
84  half0Areas[facei] = half0[facei].area(half0.points());
85  }
86 
87  // Half1
88  const cyclicPolyPatch& half1 = neighbPatch();
89  vectorField half1Areas(half1.size());
90  forAll(half1, facei)
91  {
92  half1Areas[facei] = half1[facei].area(half1.points());
93  }
94 
95  calcTransforms
96  (
97  half0,
98  half0.faceCentres(),
99  half0Areas,
100  half1.faceCentres(),
101  half1Areas
102  );
103  }
104 }
105 
106 
108 (
109  const primitivePatch& half0,
110  const pointField& half0Ctrs,
111  const vectorField& half0Areas,
112  const pointField& half1Ctrs,
113  const vectorField& half1Areas
114 )
115 {
116  if (debug && owner())
117  {
118  fileName casePath(boundaryMesh().mesh().time().path());
119  {
120  fileName nm0(casePath/name()+"_faces.obj");
121  Pout<< "cyclicPolyPatch::calcTransforms : Writing " << name()
122  << " faces to OBJ file " << nm0 << endl;
123  writeOBJ(nm0, half0, half0.points());
124  }
125  const cyclicPolyPatch& half1 = neighbPatch();
126  {
127  fileName nm1(casePath/half1.name()+"_faces.obj");
128  Pout<< "cyclicPolyPatch::calcTransforms : Writing " << half1.name()
129  << " faces to OBJ file " << nm1 << endl;
130  writeOBJ(nm1, half1, half1.points());
131  }
132  {
133  OFstream str(casePath/name()+"_to_" + half1.name() + ".obj");
134  label vertI = 0;
135  Pout<< "cyclicPolyPatch::calcTransforms :"
136  << " Writing coupled face centres as lines to " << str.name()
137  << endl;
138  forAll(half0Ctrs, i)
139  {
140  const point& p0 = half0Ctrs[i];
141  str << "v " << p0.x() << ' ' << p0.y() << ' ' << p0.z() << nl;
142  vertI++;
143  const point& p1 = half1Ctrs[i];
144  str << "v " << p1.x() << ' ' << p1.y() << ' ' << p1.z() << nl;
145  vertI++;
146  str << "l " << vertI-1 << ' ' << vertI << nl;
147  }
148  }
149  }
150 
151 
152  // Some sanity checks
153 
154  if (half0Ctrs.size() != half1Ctrs.size())
155  {
157  << "For patch " << name()
158  << " there are " << half0Ctrs.size()
159  << " face centres, for the neighbour patch " << neighbPatch().name()
160  << " there are " << half1Ctrs.size()
161  << exit(FatalError);
162  }
163 
164  if (transform() != neighbPatch().transform())
165  {
167  << "Patch " << name()
168  << " has transform type " << transformTypeNames[transform()]
169  << ", neighbour patch " << neighbPatchName()
170  << " has transform type "
171  << neighbPatch().transformTypeNames[neighbPatch().transform()]
172  << exit(FatalError);
173  }
174 
175 
176  // Calculate transformation tensors
177 
178  if (half0Ctrs.size() > 0)
179  {
180  vectorField half0Normals(half0Areas.size());
181  vectorField half1Normals(half1Areas.size());
182 
183  scalar maxAreaDiff = -great;
184  label maxAreaFacei = -1;
185 
186  forAll(half0, facei)
187  {
188  scalar magSf = mag(half0Areas[facei]);
189  scalar nbrMagSf = mag(half1Areas[facei]);
190  scalar avSf = (magSf + nbrMagSf)/2.0;
191 
192  if (magSf < rootVSmall && nbrMagSf < rootVSmall)
193  {
194  // Undetermined normal. Use dummy normal to force separation
195  // check. (note use of sqrt(vSmall) since that is how mag
196  // scales)
197  half0Normals[facei] = point(1, 0, 0);
198  half1Normals[facei] = half0Normals[facei];
199  }
200  else
201  {
202  scalar areaDiff = mag(magSf - nbrMagSf)/avSf;
203 
204  if (areaDiff > maxAreaDiff)
205  {
206  maxAreaDiff = areaDiff;
207  maxAreaFacei = facei;
208  }
209 
210  if (areaDiff > matchTolerance())
211  {
213  << "face " << facei
214  << " area does not match neighbour by "
215  << 100*areaDiff
216  << "% -- possible face ordering problem." << endl
217  << "patch:" << name()
218  << " my area:" << magSf
219  << " neighbour area:" << nbrMagSf
220  << " matching tolerance:" << matchTolerance()
221  << endl
222  << "Mesh face:" << start()+facei
223  << " fc:" << half0Ctrs[facei]
224  << endl
225  << "Neighbour fc:" << half1Ctrs[facei]
226  << endl
227  << "If you are certain your matching is correct"
228  << " you can increase the 'matchTolerance' setting"
229  << " in the patch dictionary in the boundary file."
230  << endl
231  << "Rerun with cyclic debug flag set"
232  << " for more information." << exit(FatalError);
233  }
234  else
235  {
236  half0Normals[facei] = half0Areas[facei] / magSf;
237  half1Normals[facei] = half1Areas[facei] / nbrMagSf;
238  }
239  }
240  }
241 
242 
243  // Print area match
244  if (debug)
245  {
246  Pout<< "cyclicPolyPatch::calcTransforms :"
247  << " patch:" << name()
248  << " Max area error:" << 100*maxAreaDiff << "% at face:"
249  << maxAreaFacei << " at:" << half0Ctrs[maxAreaFacei]
250  << " coupled face at:" << half1Ctrs[maxAreaFacei]
251  << endl;
252  }
253 
254 
255  // Calculate transformation tensors
256 
257  if (transform() == ROTATIONAL)
258  {
259  // Calculate using the given rotation axis and centre. Do not
260  // use calculated normals.
261  vector n0 = findFaceMaxRadius(half0Ctrs);
262  vector n1 = -findFaceMaxRadius(half1Ctrs);
263  n0 /= mag(n0) + vSmall;
264  n1 /= mag(n1) + vSmall;
265 
266  if (debug)
267  {
268  scalar theta = radToDeg(acos(n0 & n1));
269 
270  Pout<< "cyclicPolyPatch::calcTransforms :"
271  << " patch:" << name()
272  << " Specified rotation :"
273  << " n0:" << n0 << " n1:" << n1
274  << " swept angle: " << theta << " [deg]"
275  << endl;
276  }
277 
278  // Extended tensor from two local coordinate systems calculated
279  // using normal and rotation axis
280  const tensor E0
281  (
282  rotationAxis_,
283  (n0 ^ rotationAxis_),
284  n0
285  );
286  const tensor E1
287  (
288  rotationAxis_,
289  (-n1 ^ rotationAxis_),
290  -n1
291  );
292  const tensor revT(E1.T() & E0);
293 
294  const_cast<tensorField&>(forwardT()) = tensorField(1, revT.T());
295  const_cast<tensorField&>(reverseT()) = tensorField(1, revT);
296  const_cast<vectorField&>(separation()).setSize(0);
297  const_cast<boolList&>(collocated()) = boolList(1, false);
298  }
299  else
300  {
301  scalarField half0Tols
302  (
303  matchTolerance()
304  *calcFaceTol
305  (
306  half0,
307  half0.points(),
308  static_cast<const pointField&>(half0Ctrs)
309  )
310  );
311 
312  calcTransformTensors
313  (
314  static_cast<const pointField&>(half0Ctrs),
315  static_cast<const pointField&>(half1Ctrs),
316  half0Normals,
317  half1Normals,
318  half0Tols,
319  matchTolerance(),
320  transform()
321  );
322 
323 
324  if (transform() == TRANSLATIONAL)
325  {
326  if (debug)
327  {
328  Pout<< "cyclicPolyPatch::calcTransforms :"
329  << " patch:" << name()
330  << " Specified separation vector : "
331  << separationVector_ << endl;
332  }
333 
334  // Check that separation vectors are same.
335  const scalar avgTol = average(half0Tols);
336  if
337  (
338  mag(separationVector_ + neighbPatch().separationVector_)
339  > avgTol
340  )
341  {
343  << "Specified separation vector " << separationVector_
344  << " differs by that of neighbouring patch "
345  << neighbPatch().separationVector_
346  << " by more than tolerance " << avgTol << endl
347  << "patch:" << name()
348  << " neighbour:" << neighbPatchName()
349  << endl;
350  }
351 
352 
353  // Override computed transform with specified.
354  if
355  (
356  separation().size() != 1
357  || mag(separation()[0] - separationVector_) > avgTol
358  )
359  {
361  << "Specified separationVector " << separationVector_
362  << " differs from computed separation vector "
363  << separation() << endl
364  << "This probably means your geometry is not consistent"
365  << " with the specified separation and might lead"
366  << " to problems." << endl
367  << "Continuing with specified separation vector "
368  << separationVector_ << endl
369  << "patch:" << name()
370  << " neighbour:" << neighbPatchName()
371  << endl;
372  }
373 
374  // Set tensors
375  const_cast<tensorField&>(forwardT()).clear();
376  const_cast<tensorField&>(reverseT()).clear();
377  const_cast<vectorField&>(separation()) = vectorField
378  (
379  1,
380  separationVector_
381  );
382  const_cast<boolList&>(collocated()) = boolList(1, false);
383  }
384  }
385  }
386 }
387 
388 
389 void Foam::cyclicPolyPatch::getCentresAndAnchors
390 (
391  const primitivePatch& pp0,
392  const primitivePatch& pp1,
393 
394  pointField& half0Ctrs,
395  pointField& half1Ctrs,
396  pointField& anchors0,
397  scalarField& tols
398 ) const
399 {
400  // Get geometric data on both halves.
401  half0Ctrs = pp0.faceCentres();
402  anchors0 = getAnchorPoints(pp0, pp0.points(), transform());
403  half1Ctrs = pp1.faceCentres();
404 
405  if (debug)
406  {
407  Pout<< "cyclicPolyPatch::getCentresAndAnchors :"
408  << " patch:" << name() << nl
409  << "half0 untransformed faceCentres (avg) : "
410  << gAverage(half0Ctrs) << nl
411  << "half1 untransformed faceCentres (avg) : "
412  << gAverage(half1Ctrs) << endl;
413  }
414 
415  if (half0Ctrs.size())
416  {
417  switch (transform())
418  {
419  case ROTATIONAL:
420  {
421  vector n0 = findFaceMaxRadius(half0Ctrs);
422  vector n1 = -findFaceMaxRadius(half1Ctrs);
423  n0 /= mag(n0) + vSmall;
424  n1 /= mag(n1) + vSmall;
425 
426  if (debug)
427  {
428  scalar theta = radToDeg(acos(n0 & n1));
429 
430  Pout<< "cyclicPolyPatch::getCentresAndAnchors :"
431  << " patch:" << name()
432  << " Specified rotation :"
433  << " n0:" << n0 << " n1:" << n1
434  << " swept angle: " << theta << " [deg]"
435  << endl;
436  }
437 
438  // Extended tensor from two local coordinate systems calculated
439  // using normal and rotation axis
440  const tensor E0
441  (
442  rotationAxis_,
443  (n0 ^ rotationAxis_),
444  n0
445  );
446  const tensor E1
447  (
448  rotationAxis_,
449  (-n1 ^ rotationAxis_),
450  -n1
451  );
452  const tensor revT(E1.T() & E0);
453 
454  // Rotation
455  forAll(half0Ctrs, facei)
456  {
457  half0Ctrs[facei] =
459  (
460  revT,
461  half0Ctrs[facei] - rotationCentre_
462  )
463  + rotationCentre_;
464  anchors0[facei] =
466  (
467  revT,
468  anchors0[facei] - rotationCentre_
469  )
470  + rotationCentre_;
471  }
472 
473  break;
474  }
475  case TRANSLATIONAL:
476  {
477  // Transform 0 points.
478 
479  if (debug)
480  {
481  Pout<< "cyclicPolyPatch::getCentresAndAnchors :"
482  << " patch:" << name()
483  << "Specified translation : " << separationVector_
484  << endl;
485  }
486 
487  // Note: getCentresAndAnchors gets called on the slave side
488  // so separationVector is owner-slave points.
489 
490  half0Ctrs -= separationVector_;
491  anchors0 -= separationVector_;
492  break;
493  }
494  default:
495  {
496  // Assumes that cyclic is rotational. This is also the initial
497  // condition for patches without faces.
498 
499  // Determine the face with max area on both halves. These
500  // two faces are used to determine the transformation tensors
501  const label max0I = findMaxArea(pp0.points(), pp0);
502  const vector n0 = pp0[max0I].normal(pp0.points());
503 
504  const label max1I = findMaxArea(pp1.points(), pp1);
505  const vector n1 = pp1[max1I].normal(pp1.points());
506 
507  if (mag(n0 & n1) < 1-matchTolerance())
508  {
509  if (debug)
510  {
511  Pout<< "cyclicPolyPatch::getCentresAndAnchors :"
512  << " patch:" << name()
513  << " Detected rotation :"
514  << " n0:" << n0 << " n1:" << n1 << endl;
515  }
516 
517  // Rotation (around origin)
518  const tensor revT(rotationTensor(n0, -n1));
519 
520  // Rotation
521  forAll(half0Ctrs, facei)
522  {
523  half0Ctrs[facei] = Foam::transform
524  (
525  revT,
526  half0Ctrs[facei]
527  );
528  anchors0[facei] = Foam::transform
529  (
530  revT,
531  anchors0[facei]
532  );
533  }
534  }
535  else
536  {
537  // Parallel translation. Get average of all used points.
538 
539  const point ctr0(sum(pp0.localPoints())/pp0.nPoints());
540  const point ctr1(sum(pp1.localPoints())/pp1.nPoints());
541 
542  if (debug)
543  {
544  Pout<< "cyclicPolyPatch::getCentresAndAnchors :"
545  << " patch:" << name()
546  << " Detected translation :"
547  << " n0:" << n0 << " n1:" << n1
548  << " ctr0:" << ctr0 << " ctr1:" << ctr1 << endl;
549  }
550 
551  half0Ctrs += ctr1 - ctr0;
552  anchors0 += ctr1 - ctr0;
553  }
554  break;
555  }
556  }
557  }
558 
559  // Calculate typical distance per face
560  tols = matchTolerance()*calcFaceTol(pp1, pp1.points(), half1Ctrs);
561 }
562 
563 
564 Foam::vector Foam::cyclicPolyPatch::findFaceMaxRadius
565 (
566  const pointField& faceCentres
567 ) const
568 {
569  // Determine a face furthest away from the axis
570 
571  const vectorField n((faceCentres - rotationCentre_) ^ rotationAxis_);
572 
573  const scalarField magRadSqr(magSqr(n));
574 
575  label facei = findMax(magRadSqr);
576 
577  if (debug)
578  {
579  Info<< "findFaceMaxRadius(const pointField&) : patch: " << name() << nl
580  << " rotFace = " << facei << nl
581  << " point = " << faceCentres[facei] << nl
582  << " distance = " << Foam::sqrt(magRadSqr[facei])
583  << endl;
584  }
585 
586  return n[facei];
587 }
588 
589 
590 // * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * * * * //
591 
593 (
594  const word& name,
595  const label size,
596  const label start,
597  const label index,
598  const polyBoundaryMesh& bm,
599  const word& patchType,
600  const transformType transform
601 )
602 :
603  coupledPolyPatch(name, size, start, index, bm, patchType, transform),
604  neighbPatchName_(word::null),
605  neighbPatchID_(-1),
606  rotationAxis_(Zero),
607  rotationCentre_(Zero),
608  separationVector_(Zero),
609  coupledPointsPtr_(nullptr),
610  coupledEdgesPtr_(nullptr)
611 {
612  // Neighbour patch might not be valid yet so no transformation
613  // calculation possible.
614 }
615 
616 
618 (
619  const word& name,
620  const label size,
621  const label start,
622  const label index,
623  const polyBoundaryMesh& bm,
624  const word& neighbPatchName,
625  const transformType transform,
626  const vector& rotationAxis,
627  const point& rotationCentre,
628  const vector& separationVector
629 )
630 :
631  coupledPolyPatch(name, size, start, index, bm, typeName, transform),
632  neighbPatchName_(neighbPatchName),
633  neighbPatchID_(-1),
634  rotationAxis_(rotationAxis),
635  rotationCentre_(rotationCentre),
636  separationVector_(separationVector),
637  coupledPointsPtr_(nullptr),
638  coupledEdgesPtr_(nullptr)
639 {
640  // Neighbour patch might not be valid yet so no transformation
641  // calculation possible.
642 }
643 
644 
646 (
647  const word& name,
648  const dictionary& dict,
649  const label index,
650  const polyBoundaryMesh& bm,
651  const word& patchType
652 )
653 :
654  coupledPolyPatch(name, dict, index, bm, patchType),
655  neighbPatchName_(dict.lookupOrDefault("neighbourPatch", word::null)),
656  coupleGroup_(dict),
657  neighbPatchID_(-1),
658  rotationAxis_(Zero),
659  rotationCentre_(Zero),
660  separationVector_(Zero),
661  coupledPointsPtr_(nullptr),
662  coupledEdgesPtr_(nullptr)
663 {
664  if (neighbPatchName_ == word::null && !coupleGroup_.valid())
665  {
667  (
668  dict
669  ) << "No \"neighbourPatch\" provided." << endl
670  << "Is your mesh uptodate with split cyclics?" << endl
671  << "Run foamUpgradeCyclics to convert mesh and fields"
672  << " to split cyclics." << exit(FatalIOError);
673  }
674 
675  if (neighbPatchName_ == name)
676  {
678  << "Neighbour patch name " << neighbPatchName_
679  << " cannot be the same as this patch " << name
680  << exit(FatalIOError);
681  }
682 
683  switch (transform())
684  {
685  case ROTATIONAL:
686  {
687  dict.lookup("rotationAxis") >> rotationAxis_;
688  dict.lookup("rotationCentre") >> rotationCentre_;
689 
690  scalar magRot = mag(rotationAxis_);
691  if (magRot < small)
692  {
694  << "Illegal rotationAxis " << rotationAxis_ << endl
695  << "Please supply a non-zero vector."
696  << exit(FatalIOError);
697  }
698  rotationAxis_ /= magRot;
699 
700  break;
701  }
702  case TRANSLATIONAL:
703  {
704  dict.lookup("separationVector") >> separationVector_;
705  break;
706  }
707  default:
708  {
709  // no additional info required
710  }
711  }
712 
713  // Neighbour patch might not be valid yet so no transformation
714  // calculation possible.
715 }
716 
717 
719 (
720  const cyclicPolyPatch& pp,
721  const polyBoundaryMesh& bm
722 )
723 :
724  coupledPolyPatch(pp, bm),
725  neighbPatchName_(pp.neighbPatchName_),
726  coupleGroup_(pp.coupleGroup_),
727  neighbPatchID_(-1),
728  rotationAxis_(pp.rotationAxis_),
729  rotationCentre_(pp.rotationCentre_),
730  separationVector_(pp.separationVector_),
731  coupledPointsPtr_(nullptr),
732  coupledEdgesPtr_(nullptr)
733 {
734  // Neighbour patch might not be valid yet so no transformation
735  // calculation possible.
736 }
737 
738 
740 (
741  const cyclicPolyPatch& pp,
742  const polyBoundaryMesh& bm,
743  const label index,
744  const label newSize,
745  const label newStart,
746  const word& neighbName
747 )
748 :
749  coupledPolyPatch(pp, bm, index, newSize, newStart),
750  neighbPatchName_(neighbName),
751  coupleGroup_(pp.coupleGroup_),
752  neighbPatchID_(-1),
753  rotationAxis_(pp.rotationAxis_),
754  rotationCentre_(pp.rotationCentre_),
755  separationVector_(pp.separationVector_),
756  coupledPointsPtr_(nullptr),
757  coupledEdgesPtr_(nullptr)
758 {
759  if (neighbName == name())
760  {
762  << "Neighbour patch name " << neighbName
763  << " cannot be the same as this patch " << name()
764  << exit(FatalError);
765  }
766 
767  // Neighbour patch might not be valid yet so no transformation
768  // calculation possible.
769 }
770 
771 
773 (
774  const cyclicPolyPatch& pp,
775  const polyBoundaryMesh& bm,
776  const label index,
777  const labelUList& mapAddressing,
778  const label newStart
779 )
780 :
781  coupledPolyPatch(pp, bm, index, mapAddressing, newStart),
782  neighbPatchName_(pp.neighbPatchName_),
783  coupleGroup_(pp.coupleGroup_),
784  neighbPatchID_(-1),
785  rotationAxis_(pp.rotationAxis_),
786  rotationCentre_(pp.rotationCentre_),
787  separationVector_(pp.separationVector_),
788  coupledPointsPtr_(nullptr),
789  coupledEdgesPtr_(nullptr)
790 {}
791 
792 
793 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
794 
796 {
797  deleteDemandDrivenData(coupledPointsPtr_);
798  deleteDemandDrivenData(coupledEdgesPtr_);
799 }
800 
801 
802 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
803 
805 {
806  if (neighbPatchName_.empty())
807  {
808  // Try and use patchGroup to find samplePatch and sampleRegion
809  label patchID = coupleGroup_.findOtherPatchID(*this);
810 
811  neighbPatchName_ = boundaryMesh()[patchID].name();
812  }
813  return neighbPatchName_;
814 }
815 
816 
818 {
819  if (neighbPatchID_ == -1)
820  {
821  neighbPatchID_ = this->boundaryMesh().findPatchID(neighbPatchName());
822 
823  if (neighbPatchID_ == -1)
824  {
826  << "Illegal neighbourPatch name " << neighbPatchName()
827  << endl << "Valid patch names are "
828  << this->boundaryMesh().names()
829  << exit(FatalError);
830  }
831 
832  // Check that it is a cyclic
833  const cyclicPolyPatch& nbrPatch = refCast<const cyclicPolyPatch>
834  (
835  this->boundaryMesh()[neighbPatchID_]
836  );
837 
838  if (nbrPatch.neighbPatchName() != name())
839  {
841  << "Patch " << name()
842  << " specifies neighbour patch " << neighbPatchName()
843  << endl << " but that in return specifies "
844  << nbrPatch.neighbPatchName()
845  << endl;
846  }
847  }
848  return neighbPatchID_;
849 }
850 
851 
853 {
854  if (!parallel())
855  {
856  if (transform() == ROTATIONAL)
857  {
858  l =
859  Foam::transform(forwardT(), l-rotationCentre_)
860  + rotationCentre_;
861  }
862  else
863  {
864  l = Foam::transform(forwardT(), l);
865  }
866  }
867  else if (separated())
868  {
869  // transformPosition gets called on the receiving side,
870  // separation gets calculated on the sending side so subtract.
871 
872  const vectorField& s = separation();
873  if (s.size() == 1)
874  {
875  forAll(l, i)
876  {
877  l[i] -= s[0];
878  }
879  }
880  else
881  {
882  l -= s;
883  }
884  }
885 }
886 
887 
889 {
890  if (!parallel())
891  {
892  const tensor& T =
893  (
894  forwardT().size() == 1
895  ? forwardT()[0]
896  : forwardT()[facei]
897  );
898 
899  if (transform() == ROTATIONAL)
900  {
901  l = Foam::transform(T, l-rotationCentre_) + rotationCentre_;
902  }
903  else
904  {
905  l = Foam::transform(T, l);
906  }
907  }
908  else if (separated())
909  {
910  const vector& s =
911  (
912  separation().size() == 1
913  ? separation()[0]
914  : separation()[facei]
915  );
916 
917  l -= s;
918  }
919 }
920 
921 
923 {
925 }
926 
927 
929 (
930  const primitivePatch& referPatch,
931  pointField& nbrCtrs,
932  vectorField& nbrAreas,
933  pointField& nbrCc
934 )
935 {}
936 
937 
939 (
940  const primitivePatch& referPatch,
941  const pointField& thisCtrs,
942  const vectorField& thisAreas,
943  const pointField& thisCc,
944  const pointField& nbrCtrs,
945  const vectorField& nbrAreas,
946  const pointField& nbrCc
947 )
948 {
949  calcTransforms
950  (
951  referPatch,
952  thisCtrs,
953  thisAreas,
954  nbrCtrs,
955  nbrAreas
956  );
957 }
958 
959 
961 {
962  calcGeometry
963  (
964  *this,
965  faceCentres(),
966  faceAreas(),
967  faceCellCentres(),
968  neighbPatch().faceCentres(),
969  neighbPatch().faceAreas(),
970  neighbPatch().faceCellCentres()
971  );
972 }
973 
974 
976 (
977  PstreamBuffers& pBufs,
978  const pointField& p
979 )
980 {
981  polyPatch::initMovePoints(pBufs, p);
982 }
983 
984 
986 (
987  PstreamBuffers& pBufs,
988  const pointField& p
989 )
990 {
991  polyPatch::movePoints(pBufs, p);
992  calcTransforms();
993 }
994 
995 
997 {
999 }
1000 
1001 
1003 {
1004  polyPatch::updateMesh(pBufs);
1005  deleteDemandDrivenData(coupledPointsPtr_);
1006  deleteDemandDrivenData(coupledEdgesPtr_);
1007 }
1008 
1009 
1011 {
1012  if (!coupledPointsPtr_)
1013  {
1014  const faceList& nbrLocalFaces = neighbPatch().localFaces();
1015  const labelList& nbrMeshPoints = neighbPatch().meshPoints();
1016 
1017  // Now all we know is that relative face index in *this is same
1018  // as coupled face in nbrPatch and also that the 0th vertex
1019  // corresponds.
1020 
1021  // From local point to nbrPatch or -1.
1022  labelList coupledPoint(nPoints(), -1);
1023 
1024  forAll(*this, patchFacei)
1025  {
1026  const face& fA = localFaces()[patchFacei];
1027  const face& fB = nbrLocalFaces[patchFacei];
1028 
1029  forAll(fA, indexA)
1030  {
1031  label patchPointA = fA[indexA];
1032 
1033  if (coupledPoint[patchPointA] == -1)
1034  {
1035  label indexB = (fB.size() - indexA) % fB.size();
1036 
1037  // Filter out points on wedge axis
1038  if (meshPoints()[patchPointA] != nbrMeshPoints[fB[indexB]])
1039  {
1040  coupledPoint[patchPointA] = fB[indexB];
1041  }
1042  }
1043  }
1044  }
1045 
1046  coupledPointsPtr_ = new edgeList(nPoints());
1047  edgeList& connected = *coupledPointsPtr_;
1048 
1049  // Extract coupled points.
1050  label connectedI = 0;
1051 
1052  forAll(coupledPoint, i)
1053  {
1054  if (coupledPoint[i] != -1)
1055  {
1056  connected[connectedI++] = edge(i, coupledPoint[i]);
1057  }
1058  }
1059 
1060  connected.setSize(connectedI);
1061 
1062  if (debug)
1063  {
1064  OFstream str
1065  (
1066  boundaryMesh().mesh().time().path()
1067  /name() + "_coupledPoints.obj"
1068  );
1069  label vertI = 0;
1070 
1071  Pout<< "Writing file " << str.name() << " with coordinates of "
1072  << "coupled points" << endl;
1073 
1074  forAll(connected, i)
1075  {
1076  const point& a = points()[meshPoints()[connected[i][0]]];
1077  const point& b = points()[nbrMeshPoints[connected[i][1]]];
1078 
1079  str<< "v " << a.x() << ' ' << a.y() << ' ' << a.z() << nl;
1080  str<< "v " << b.x() << ' ' << b.y() << ' ' << b.z() << nl;
1081  vertI += 2;
1082 
1083  str<< "l " << vertI-1 << ' ' << vertI << nl;
1084  }
1085  }
1086  }
1087  return *coupledPointsPtr_;
1088 }
1089 
1090 
1092 {
1093  if (!coupledEdgesPtr_)
1094  {
1095  const edgeList& pointCouples = coupledPoints();
1096 
1097  // Build map from points on *this (A) to points on neighbourpatch (B)
1098  Map<label> aToB(2*pointCouples.size());
1099 
1100  forAll(pointCouples, i)
1101  {
1102  const edge& e = pointCouples[i];
1103 
1104  aToB.insert(e[0], e[1]);
1105  }
1106 
1107  // Map from edge on A to points (in B indices)
1108  EdgeMap<label> edgeMap(nEdges());
1109 
1110  forAll(*this, patchFacei)
1111  {
1112  const labelList& fEdges = faceEdges()[patchFacei];
1113 
1114  forAll(fEdges, i)
1115  {
1116  label edgeI = fEdges[i];
1117 
1118  const edge& e = edges()[edgeI];
1119 
1120  // Convert edge end points to corresponding points on B side.
1121  Map<label>::const_iterator fnd0 = aToB.find(e[0]);
1122  if (fnd0 != aToB.end())
1123  {
1124  Map<label>::const_iterator fnd1 = aToB.find(e[1]);
1125  if (fnd1 != aToB.end())
1126  {
1127  edgeMap.insert(edge(fnd0(), fnd1()), edgeI);
1128  }
1129  }
1130  }
1131  }
1132 
1133  // Use the edgeMap to get the edges on the B side.
1134 
1135  const cyclicPolyPatch& neighbPatch = this->neighbPatch();
1136  const labelList& nbrMp = neighbPatch.meshPoints();
1137  const labelList& mp = meshPoints();
1138 
1139 
1140 
1141  coupledEdgesPtr_ = new edgeList(edgeMap.size());
1142  edgeList& coupledEdges = *coupledEdgesPtr_;
1143  label coupleI = 0;
1144 
1145  forAll(neighbPatch, patchFacei)
1146  {
1147  const labelList& fEdges = neighbPatch.faceEdges()[patchFacei];
1148 
1149  forAll(fEdges, i)
1150  {
1151  label edgeI = fEdges[i];
1152 
1153  const edge& e = neighbPatch.edges()[edgeI];
1154 
1155  // Look up A edge from HashTable.
1156  EdgeMap<label>::iterator iter = edgeMap.find(e);
1157 
1158  if (iter != edgeMap.end())
1159  {
1160  label edgeA = iter();
1161  const edge& eA = edges()[edgeA];
1162 
1163  // Store correspondence. Filter out edges on wedge axis.
1164  if
1165  (
1166  edge(mp[eA[0]], mp[eA[1]])
1167  != edge(nbrMp[e[0]], nbrMp[e[1]])
1168  )
1169  {
1170  coupledEdges[coupleI++] = edge(edgeA, edgeI);
1171  }
1172 
1173  // Remove so we build unique list
1174  edgeMap.erase(iter);
1175  }
1176  }
1177  }
1178  coupledEdges.setSize(coupleI);
1179 
1180 
1181  // Some checks
1182 
1183  forAll(coupledEdges, i)
1184  {
1185  const edge& e = coupledEdges[i];
1186 
1187  if (e[0] < 0 || e[1] < 0)
1188  {
1190  << "Problem : at position " << i
1191  << " illegal couple:" << e
1192  << abort(FatalError);
1193  }
1194  }
1195 
1196  if (debug)
1197  {
1198  OFstream str
1199  (
1200  boundaryMesh().mesh().time().path()
1201  /name() + "_coupledEdges.obj"
1202  );
1203  label vertI = 0;
1204 
1205  Pout<< "Writing file " << str.name() << " with centres of "
1206  << "coupled edges" << endl;
1207 
1208  forAll(coupledEdges, i)
1209  {
1210  const edge& e = coupledEdges[i];
1211 
1212  const point& a = edges()[e[0]].centre(localPoints());
1213  const point& b = neighbPatch.edges()[e[1]].centre
1214  (
1215  neighbPatch.localPoints()
1216  );
1217 
1218  str<< "v " << a.x() << ' ' << a.y() << ' ' << a.z() << nl;
1219  str<< "v " << b.x() << ' ' << b.y() << ' ' << b.z() << nl;
1220  vertI += 2;
1221 
1222  str<< "l " << vertI-1 << ' ' << vertI << nl;
1223  }
1224  }
1225  }
1226  return *coupledEdgesPtr_;
1227 }
1228 
1229 
1232  PstreamBuffers&,
1233  const primitivePatch& pp
1234 ) const
1235 {
1236  if (owner())
1237  {
1238  // Save patch for use in non-owner side ordering. Equivalent to
1239  // processorPolyPatch using OPstream.
1240  ownerPatchPtr_.reset
1241  (
1242  new primitivePatch
1243  (
1244  pp,
1245  pp.points()
1246  )
1247  );
1248  }
1249 }
1250 
1251 
1254  PstreamBuffers& pBufs,
1255  const primitivePatch& pp,
1256  labelList& faceMap,
1257  labelList& rotation
1258 ) const
1259 {
1260  if (debug)
1261  {
1262  Pout<< "order : of " << pp.size()
1263  << " faces of patch:" << name()
1264  << " neighbour:" << neighbPatchName()
1265  << endl;
1266  }
1267  faceMap.setSize(pp.size());
1268  faceMap = -1;
1269 
1270  rotation.setSize(pp.size());
1271  rotation = 0;
1272 
1273  if (transform() == NOORDERING)
1274  {
1275  // No faces, nothing to change.
1276  return false;
1277  }
1278 
1279  if (owner())
1280  {
1281  // Do nothing (i.e. identical mapping, zero rotation).
1282  // See comment at top.
1283  forAll(faceMap, patchFacei)
1284  {
1285  faceMap[patchFacei] = patchFacei;
1286  }
1287 
1288  return false;
1289  }
1290  else
1291  {
1292  // Get stored geometry from initOrder invocation of owner.
1293  const primitivePatch& pp0 = neighbPatch().ownerPatchPtr_();
1294 
1295  // Get geometric quantities
1296  pointField half0Ctrs, half1Ctrs, anchors0;
1297  scalarField tols;
1298  getCentresAndAnchors
1299  (
1300  pp0,
1301  pp,
1302 
1303  half0Ctrs,
1304  half1Ctrs,
1305  anchors0,
1306  tols
1307  );
1308 
1309  if (debug)
1310  {
1311  Pout<< "half0 transformed faceCentres (avg) : "
1312  << gAverage(half0Ctrs) << nl
1313  << "half1 untransformed faceCentres (avg) : "
1314  << gAverage(half1Ctrs) << endl;
1315  }
1316 
1317  // Geometric match of face centre vectors
1318  bool matchedAll = matchPoints
1319  (
1320  half1Ctrs,
1321  half0Ctrs,
1322  tols,
1323  true,
1324  faceMap
1325  );
1326 
1327  if (!matchedAll || debug)
1328  {
1329  // Dump halves
1330  fileName nm0
1331  (
1332  boundaryMesh().mesh().time().path()
1333  /neighbPatch().name()+"_faces.obj"
1334  );
1335  Pout<< "cyclicPolyPatch::order : Writing neighbour"
1336  << " faces to OBJ file " << nm0 << endl;
1337  writeOBJ(nm0, pp0, pp0.points());
1338 
1339  fileName nm1
1340  (
1341  boundaryMesh().mesh().time().path()
1342  /name()+"_faces.obj"
1343  );
1344  Pout<< "cyclicPolyPatch::order : Writing my"
1345  << " faces to OBJ file " << nm1 << endl;
1346  writeOBJ(nm1, pp, pp.points());
1347 
1348  OFstream ccStr
1349  (
1350  boundaryMesh().mesh().time().path()
1351  /name() + "_faceCentres.obj"
1352  );
1353  Pout<< "cyclicPolyPatch::order : "
1354  << "Dumping currently found cyclic match as lines between"
1355  << " corresponding face centres to file " << ccStr.name()
1356  << endl;
1357 
1358  // Recalculate untransformed face centres
1359  // pointField rawHalf0Ctrs =
1360  // calcFaceCentres(half0Faces, pp.points());
1361  label vertI = 0;
1362 
1363  forAll(half1Ctrs, i)
1364  {
1365  if (faceMap[i] != -1)
1366  {
1367  // Write edge between c1 and c0
1368  const point& c0 = half0Ctrs[faceMap[i]];
1369  const point& c1 = half1Ctrs[i];
1370  writeOBJ(ccStr, c0, c1, vertI);
1371  }
1372  }
1373  }
1374 
1375  if (!matchedAll)
1376  {
1378  << "Patch:" << name() << " : "
1379  << "Cannot match vectors to faces on both sides of patch"
1380  << endl
1381  << " Perhaps your faces do not match?"
1382  << " The obj files written contain the current match." << endl
1383  << " Continuing with incorrect face ordering from now on!"
1384  << endl;
1385 
1386  return false;
1387  }
1388 
1389 
1390  // Set rotation.
1391  forAll(faceMap, oldFacei)
1392  {
1393  // The face f will be at newFacei (after morphing) and we want its
1394  // anchorPoint (= f[0]) to align with the anchorpoint for the
1395  // corresponding face on the other side.
1396 
1397  label newFacei = faceMap[oldFacei];
1398 
1399  const point& wantedAnchor = anchors0[newFacei];
1400 
1401  rotation[newFacei] = getRotation
1402  (
1403  pp.points(),
1404  pp[oldFacei],
1405  wantedAnchor,
1406  tols[oldFacei]
1407  );
1408 
1409  if (rotation[newFacei] == -1)
1410  {
1412  << "in patch " << name()
1413  << " : "
1414  << "Cannot find point on face " << pp[oldFacei]
1415  << " with vertices "
1416  << IndirectList<point>(pp.points(), pp[oldFacei])()
1417  << " that matches point " << wantedAnchor
1418  << " when matching the halves of processor patch " << name()
1419  << "Continuing with incorrect face ordering from now on!"
1420  << endl;
1421 
1422  return false;
1423  }
1424  }
1425 
1426  ownerPatchPtr_.clear();
1427 
1428  // Return false if no change necessary, true otherwise.
1429 
1430  forAll(faceMap, facei)
1431  {
1432  if (faceMap[facei] != facei || rotation[facei] != 0)
1433  {
1434  return true;
1435  }
1436  }
1437 
1438  return false;
1439  }
1440 }
1441 
1442 
1444 {
1446  if (!neighbPatchName_.empty())
1447  {
1448  os.writeKeyword("neighbourPatch") << neighbPatchName_
1449  << token::END_STATEMENT << nl;
1450  }
1451  coupleGroup_.write(os);
1452  switch (transform())
1453  {
1454  case ROTATIONAL:
1455  {
1456  os.writeKeyword("rotationAxis") << rotationAxis_
1457  << token::END_STATEMENT << nl;
1458  os.writeKeyword("rotationCentre") << rotationCentre_
1459  << token::END_STATEMENT << nl;
1460  break;
1461  }
1462  case TRANSLATIONAL:
1463  {
1464  os.writeKeyword("separationVector") << separationVector_
1465  << token::END_STATEMENT << nl;
1466  break;
1467  }
1468  case NOORDERING:
1469  {
1470  break;
1471  }
1472  default:
1473  {
1474  // no additional info to write
1475  }
1476  }
1477 }
1478 
1479 
1480 // ************************************************************************* //
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
virtual const fileName & name() const
Return the name of the stream.
Definition: OSstream.H:88
label nPoints() const
Return number of points supporting patch faces.
static bool valid(char)
Is this character valid for a word.
Definition: wordI.H:115
dimensionedScalar acos(const dimensionedScalar &ds)
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
scalar radToDeg(const scalar rad)
Conversion from radians to degrees.
A class for handling file names.
Definition: fileName.H:69
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const edgeList & coupledPoints() const
Return connected points (from patch local to neighbour patch local)
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:110
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
virtual void initUpdateMesh(PstreamBuffers &)
Initialise the update of the patch topology.
Definition: polyPatch.H:115
tUEqn clear()
Unit conversion functions.
Output to file stream.
Definition: OFstream.H:82
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
Tensor< Cmpt > T() const
Return transpose.
Definition: TensorI.H:321
List< face > faceList
Definition: faceListFwd.H:43
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface...
Definition: boundaryMesh.H:59
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
virtual void initMovePoints(PstreamBuffers &, const pointField &)
Initialise the patches for moving points.
cyclicPolyPatch(const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm, const word &patchType, const transformType transform=UNKNOWN)
Construct from components.
const Field< PointType > & faceCentres() const
Return face centres for patch.
The coupledPolyPatch is an abstract base class for patches that couple regions of the computational d...
const Cmpt & z() const
Definition: VectorI.H:87
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of point.
Definition: meshTools.C:203
tensor rotationTensor(const vector &n1, const vector &n2)
Rotational transformation tensor from vector n1 to n2.
Definition: transform.H:47
virtual void initMovePoints(PstreamBuffers &, const pointField &)
Initialise the patches for moving points.
Definition: polyPatch.H:108
const edgeList & coupledEdges() const
Return connected edges (from patch local to neighbour patch local).
const fileName & name() const
Return the name of the stream.
Definition: OFstream.H:120
virtual bool order(PstreamBuffers &, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for primitivePatch.
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
Determine correspondence between points. See below.
Macros for easy insertion into run-time selection tables.
const labelList & meshPoints() const
Return labelList of mesh points in patch. They are constructed.
const Cmpt & y() const
Definition: VectorI.H:81
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
bool erase(const iterator &)
Erase a hashedEntry specified by given iterator.
Definition: HashTable.C:371
virtual void movePoints(PstreamBuffers &, const pointField &p)
Correct patches after moving points.
Definition: polyPatch.C:57
#define SeriousErrorInFunction
Report an error message using Foam::SeriousError.
List< bool > boolList
Bool container classes.
Definition: boolList.H:50
bool insert(const Key &, const T &newElmt)
Insert a new hashedEntry.
Definition: HashTableI.H:80
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
Definition: HashTable.C:142
virtual void initUpdateMesh(PstreamBuffers &)
Initialise the update of the patch topology.
const Vector< Cmpt > & centre(const Foam::List< Vector< Cmpt >> &) const
Return *this (used for point which is a typedef to Vector<scalar>.
Definition: VectorI.H:116
points setSize(newPointi)
Spatial transformation functions for primitive fields.
A list of faces which address into the list of points.
bool matchPoints(const UList< point > &pts0, const UList< point > &pts1, const UList< scalar > &matchDistance, const bool verbose, labelList &from0To1, const point &origin=point::zero)
Determine correspondence between pointFields. Gets passed.
Definition: matchPoints.C:33
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))
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
dynamicFvMesh & mesh
const pointField & points
List< edge > edgeList
Definition: edgeList.H:38
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:27
A class for handling words, derived from string.
Definition: word.H:59
label nPoints
Cyclic plane patch.
word name() const
Return file name (part beyond last /)
Definition: fileName.C:179
const Field< PointType > & points() const
Return reference to global points.
static const word null
An empty word.
Definition: word.H:77
virtual void movePoints(PstreamBuffers &, const pointField &)
Correct patches after moving points.
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
static const zero Zero
Definition: zero.H:97
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
Definition: polyPatch.C:62
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)
errorManip< error > abort(error &err)
Definition: errorManip.H:131
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:61
const Cmpt & x() const
Definition: VectorI.H:75
Foam::polyBoundaryMesh.
label findMax(const ListType &, const label start=0)
Find index of max element (and larger than given element).
dimensioned< scalar > magSqr(const dimensioned< Type > &)
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
const Field< PointType > & localPoints() const
Return pointField of points in patch.
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
static const char nl
Definition: Ostream.H:265
fileName path(UMean.rootPath()/UMean.caseName()/"graphs"/UMean.instance())
defineTypeNameAndDebug(combustionModel, 0)
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:54
Buffers for inter-processor communications streams (UOPstream, UIPstream).
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
Definition: polyPatch.H:100
Field< tensor > tensorField
Specialisation of Field<T> for tensor.
void setSize(const label)
Reset size of List.
Definition: List.C:281
Template functions to aid in the implementation of demand driven data.
virtual void transformPosition(pointField &l) const
Transform a patch-based position from other side to this side.
virtual ~cyclicPolyPatch()
Destructor.
vector point
Point is a vector.
Definition: point.H:41
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
#define WarningInFunction
Report a warning using Foam::Warning.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
Type gAverage(const FieldField< Field, Type > &f)
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
virtual label neighbPatchID() const
Neighbour patchID.
const word & neighbPatchName() const
Neighbour patch name.
const dimensionedScalar c1
First radiation constant: default SI units: [W/m2].
messageStream Info
const labelListList & faceEdges() const
Return face-edge addressing.
dimensioned< scalar > mag(const dimensioned< Type > &)
label n
virtual void initOrder(PstreamBuffers &, const primitivePatch &) const
Initialize ordering for primitivePatch. Does not.
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:98
void deleteDemandDrivenData(DataPtr &dataPtr)
A List with indirect addressing.
Definition: IndirectList.H:102
const dimensionedScalar mp
Proton mass.
virtual void calcTransforms()
Recalculate the transformation tensors.
Namespace for OpenFOAM.
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:477
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:576
IOerror FatalIOError
const vectorField::subField faceCentres() const
Return face centres.
Definition: polyPatch.C:284