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