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