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