PointEdgeWave.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-2020 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 "PointEdgeWave.H"
27 #include "polyMesh.H"
28 #include "processorPolyPatch.H"
29 #include "cyclicPolyPatch.H"
30 #include "OPstream.H"
31 #include "IPstream.H"
33 #include "debug.H"
34 #include "typeInfo.H"
35 #include "globalMeshData.H"
36 #include "pointFields.H"
37 
38 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
39 
40 template<class Type, class TrackingData>
42 
43 template<class Type, class TrackingData>
45 
46 namespace Foam
47 {
48  //- Reduction class. If x and y are not equal assign value.
49  template<class Type, class TrackingData>
51  {
52  TrackingData& td_;
53 
54  public:
55  combineEqOp(TrackingData& td)
56  :
57  td_(td)
58  {}
59 
60  void operator()(Type& x, const Type& y) const
61  {
62  if (!x.valid(td_) && y.valid(td_))
63  {
64  x = y;
65  }
66  }
67  };
68 }
69 
70 
71 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
72 
73 // Handle leaving domain. Implementation referred to Type
74 template<class Type, class TrackingData>
76 (
77  const polyPatch& patch,
78  const labelList& patchPointLabels,
79  List<Type>& pointInfo
80 ) const
81 {
82  const labelList& meshPoints = patch.meshPoints();
83 
84  forAll(patchPointLabels, i)
85  {
86  label patchPointi = patchPointLabels[i];
87 
88  const point& pt = patch.points()[meshPoints[patchPointi]];
89 
90  pointInfo[i].leaveDomain(patch, patchPointi, pt, td_);
91  }
92 }
93 
94 
95 // Handle entering domain. Implementation referred to Type
96 template<class Type, class TrackingData>
98 (
99  const polyPatch& patch,
100  const labelList& patchPointLabels,
101  List<Type>& pointInfo
102 ) const
103 {
104  const labelList& meshPoints = patch.meshPoints();
105 
106  forAll(patchPointLabels, i)
107  {
108  label patchPointi = patchPointLabels[i];
109 
110  const point& pt = patch.points()[meshPoints[patchPointi]];
111 
112  pointInfo[i].enterDomain(patch, patchPointi, pt, td_);
113  }
114 }
115 
116 
117 // Transform. Implementation referred to Type
118 template<class Type, class TrackingData>
120 (
121  const polyPatch& patch,
122  const tensor& rotTensor,
123  List<Type>& pointInfo
124 ) const
125 {
126  forAll(pointInfo, i)
127  {
128  pointInfo[i].transform(rotTensor, td_);
129  }
130 }
131 
132 
133 // Update info for pointi, at position pt, with information from
134 // neighbouring edge.
135 // Updates:
136 // - changedPoint_, changedPoints_, nChangedPoints_,
137 // - statistics: nEvals_, nUnvisitedPoints_
138 template<class Type, class TrackingData>
140 (
141  const label pointi,
142  const label neighbourEdgeI,
143  const Type& neighbourInfo,
144  Type& pointInfo
145 )
146 {
147  nEvals_++;
148 
149  bool wasValid = pointInfo.valid(td_);
150 
151  bool propagate =
152  pointInfo.updatePoint
153  (
154  mesh_,
155  pointi,
156  neighbourEdgeI,
157  neighbourInfo,
158  propagationTol_,
159  td_
160  );
161 
162  if (propagate)
163  {
164  if (!changedPoint_[pointi])
165  {
166  changedPoint_[pointi] = true;
167  changedPoints_[nChangedPoints_++] = pointi;
168  }
169  }
170 
171  if (!wasValid && pointInfo.valid(td_))
172  {
173  --nUnvisitedPoints_;
174  }
175 
176  return propagate;
177 }
178 
179 
180 // Update info for pointi, at position pt, with information from
181 // same point.
182 // Updates:
183 // - changedPoint_, changedPoints_, nChangedPoints_,
184 // - statistics: nEvals_, nUnvisitedPoints_
185 template<class Type, class TrackingData>
187 (
188  const label pointi,
189  const Type& neighbourInfo,
190  Type& pointInfo
191 )
192 {
193  nEvals_++;
194 
195  bool wasValid = pointInfo.valid(td_);
196 
197  bool propagate =
198  pointInfo.updatePoint
199  (
200  mesh_,
201  pointi,
202  neighbourInfo,
203  propagationTol_,
204  td_
205  );
206 
207  if (propagate)
208  {
209  if (!changedPoint_[pointi])
210  {
211  changedPoint_[pointi] = true;
212  changedPoints_[nChangedPoints_++] = pointi;
213  }
214  }
215 
216  if (!wasValid && pointInfo.valid(td_))
217  {
218  --nUnvisitedPoints_;
219  }
220 
221  return propagate;
222 }
223 
224 
225 // Update info for edgeI, at position pt, with information from
226 // neighbouring point.
227 // Updates:
228 // - changedEdge_, changedEdges_, nChangedEdges_,
229 // - statistics: nEvals_, nUnvisitedEdge_
230 template<class Type, class TrackingData>
232 (
233  const label edgeI,
234  const label neighbourPointi,
235  const Type& neighbourInfo,
236  Type& edgeInfo
237 )
238 {
239  nEvals_++;
240 
241  bool wasValid = edgeInfo.valid(td_);
242 
243  bool propagate =
244  edgeInfo.updateEdge
245  (
246  mesh_,
247  edgeI,
248  neighbourPointi,
249  neighbourInfo,
250  propagationTol_,
251  td_
252  );
253 
254  if (propagate)
255  {
256  if (!changedEdge_[edgeI])
257  {
258  changedEdge_[edgeI] = true;
259  changedEdges_[nChangedEdges_++] = edgeI;
260  }
261  }
262 
263  if (!wasValid && edgeInfo.valid(td_))
264  {
265  --nUnvisitedEdges_;
266  }
267 
268  return propagate;
269 }
270 
271 
272 // Check if patches of given type name are present
273 template<class Type, class TrackingData>
274 template<class PatchType>
276 {
277  label nPatches = 0;
278 
279  forAll(mesh_.boundaryMesh(), patchi)
280  {
281  if (isA<PatchType>(mesh_.boundaryMesh()[patchi]))
282  {
283  nPatches++;
284  }
285  }
286  return nPatches;
287 }
288 
289 
290 // Transfer all the information to/from neighbouring processors
291 template<class Type, class TrackingData>
293 {
294  // 1. Send all point info on processor patches.
295 
297 
298  DynamicList<Type> patchInfo;
299  DynamicList<label> thisPoints;
300  DynamicList<label> nbrPoints;
301 
302  forAll(mesh_.globalData().processorPatches(), i)
303  {
304  label patchi = mesh_.globalData().processorPatches()[i];
305  const processorPolyPatch& procPatch =
306  refCast<const processorPolyPatch>(mesh_.boundaryMesh()[patchi]);
307 
308  patchInfo.clear();
309  patchInfo.reserve(procPatch.nPoints());
310  thisPoints.clear();
311  thisPoints.reserve(procPatch.nPoints());
312  nbrPoints.clear();
313  nbrPoints.reserve(procPatch.nPoints());
314 
315  // Get all changed points in reverse order
316  const labelList& neighbPoints = procPatch.nbrPoints();
317  forAll(neighbPoints, thisPointi)
318  {
319  label meshPointi = procPatch.meshPoints()[thisPointi];
320  if (changedPoint_[meshPointi])
321  {
322  patchInfo.append(allPointInfo_[meshPointi]);
323  thisPoints.append(thisPointi);
324  nbrPoints.append(neighbPoints[thisPointi]);
325  }
326  }
327 
328  // Adapt for leaving domain
329  leaveDomain(procPatch, thisPoints, patchInfo);
330 
331  // if (debug)
332  //{
333  // Pout<< "Processor patch " << patchi << ' ' << procPatch.name()
334  // << " communicating with " << procPatch.neighbProcNo()
335  // << " Sending:" << patchInfo.size() << endl;
336  //}
337 
338  UOPstream toNeighbour(procPatch.neighbProcNo(), pBufs);
339  toNeighbour << nbrPoints << patchInfo;
340  }
341 
342 
343  pBufs.finishedSends();
344 
345  //
346  // 2. Receive all point info on processor patches.
347  //
348 
349  forAll(mesh_.globalData().processorPatches(), i)
350  {
351  label patchi = mesh_.globalData().processorPatches()[i];
352  const processorPolyPatch& procPatch =
353  refCast<const processorPolyPatch>(mesh_.boundaryMesh()[patchi]);
354 
355  List<Type> patchInfo;
356  labelList patchPoints;
357 
358  {
359  UIPstream fromNeighbour(procPatch.neighbProcNo(), pBufs);
360  fromNeighbour >> patchPoints >> patchInfo;
361  }
362 
363  // if (debug)
364  //{
365  // Pout<< "Processor patch " << patchi << ' ' << procPatch.name()
366  // << " communicating with " << procPatch.neighbProcNo()
367  // << " Received:" << patchInfo.size() << endl;
368  //}
369 
370  // Apply transform to received data for non-parallel planes
371  if (procPatch.transform().transforms())
372  {
373  transform(procPatch, procPatch.transform().T(), patchInfo);
374  }
375 
376  // Adapt for entering domain
377  enterDomain(procPatch, patchPoints, patchInfo);
378 
379  // Merge received info
380  const labelList& meshPoints = procPatch.meshPoints();
381  forAll(patchInfo, i)
382  {
383  label meshPointi = meshPoints[patchPoints[i]];
384 
385  if (!allPointInfo_[meshPointi].equal(patchInfo[i], td_))
386  {
387  updatePoint
388  (
389  meshPointi,
390  patchInfo[i],
391  allPointInfo_[meshPointi]
392  );
393  }
394  }
395  }
396 
397  // Collocated points should be handled by face based transfer
398  // (since that is how connectivity is worked out)
399  // They are also explicitly equalised in handleCollocatedPoints to
400  // guarantee identical values.
401 }
402 
403 
404 template<class Type, class TrackingData>
406 {
407  // 1. Send all point info on cyclic patches.
408 
409  DynamicList<Type> nbrInfo;
410  DynamicList<label> nbrPoints;
411  DynamicList<label> thisPoints;
412 
413  forAll(mesh_.boundaryMesh(), patchi)
414  {
415  const polyPatch& patch = mesh_.boundaryMesh()[patchi];
416 
417  if (isA<cyclicPolyPatch>(patch))
418  {
419  const cyclicPolyPatch& cycPatch =
420  refCast<const cyclicPolyPatch>(patch);
421 
422  nbrInfo.clear();
423  nbrInfo.reserve(cycPatch.nPoints());
424  nbrPoints.clear();
425  nbrPoints.reserve(cycPatch.nPoints());
426  thisPoints.clear();
427  thisPoints.reserve(cycPatch.nPoints());
428 
429  // Collect nbrPatch points that have changed
430  {
431  const cyclicPolyPatch& nbrPatch = cycPatch.nbrPatch();
432  const edgeList& pairs = cycPatch.coupledPoints();
433  const labelList& meshPoints = nbrPatch.meshPoints();
434 
435  forAll(pairs, pairI)
436  {
437  label thisPointi = pairs[pairI][0];
438  label nbrPointi = pairs[pairI][1];
439  label meshPointi = meshPoints[nbrPointi];
440 
441  if (changedPoint_[meshPointi])
442  {
443  nbrInfo.append(allPointInfo_[meshPointi]);
444  nbrPoints.append(nbrPointi);
445  thisPoints.append(thisPointi);
446  }
447  }
448 
449  // nbr : adapt for leaving domain
450  leaveDomain(nbrPatch, nbrPoints, nbrInfo);
451  }
452 
453 
454  // Apply rotation for non-parallel planes
455 
456  if (cycPatch.transform().transforms())
457  {
458  // received data from half1
459  transform(cycPatch, cycPatch.transform().T(), nbrInfo);
460  }
461 
462  // if (debug)
463  //{
464  // Pout<< "Cyclic patch " << patchi << ' ' << patch.name()
465  // << " Changed : " << nbrInfo.size()
466  // << endl;
467  //}
468 
469  // Adapt for entering domain
470  enterDomain(cycPatch, thisPoints, nbrInfo);
471 
472  // Merge received info
473  const labelList& meshPoints = cycPatch.meshPoints();
474  forAll(nbrInfo, i)
475  {
476  label meshPointi = meshPoints[thisPoints[i]];
477 
478  if (!allPointInfo_[meshPointi].equal(nbrInfo[i], td_))
479  {
480  updatePoint
481  (
482  meshPointi,
483  nbrInfo[i],
484  allPointInfo_[meshPointi]
485  );
486  }
487  }
488  }
489  }
490 }
491 
492 
493 // Guarantee collocated points have same information.
494 // Return number of points changed.
495 template<class Type, class TrackingData>
497 {
498  // Transfer onto coupled patch
499  const globalMeshData& gmd = mesh_.globalData();
500  const indirectPrimitivePatch& cpp = gmd.coupledPatch();
501  const labelList& meshPoints = cpp.meshPoints();
502 
503  const mapDistribute& slavesMap = gmd.globalPointSlavesMap();
504  const labelListList& slaves = gmd.globalPointSlaves();
505 
506  List<Type> elems(slavesMap.constructSize());
507  forAll(meshPoints, pointi)
508  {
509  elems[pointi] = allPointInfo_[meshPoints[pointi]];
510  }
511 
512  // Pull slave data onto master (which might or might not have any
513  // initialised points). No need to update transformed slots.
514  slavesMap.distribute(elems, false);
515 
516  // Combine master data with slave data
518 
519  forAll(slaves, pointi)
520  {
521  Type& elem = elems[pointi];
522 
523  const labelList& slavePoints = slaves[pointi];
524 
525  // Combine master with untransformed slave data
526  forAll(slavePoints, j)
527  {
528  cop(elem, elems[slavePoints[j]]);
529  }
530 
531  // Copy result back to slave slots
532  forAll(slavePoints, j)
533  {
534  elems[slavePoints[j]] = elem;
535  }
536  }
537 
538  // Push slave-slot data back to slaves
539  slavesMap.reverseDistribute(elems.size(), elems, false);
540 
541  // Extract back onto mesh
542  forAll(meshPoints, pointi)
543  {
544  if (elems[pointi].valid(td_))
545  {
546  label meshPointi = meshPoints[pointi];
547 
548  Type& elem = allPointInfo_[meshPointi];
549 
550  bool wasValid = elem.valid(td_);
551 
552  // Like updatePoint but bypass Type::updatePoint with its tolerance
553  // checking
554  // if (!elem.valid(td_) || !elem.equal(elems[pointi], td_))
555  if (!elem.equal(elems[pointi], td_))
556  {
557  nEvals_++;
558  elem = elems[pointi];
559 
560  // See if element now valid
561  if (!wasValid && elem.valid(td_))
562  {
563  --nUnvisitedPoints_;
564  }
565 
566  // Update database of changed points
567  if (!changedPoint_[meshPointi])
568  {
569  changedPoint_[meshPointi] = true;
570  changedPoints_[nChangedPoints_++] = meshPointi;
571  }
572  }
573  }
574  }
575 
576  // Sum nChangedPoints over all procs
577  label totNChanged = nChangedPoints_;
578 
579  reduce(totNChanged, sumOp<label>());
580 
581  return totNChanged;
582 }
583 
584 
585 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
586 
587 // Iterate, propagating changedPointsInfo across mesh, until no change (or
588 // maxIter reached). Initial point values specified.
589 template<class Type, class TrackingData>
591 (
592  const polyMesh& mesh,
593  const labelList& changedPoints,
594  const List<Type>& changedPointsInfo,
595 
596  UList<Type>& allPointInfo,
597  UList<Type>& allEdgeInfo,
598 
599  const label maxIter,
600  TrackingData& td
601 )
602 :
603  mesh_(mesh),
604  allPointInfo_(allPointInfo),
605  allEdgeInfo_(allEdgeInfo),
606  td_(td),
607  changedPoint_(mesh_.nPoints(), false),
608  changedPoints_(mesh_.nPoints()),
609  nChangedPoints_(0),
610  changedEdge_(mesh_.nEdges(), false),
611  changedEdges_(mesh_.nEdges()),
612  nChangedEdges_(0),
613  nCyclicPatches_(countPatchType<cyclicPolyPatch>()),
614  nEvals_(0),
615  nUnvisitedPoints_(mesh_.nPoints()),
616  nUnvisitedEdges_(mesh_.nEdges())
617 {
618  if (allPointInfo_.size() != mesh_.nPoints())
619  {
621  << "size of pointInfo work array is not equal to the number"
622  << " of points in the mesh" << endl
623  << " pointInfo :" << allPointInfo_.size() << endl
624  << " mesh.nPoints:" << mesh_.nPoints()
625  << exit(FatalError);
626  }
627  if (allEdgeInfo_.size() != mesh_.nEdges())
628  {
630  << "size of edgeInfo work array is not equal to the number"
631  << " of edges in the mesh" << endl
632  << " edgeInfo :" << allEdgeInfo_.size() << endl
633  << " mesh.nEdges:" << mesh_.nEdges()
634  << exit(FatalError);
635  }
636 
637 
638  // Set from initial changed points data
639  setPointInfo(changedPoints, changedPointsInfo);
640 
641  if (debug)
642  {
643  Info<< typeName << ": Seed points : "
644  << returnReduce(nChangedPoints_, sumOp<label>()) << endl;
645  }
646 
647  // Iterate until nothing changes
648  label iter = iterate(maxIter);
649 
650  if ((maxIter > 0) && (iter >= maxIter))
651  {
653  << "Maximum number of iterations reached. Increase maxIter." << endl
654  << " maxIter:" << maxIter << endl
655  << " nChangedPoints:" << nChangedPoints_ << endl
656  << " nChangedEdges:" << nChangedEdges_ << endl
657  << exit(FatalError);
658  }
659 }
660 
661 
662 template<class Type, class TrackingData>
664 (
665  const polyMesh& mesh,
666  UList<Type>& allPointInfo,
667  UList<Type>& allEdgeInfo,
668  TrackingData& td
669 )
670 :
671  mesh_(mesh),
672  allPointInfo_(allPointInfo),
673  allEdgeInfo_(allEdgeInfo),
674  td_(td),
675  changedPoint_(mesh_.nPoints(), false),
676  changedPoints_(mesh_.nPoints()),
677  nChangedPoints_(0),
678  changedEdge_(mesh_.nEdges(), false),
679  changedEdges_(mesh_.nEdges()),
680  nChangedEdges_(0),
681  nCyclicPatches_(countPatchType<cyclicPolyPatch>()),
682  nEvals_(0),
683  nUnvisitedPoints_(mesh_.nPoints()),
684  nUnvisitedEdges_(mesh_.nEdges())
685 {}
686 
687 
688 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
689 
690 template<class Type, class TrackingData>
692 {}
693 
694 
695 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
696 
697 
698 template<class Type, class TrackingData>
700 {
701  return nUnvisitedPoints_;
702 }
703 
704 
705 template<class Type, class TrackingData>
707 {
708  return nUnvisitedEdges_;
709 }
710 
711 
712 // Copy point information into member data
713 template<class Type, class TrackingData>
715 (
716  const labelList& changedPoints,
717  const List<Type>& changedPointsInfo
718 )
719 {
720  forAll(changedPoints, changedPointi)
721  {
722  label pointi = changedPoints[changedPointi];
723 
724  bool wasValid = allPointInfo_[pointi].valid(td_);
725 
726  // Copy info for pointi
727  allPointInfo_[pointi] = changedPointsInfo[changedPointi];
728 
729  // Maintain count of unset points
730  if (!wasValid && allPointInfo_[pointi].valid(td_))
731  {
732  --nUnvisitedPoints_;
733  }
734 
735  // Mark pointi as changed, both on list and on point itself.
736 
737  if (!changedPoint_[pointi])
738  {
739  changedPoint_[pointi] = true;
740  changedPoints_[nChangedPoints_++] = pointi;
741  }
742  }
743 
744  // Sync
745  handleCollocatedPoints();
746 }
747 
748 
749 // Propagate information from edge to point. Return number of points changed.
750 template<class Type, class TrackingData>
752 {
753  for
754  (
755  label changedEdgeI = 0;
756  changedEdgeI < nChangedEdges_;
757  changedEdgeI++
758  )
759  {
760  label edgeI = changedEdges_[changedEdgeI];
761 
762  if (!changedEdge_[edgeI])
763  {
765  << "edge " << edgeI
766  << " not marked as having been changed" << nl
767  << "This might be caused by multiple occurrences of the same"
768  << " seed point." << abort(FatalError);
769  }
770 
771 
772  const Type& neighbourWallInfo = allEdgeInfo_[edgeI];
773 
774  // Evaluate all connected points (= edge endpoints)
775  const edge& e = mesh_.edges()[edgeI];
776 
777  forAll(e, eI)
778  {
779  Type& currentWallInfo = allPointInfo_[e[eI]];
780 
781  if (!currentWallInfo.equal(neighbourWallInfo, td_))
782  {
783  updatePoint
784  (
785  e[eI],
786  edgeI,
787  neighbourWallInfo,
788  currentWallInfo
789  );
790  }
791  }
792 
793  // Reset status of edge
794  changedEdge_[edgeI] = false;
795  }
796 
797  // Handled all changed edges by now
798  nChangedEdges_ = 0;
799 
800  if (nCyclicPatches_ > 0)
801  {
802  // Transfer changed points across cyclic halves
803  handleCyclicPatches();
804  }
805  if (Pstream::parRun())
806  {
807  // Transfer changed points from neighbouring processors.
808  handleProcPatches();
809  }
810 
811  // if (debug)
812  //{
813  // Pout<< "Changed points : " << nChangedPoints_ << endl;
814  //}
815 
816  // Sum nChangedPoints over all procs
817  label totNChanged = nChangedPoints_;
818 
819  reduce(totNChanged, sumOp<label>());
820 
821  return totNChanged;
822 }
823 
824 
825 // Propagate information from point to edge. Return number of edges changed.
826 template<class Type, class TrackingData>
828 {
829  const labelListList& pointEdges = mesh_.pointEdges();
830 
831  for
832  (
833  label changedPointi = 0;
834  changedPointi < nChangedPoints_;
835  changedPointi++
836  )
837  {
838  label pointi = changedPoints_[changedPointi];
839 
840  if (!changedPoint_[pointi])
841  {
843  << "Point " << pointi
844  << " not marked as having been changed" << nl
845  << "This might be caused by multiple occurrences of the same"
846  << " seed point." << abort(FatalError);
847  }
848 
849  const Type& neighbourWallInfo = allPointInfo_[pointi];
850 
851  // Evaluate all connected edges
852 
853  const labelList& edgeLabels = pointEdges[pointi];
854  forAll(edgeLabels, edgeLabelI)
855  {
856  label edgeI = edgeLabels[edgeLabelI];
857 
858  Type& currentWallInfo = allEdgeInfo_[edgeI];
859 
860  if (!currentWallInfo.equal(neighbourWallInfo, td_))
861  {
862  updateEdge
863  (
864  edgeI,
865  pointi,
866  neighbourWallInfo,
867  currentWallInfo
868  );
869  }
870  }
871 
872  // Reset status of point
873  changedPoint_[pointi] = false;
874  }
875 
876  // Handled all changed points by now
877  nChangedPoints_ = 0;
878 
879  // if (debug)
880  //{
881  // Pout<< "Changed edges : " << nChangedEdges_ << endl;
882  //}
883 
884  // Sum nChangedPoints over all procs
885  label totNChanged = nChangedEdges_;
886 
887  reduce(totNChanged, sumOp<label>());
888 
889  return totNChanged;
890 }
891 
892 
893 // Iterate
894 template<class Type, class TrackingData>
896 (
897  const label maxIter
898 )
899 {
900  if (nCyclicPatches_ > 0)
901  {
902  // Transfer changed points across cyclic halves
903  handleCyclicPatches();
904  }
905  if (Pstream::parRun())
906  {
907  // Transfer changed points from neighbouring processors.
908  handleProcPatches();
909  }
910 
911  nEvals_ = 0;
912 
913  label iter = 0;
914 
915  while (iter < maxIter)
916  {
917  while (iter < maxIter)
918  {
919  if (debug)
920  {
921  Info<< typeName << ": Iteration " << iter << endl;
922  }
923 
924  label nEdges = pointToEdge();
925 
926  if (debug)
927  {
928  Info<< typeName << ": Total changed edges : "
929  << nEdges << endl;
930  }
931 
932  if (nEdges == 0)
933  {
934  break;
935  }
936 
937  label nPoints = edgeToPoint();
938 
939  if (debug)
940  {
941  Info<< typeName << ": Total changed points : "
942  << nPoints << nl
943  << typeName << ": Total evaluations : "
944  << returnReduce(nEvals_, sumOp<label>()) << nl
945  << typeName << ": Remaining unvisited points: "
946  << returnReduce(nUnvisitedPoints_, sumOp<label>()) << nl
947  << typeName << ": Remaining unvisited edges : "
948  << returnReduce(nUnvisitedEdges_, sumOp<label>()) << nl
949  << endl;
950  }
951 
952  if (nPoints == 0)
953  {
954  break;
955  }
956 
957  iter++;
958  }
959 
960 
961  // Enforce collocated points are exactly equal. This might still mean
962  // non-collocated points are not equal though. WIP.
963  label nPoints = handleCollocatedPoints();
964  if (debug)
965  {
966  Info<< typeName << ": Collocated point sync : "
967  << nPoints << nl << endl;
968  }
969 
970  if (nPoints == 0)
971  {
972  break;
973  }
974  }
975 
976  return iter;
977 }
978 
979 
980 // ************************************************************************* //
label nPatches
Definition: readKivaGrid.H:402
label nPoints() const
Return number of points supporting patch faces.
const labelListList & globalPointSlaves() const
Reduction class. If x and y are not equal assign value.
Definition: PointEdgeWave.C:50
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
const labelList & nbrPoints() const
Return neighbour point labels. WIP.
void finishedSends(const bool block=true)
Mark all sends as having been done. This will start receives.
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)
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
void distribute(List< T > &fld, const bool dummyTransform=true, const int tag=UPstream::msgType()) const
Distribute data using default commsType.
bool transforms() const
Return true if the transformer transforms a type.
Definition: transformerI.H:133
label pointToEdge()
Propagate from point to edge. Returns total number of edges.
Various mesh related information for a parallel run. Upon construction, constructs all info using par...
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:291
label iterate(const label maxIter)
Iterate until no changes or maxIter reached. Returns actual.
PointEdgeWave(const polyMesh &mesh, const labelList &initialPoints, const List< Type > &initialPointsInfo, UList< Type > &allPointInfo, UList< Type > &allEdgeInfo, const label maxIter, TrackingData &td=dummyTrackData_)
Construct from mesh, list of changed points with the Type.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const mapDistribute & globalPointSlavesMap() const
label edgeToPoint()
Propagate from edge to point. Returns total number of points.
virtual const transformer & transform() const
Return transformation between the coupled patches.
~PointEdgeWave()
Destructor.
void setPointInfo(const labelList &changedPoints, const List< Type > &changedPointsInfo)
Copy initial data into allPointInfo_.
const cyclicPolyPatch & nbrPatch() const
combineEqOp(TrackingData &td)
Definition: PointEdgeWave.C:55
Combination-Reduction operation for a parallel run. The information from all nodes is collected on th...
void reserve(const label)
Reserve allocation space for at least this size.
Definition: DynamicListI.H:152
const labelList & meshPoints() const
Return labelList of mesh points in patch. They are constructed.
Neighbour processor patch.
Input inter-processor communications stream operating on external buffer.
Definition: UIPstream.H:53
scalar y
A list of faces which address into the list of points.
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
Definition: DynamicList.H:56
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
label nPoints
Cyclic plane patch.
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:296
const Field< PointType > & points() const
Return reference to global points.
int neighbProcNo() const
Return neighbour processor number.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Output inter-processor communications stream operating on external buffer.
Definition: UOPstream.H:54
const tensor & T() const
Return the transformation tensor.
Definition: transformerI.H:94
static const char nl
Definition: Ostream.H:260
void operator()(Type &x, const Type &y) const
Definition: PointEdgeWave.C:60
Buffers for inter-processor communications streams (UOPstream, UIPstream).
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
void reverseDistribute(const label constructSize, List< T > &, const bool dummyTransform=true, const int tag=UPstream::msgType()) const
Reverse distribute data using default commsType.
const indirectPrimitivePatch & coupledPatch() const
Return patch of all coupled faces.
label constructSize() const
Constructed data size.
virtual const transformer & transform() const
Return null transform between processor patches.
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
label patchi
Class containing processor-to-processor mapping information.
bool equal(const T &s1, const T &s2)
Definition: doubleFloat.H:62
messageStream Info
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
label getUnsetPoints() const
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
void clear()
Clear the addressed list, i.e. set the size to zero.
Definition: DynamicListI.H:236
Wave propagation of information through grid. Every iteration information goes through one layer of e...
Definition: PointEdgeWave.H:85
label getUnsetEdges() const
Get number of unvisited edges, i.e. edges that were not (yet)
Namespace for OpenFOAM.
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:477