FaceCellWave.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 "FaceCellWave.H"
27 #include "polyMesh.H"
28 #include "processorPolyPatch.H"
29 #include "cyclicPolyPatch.H"
30 #include "OPstream.H"
31 #include "IPstream.H"
32 #include "PstreamReduceOps.H"
33 #include "debug.H"
34 #include "typeInfo.H"
35 #include "SubField.H"
36 #include "globalMeshData.H"
37 
38 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
39 
40 template<class Type, class TrackingData>
42 
43 template<class Type, class TrackingData>
45 
46 template<class Type, class TrackingData>
48 
49 
50 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
51 
52 template<class Type, class TrackingData>
54 (
55  const label celli,
56  const label neighbourFacei,
57  const Type& neighbourInfo,
58  const scalar tol,
59  Type& cellInfo
60 )
61 {
62  // Update info for celli, at position pt, with information from
63  // neighbouring face/cell.
64  // Updates:
65  // - changedCell_, changedCells_
66  // - statistics: nEvals_, nUnvisitedCells_
67 
68  nEvals_++;
69 
70  bool wasValid = cellInfo.valid(td_);
71 
72  bool propagate =
74  (
75  mesh_,
76  celli,
77  neighbourFacei,
78  neighbourInfo,
79  tol,
80  td_
81  );
82 
83  if (propagate)
84  {
85  if (!changedCell_[celli])
86  {
87  changedCell_[celli] = true;
88  changedCells_.append(celli);
89  }
90  }
91 
92  if (!wasValid && cellInfo.valid(td_))
93  {
94  --nUnvisitedCells_;
95  }
96 
97  return propagate;
98 }
99 
100 
101 template<class Type, class TrackingData>
103 (
104  const label facei,
105  const label neighbourCelli,
106  const Type& neighbourInfo,
107  const scalar tol,
108  Type& faceInfo
109 )
110 {
111  // Update info for facei, at position pt, with information from
112  // neighbouring face/cell.
113  // Updates:
114  // - changedFace_, changedFaces_,
115  // - statistics: nEvals_, nUnvisitedFaces_
116 
117  nEvals_++;
118 
119  bool wasValid = faceInfo.valid(td_);
120 
121  bool propagate =
122  faceInfo.updateFace
123  (
124  mesh_,
125  facei,
126  neighbourCelli,
127  neighbourInfo,
128  tol,
129  td_
130  );
131 
132  if (propagate)
133  {
134  if (!changedFace_[facei])
135  {
136  changedFace_[facei] = true;
137  changedFaces_.append(facei);
138  }
139  }
140 
141  if (!wasValid && faceInfo.valid(td_))
142  {
143  --nUnvisitedFaces_;
144  }
145 
146  return propagate;
147 }
148 
149 
150 template<class Type, class TrackingData>
152 (
153  const label facei,
154  const Type& neighbourInfo,
155  const scalar tol,
156  Type& faceInfo
157 )
158 {
159  // Update info for facei, at position pt, with information from
160  // same face.
161  // Updates:
162  // - changedFace_, changedFaces_,
163  // - statistics: nEvals_, nUnvisitedFaces_
164 
165  nEvals_++;
166 
167  bool wasValid = faceInfo.valid(td_);
168 
169  bool propagate =
170  faceInfo.updateFace
171  (
172  mesh_,
173  facei,
174  neighbourInfo,
175  tol,
176  td_
177  );
178 
179  if (propagate)
180  {
181  if (!changedFace_[facei])
182  {
183  changedFace_[facei] = true;
184  changedFaces_.append(facei);
185  }
186  }
187 
188  if (!wasValid && faceInfo.valid(td_))
189  {
190  --nUnvisitedFaces_;
191  }
192 
193  return propagate;
194 }
195 
196 
197 template<class Type, class TrackingData>
199 (
200  const polyPatch& patch
201 ) const
202 {
203  // For debugging: check status on both sides of cyclic
204 
205  const cyclicPolyPatch& nbrPatch =
206  refCast<const cyclicPolyPatch>(patch).nbrPatch();
207 
208  forAll(patch, patchFacei)
209  {
210  label i1 = patch.start() + patchFacei;
211  label i2 = nbrPatch.start() + patchFacei;
212 
213  if
214  (
215  !allFaceInfo_[i1].sameGeometry
216  (
217  mesh_,
218  allFaceInfo_[i2],
219  geomTol_,
220  td_
221  )
222  )
223  {
225  << " faceInfo:" << allFaceInfo_[i1]
226  << " otherfaceInfo:" << allFaceInfo_[i2]
227  << abort(FatalError);
228  }
229 
230  if (changedFace_[i1] != changedFace_[i2])
231  {
233  << " faceInfo:" << allFaceInfo_[i1]
234  << " otherfaceInfo:" << allFaceInfo_[i2]
235  << " changedFace:" << changedFace_[i1]
236  << " otherchangedFace:" << changedFace_[i2]
237  << abort(FatalError);
238  }
239  }
240 }
241 
242 
243 template<class Type, class TrackingData>
244 template<class PatchType>
246 {
247  forAll(mesh_.boundaryMesh(), patchi)
248  {
249  if (isA<PatchType>(mesh_.boundaryMesh()[patchi]))
250  {
251  return true;
252  }
253  }
254  return false;
255 }
256 
257 
258 template<class Type, class TrackingData>
260 (
261  const labelList& changedFaces,
262  const List<Type>& changedFacesInfo
263 )
264 {
265  forAll(changedFaces, changedFacei)
266  {
267  label facei = changedFaces[changedFacei];
268 
269  bool wasValid = allFaceInfo_[facei].valid(td_);
270 
271  // Copy info for facei
272  allFaceInfo_[facei] = changedFacesInfo[changedFacei];
273 
274  // Maintain count of unset faces
275  if (!wasValid && allFaceInfo_[facei].valid(td_))
276  {
277  --nUnvisitedFaces_;
278  }
279 
280  // Mark facei as changed, both on list and on face itself.
281 
282  changedFace_[facei] = true;
283  changedFaces_.append(facei);
284  }
285 }
286 
287 
288 template<class Type, class TrackingData>
290 (
291  const polyPatch& patch,
292  const label nFaces,
293  const labelList& changedFaces,
294  const List<Type>& changedFacesInfo
295 )
296 {
297  // Merge face information into member data
298 
299  for (label changedFacei = 0; changedFacei < nFaces; changedFacei++)
300  {
301  const Type& neighbourWallInfo = changedFacesInfo[changedFacei];
302  label patchFacei = changedFaces[changedFacei];
303 
304  label meshFacei = patch.start() + patchFacei;
305 
306  Type& currentWallInfo = allFaceInfo_[meshFacei];
307 
308  if (!currentWallInfo.equal(neighbourWallInfo, td_))
309  {
310  updateFace
311  (
312  meshFacei,
313  neighbourWallInfo,
314  propagationTol_,
315  currentWallInfo
316  );
317  }
318  }
319 }
320 
321 
322 template<class Type, class TrackingData>
324 (
325  const polyPatch& patch,
326  const label startFacei,
327  const label nFaces,
328  labelList& changedPatchFaces,
329  List<Type>& changedPatchFacesInfo
330 ) const
331 {
332  // Construct compact patchFace change arrays for a (slice of a) single
333  // patch. changedPatchFaces in local patch numbering.
334  // Return length of arrays.
335  label nChangedPatchFaces = 0;
336 
337  for (label i = 0; i < nFaces; i++)
338  {
339  label patchFacei = i + startFacei;
340 
341  label meshFacei = patch.start() + patchFacei;
342 
343  if (changedFace_[meshFacei])
344  {
345  changedPatchFaces[nChangedPatchFaces] = patchFacei;
346  changedPatchFacesInfo[nChangedPatchFaces] = allFaceInfo_[meshFacei];
347  nChangedPatchFaces++;
348  }
349  }
350  return nChangedPatchFaces;
351 }
352 
353 
354 template<class Type, class TrackingData>
356 (
357  const polyPatch& patch,
358  const label nFaces,
359  const labelList& faceLabels,
360  const transformer& transform,
361  List<Type>& faceInfo
362 )
363 {
364  for (label i = 0; i < nFaces; i++)
365  {
366  faceInfo[i].transform(patch, faceLabels[i], transform, td_);
367  }
368 }
369 
370 
371 template<class Type, class TrackingData>
373 {
374  // Transfer all the information to/from neighbouring processors
375 
376  const globalMeshData& pData = mesh_.globalData();
377 
378  // Which patches are processor patches
379  const labelList& procPatches = pData.processorPatches();
380 
381  // Send all
382 
383  PstreamBuffers pBufs(Pstream::commsTypes::nonBlocking);
384 
385  forAll(procPatches, i)
386  {
387  label patchi = procPatches[i];
388 
389  const processorPolyPatch& procPatch =
390  refCast<const processorPolyPatch>(mesh_.boundaryMesh()[patchi]);
391 
392  // Allocate buffers
393  label nSendFaces;
394  labelList sendFaces(procPatch.size());
395  List<Type> sendFacesInfo(procPatch.size());
396 
397  // Determine which faces changed on current patch
398  nSendFaces = getChangedPatchFaces
399  (
400  procPatch,
401  0,
402  procPatch.size(),
403  sendFaces,
404  sendFacesInfo
405  );
406 
407  if (debug & 2)
408  {
409  Pout<< " Processor patch " << patchi << ' ' << procPatch.name()
410  << " communicating with " << procPatch.neighbProcNo()
411  << " Sending:" << nSendFaces
412  << endl;
413  }
414 
415  // Send
416  UOPstream toNeighbour(procPatch.neighbProcNo(), pBufs);
417  toNeighbour
418  << SubList<label>(sendFaces, nSendFaces)
419  << SubList<Type>(sendFacesInfo, nSendFaces);
420  }
421 
422  pBufs.finishedSends();
423 
424  // Receive all
425 
426  forAll(procPatches, i)
427  {
428  label patchi = procPatches[i];
429 
430  const processorPolyPatch& procPatch =
431  refCast<const processorPolyPatch>(mesh_.boundaryMesh()[patchi]);
432 
433  // Allocate buffers
434  labelList receiveFaces;
435  List<Type> receiveFacesInfo;
436 
437  // Receive
438  {
439  UIPstream fromNeighbour(procPatch.neighbProcNo(), pBufs);
440  fromNeighbour >> receiveFaces >> receiveFacesInfo;
441  }
442 
443  if (debug & 2)
444  {
445  Pout<< " Processor patch " << patchi << ' ' << procPatch.name()
446  << " communicating with " << procPatch.neighbProcNo()
447  << " Receiving:" << receiveFaces.size()
448  << endl;
449  }
450 
451  // Transform info across the interface
452  transform
453  (
454  procPatch,
455  receiveFaces.size(),
456  receiveFaces,
457  procPatch.transform(),
458  receiveFacesInfo
459  );
460 
461  // Merge received info
462  mergeFaceInfo
463  (
464  procPatch,
465  receiveFaces.size(),
466  receiveFaces,
467  receiveFacesInfo
468  );
469  }
470 }
471 
472 
473 template<class Type, class TrackingData>
475 {
476  // Transfer information across cyclic halves.
477 
478  forAll(mesh_.boundaryMesh(), patchi)
479  {
480  const polyPatch& patch = mesh_.boundaryMesh()[patchi];
481 
482  if (isA<cyclicPolyPatch>(patch))
483  {
484  const cyclicPolyPatch& cycPatch =
485  refCast<const cyclicPolyPatch>(patch);
486  const cyclicPolyPatch& nbrPatch = cycPatch.nbrPatch();
487 
488  // Allocate buffers
489  label nReceiveFaces;
490  labelList receiveFaces(patch.size());
491  List<Type> receiveFacesInfo(patch.size());
492 
493  // Determine which faces changed
494  nReceiveFaces = getChangedPatchFaces
495  (
496  nbrPatch,
497  0,
498  nbrPatch.size(),
499  receiveFaces,
500  receiveFacesInfo
501  );
502 
503  if (debug & 2)
504  {
505  Pout<< " Cyclic patch " << patchi << ' ' << cycPatch.name()
506  << " Changed : " << nReceiveFaces
507  << endl;
508  }
509 
510  // Transform info across the interface
511  transform
512  (
513  cycPatch,
514  nReceiveFaces,
515  receiveFaces,
516  cycPatch.transform(),
517  receiveFacesInfo
518  );
519 
520  // Merge into global storage
521  mergeFaceInfo
522  (
523  cycPatch,
524  nReceiveFaces,
525  receiveFaces,
526  receiveFacesInfo
527  );
528 
529  if (debug)
530  {
531  checkCyclic(cycPatch);
532  }
533  }
534  }
535 }
536 
537 
538 template<class Type, class TrackingData>
540 {
541  // Collect changed information
542 
543  DynamicList<label> f0Baffle(explicitConnections_.size());
544  DynamicList<Type> f0Info(explicitConnections_.size());
545 
546  DynamicList<label> f1Baffle(explicitConnections_.size());
547  DynamicList<Type> f1Info(explicitConnections_.size());
548 
549  forAll(explicitConnections_, connI)
550  {
551  const labelPair& baffle = explicitConnections_[connI];
552 
553  label f0 = baffle[0];
554  if (changedFace_[f0])
555  {
556  f0Baffle.append(connI);
557  f0Info.append(allFaceInfo_[f0]);
558  }
559 
560  label f1 = baffle[1];
561  if (changedFace_[f1])
562  {
563  f1Baffle.append(connI);
564  f1Info.append(allFaceInfo_[f1]);
565  }
566  }
567 
568 
569  // Update other side with changed information
570 
571  forAll(f1Info, index)
572  {
573  const labelPair& baffle = explicitConnections_[f1Baffle[index]];
574 
575  label f0 = baffle[0];
576  Type& currentWallInfo = allFaceInfo_[f0];
577 
578  if (!currentWallInfo.equal(f1Info[index], td_))
579  {
580  updateFace
581  (
582  f0,
583  f1Info[index],
584  propagationTol_,
585  currentWallInfo
586  );
587  }
588  }
589 
590  forAll(f0Info, index)
591  {
592  const labelPair& baffle = explicitConnections_[f0Baffle[index]];
593 
594  label f1 = baffle[1];
595  Type& currentWallInfo = allFaceInfo_[f1];
596 
597  if (!currentWallInfo.equal(f0Info[index], td_))
598  {
599  updateFace
600  (
601  f1,
602  f0Info[index],
603  propagationTol_,
604  currentWallInfo
605  );
606  }
607  }
608 }
609 
610 
611 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
612 
613 template<class Type, class TrackingData>
615 (
616  const polyMesh& mesh,
617  UList<Type>& allFaceInfo,
618  UList<Type>& allCellInfo,
619  TrackingData& td
620 )
621 :
622  mesh_(mesh),
623  explicitConnections_(0),
624  allFaceInfo_(allFaceInfo),
625  allCellInfo_(allCellInfo),
626  td_(td),
627  changedFace_(mesh_.nFaces(), false),
628  changedFaces_(mesh_.nFaces()),
629  changedCell_(mesh_.nCells(), false),
630  changedCells_(mesh_.nCells()),
631  hasCyclicPatches_(hasPatch<cyclicPolyPatch>()),
632  nEvals_(0),
633  nUnvisitedCells_(mesh_.nCells()),
634  nUnvisitedFaces_(mesh_.nFaces())
635 {
636  if
637  (
639  || allCellInfo.size() != mesh_.nCells()
640  )
641  {
643  << "face and cell storage not the size of mesh faces, cells:"
644  << endl
645  << " allFaceInfo :" << allFaceInfo.size() << endl
646  << " mesh_.nFaces():" << mesh_.nFaces() << endl
647  << " allCellInfo :" << allCellInfo.size() << endl
648  << " mesh_.nCells():" << mesh_.nCells()
649  << exit(FatalError);
650  }
651 }
652 
653 
654 template<class Type, class TrackingData>
656 (
657  const polyMesh& mesh,
658  const labelList& changedFaces,
659  const List<Type>& changedFacesInfo,
660  UList<Type>& allFaceInfo,
661  UList<Type>& allCellInfo,
662  const label maxIter,
663  TrackingData& td
664 )
665 :
666  mesh_(mesh),
667  explicitConnections_(0),
668  allFaceInfo_(allFaceInfo),
669  allCellInfo_(allCellInfo),
670  td_(td),
671  changedFace_(mesh_.nFaces(), false),
672  changedFaces_(mesh_.nFaces()),
673  changedCell_(mesh_.nCells(), false),
674  changedCells_(mesh_.nCells()),
675  hasCyclicPatches_(hasPatch<cyclicPolyPatch>()),
676  nEvals_(0),
677  nUnvisitedCells_(mesh_.nCells()),
678  nUnvisitedFaces_(mesh_.nFaces())
679 {
680  if
681  (
683  || allCellInfo.size() != mesh_.nCells()
684  )
685  {
687  << "face and cell storage not the size of mesh faces, cells:"
688  << endl
689  << " allFaceInfo :" << allFaceInfo.size() << endl
690  << " mesh_.nFaces():" << mesh_.nFaces() << endl
691  << " allCellInfo :" << allCellInfo.size() << endl
692  << " mesh_.nCells():" << mesh_.nCells()
693  << exit(FatalError);
694  }
695 
696  // Copy initial changed faces data
697  setFaceInfo(changedFaces, changedFacesInfo);
698 
699  // Iterate until nothing changes
700  label iter = iterate(maxIter);
701 
702  if ((maxIter > 0) && (iter >= maxIter))
703  {
705  << "Maximum number of iterations reached. Increase maxIter." << endl
706  << " maxIter:" << maxIter << endl
707  << " nChangedCells:" << changedCells_.size() << endl
708  << " nChangedFaces:" << changedFaces_.size() << endl
709  << exit(FatalError);
710  }
711 }
712 
713 
714 template<class Type, class TrackingData>
716 (
717  const polyMesh& mesh,
718  const labelPairList& explicitConnections,
719  const labelList& changedFaces,
720  const List<Type>& changedFacesInfo,
721  UList<Type>& allFaceInfo,
722  UList<Type>& allCellInfo,
723  const label maxIter,
724  TrackingData& td
725 )
726 :
727  mesh_(mesh),
728  explicitConnections_(explicitConnections),
729  allFaceInfo_(allFaceInfo),
730  allCellInfo_(allCellInfo),
731  td_(td),
732  changedFace_(mesh_.nFaces(), false),
733  changedFaces_(mesh_.nFaces()),
734  changedCell_(mesh_.nCells(), false),
735  changedCells_(mesh_.nCells()),
736  hasCyclicPatches_(hasPatch<cyclicPolyPatch>()),
737  nEvals_(0),
738  nUnvisitedCells_(mesh_.nCells()),
739  nUnvisitedFaces_(mesh_.nFaces())
740 {
741  if
742  (
744  || allCellInfo.size() != mesh_.nCells()
745  )
746  {
748  << "face and cell storage not the size of mesh faces, cells:"
749  << endl
750  << " allFaceInfo :" << allFaceInfo.size() << endl
751  << " mesh_.nFaces():" << mesh_.nFaces() << endl
752  << " allCellInfo :" << allCellInfo.size() << endl
753  << " mesh_.nCells():" << mesh_.nCells()
754  << exit(FatalError);
755  }
756 
757  // Copy initial changed faces data
758  setFaceInfo(changedFaces, changedFacesInfo);
759 
760  // Iterate until nothing changes
761  label iter = iterate(maxIter);
762 
763  if ((maxIter > 0) && (iter >= maxIter))
764  {
766  << "Maximum number of iterations reached. Increase maxIter." << endl
767  << " maxIter:" << maxIter << endl
768  << " nChangedCells:" << changedCells_.size() << endl
769  << " nChangedFaces:" << changedFaces_.size() << endl
770  << exit(FatalError);
771  }
772 }
773 
774 
775 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
776 
777 template<class Type, class TrackingData>
779 {
780  return nUnvisitedCells_;
781 }
782 
783 
784 template<class Type, class TrackingData>
786 {
787  return nUnvisitedFaces_;
788 }
789 
790 
791 template<class Type, class TrackingData>
793 {
794  // Propagate face to cell
795 
796  const labelList& owner = mesh_.faceOwner();
797  const labelList& neighbour = mesh_.faceNeighbour();
798  label nInternalFaces = mesh_.nInternalFaces();
799 
800  forAll(changedFaces_, changedFacei)
801  {
802  label facei = changedFaces_[changedFacei];
803  if (!changedFace_[facei])
804  {
806  << "Face " << facei
807  << " not marked as having been changed"
808  << abort(FatalError);
809  }
810 
811 
812  const Type& neighbourWallInfo = allFaceInfo_[facei];
813 
814  // Evaluate all connected cells
815 
816  // Owner
817  label celli = owner[facei];
818  Type& currentWallInfo = allCellInfo_[celli];
819 
820  if (!currentWallInfo.equal(neighbourWallInfo, td_))
821  {
822  updateCell
823  (
824  celli,
825  facei,
826  neighbourWallInfo,
827  propagationTol_,
828  currentWallInfo
829  );
830  }
831 
832  // Neighbour.
833  if (facei < nInternalFaces)
834  {
835  celli = neighbour[facei];
836  Type& currentWallInfo2 = allCellInfo_[celli];
837 
838  if (!currentWallInfo2.equal(neighbourWallInfo, td_))
839  {
840  updateCell
841  (
842  celli,
843  facei,
844  neighbourWallInfo,
845  propagationTol_,
846  currentWallInfo2
847  );
848  }
849  }
850 
851  // Reset status of face
852  changedFace_[facei] = false;
853  }
854 
855  // Handled all changed faces by now
856  changedFaces_.clear();
857 
858  if (debug & 2)
859  {
860  Pout<< " Changed cells : " << changedCells_.size() << endl;
861  }
862 
863  // Sum changedCells over all procs
864  label totNChanged = changedCells_.size();
865 
866  reduce(totNChanged, sumOp<label>());
867 
868  return totNChanged;
869 }
870 
871 
872 template<class Type, class TrackingData>
874 {
875  // Propagate cell to face
876 
877  const cellList& cells = mesh_.cells();
878 
879  forAll(changedCells_, changedCelli)
880  {
881  label celli = changedCells_[changedCelli];
882  if (!changedCell_[celli])
883  {
885  << "Cell " << celli << " not marked as having been changed"
886  << abort(FatalError);
887  }
888 
889  const Type& neighbourWallInfo = allCellInfo_[celli];
890 
891  // Evaluate all connected faces
892 
893  const labelList& faceLabels = cells[celli];
894  forAll(faceLabels, faceLabelI)
895  {
896  label facei = faceLabels[faceLabelI];
897  Type& currentWallInfo = allFaceInfo_[facei];
898 
899  if (!currentWallInfo.equal(neighbourWallInfo, td_))
900  {
901  updateFace
902  (
903  facei,
904  celli,
905  neighbourWallInfo,
906  propagationTol_,
907  currentWallInfo
908  );
909  }
910  }
911 
912  // Reset status of cell
913  changedCell_[celli] = false;
914  }
915 
916  // Handled all changed cells by now
917  changedCells_.clear();
918 
919 
920  // Transfer across any explicitly provided internal connections
921  handleExplicitConnections();
922 
923  if (hasCyclicPatches_)
924  {
925  // Transfer changed faces across cyclic halves
926  handleCyclicPatches();
927  }
928 
929  if (Pstream::parRun())
930  {
931  // Transfer changed faces from neighbouring processors.
932  handleProcPatches();
933  }
934 
935  if (debug & 2)
936  {
937  Pout<< " Changed faces : " << changedFaces_.size() << endl;
938  }
939 
940  // Sum nChangedFaces over all procs
941  label totNChanged = changedFaces_.size();
942 
943  reduce(totNChanged, sumOp<label>());
944 
945  return totNChanged;
946 }
947 
948 
949 // Iterate
950 template<class Type, class TrackingData>
952 {
953  if (hasCyclicPatches_)
954  {
955  // Transfer changed faces across cyclic halves
956  handleCyclicPatches();
957  }
958 
959  if (Pstream::parRun())
960  {
961  // Transfer changed faces from neighbouring processors.
962  handleProcPatches();
963  }
964 
965  label iter = 0;
966 
967  while (iter < maxIter)
968  {
969  if (debug)
970  {
971  Info<< " Iteration " << iter << endl;
972  }
973 
974  nEvals_ = 0;
975 
976  label nCells = faceToCell();
977 
978  if (debug)
979  {
980  Info<< " Total changed cells : " << nCells << endl;
981  }
982 
983  if (nCells == 0)
984  {
985  break;
986  }
987 
988  label nFaces = cellToFace();
989 
990  if (debug)
991  {
992  Info<< " Total changed faces : " << nFaces << nl
993  << " Total evaluations : " << nEvals_ << nl
994  << " Remaining unvisited cells: " << nUnvisitedCells_ << nl
995  << " Remaining unvisited faces: " << nUnvisitedFaces_ << endl;
996  }
997 
998  if (nFaces == 0)
999  {
1000  break;
1001  }
1002 
1003  ++iter;
1004  }
1005 
1006  return iter;
1007 }
1008 
1009 
1010 // ************************************************************************* //
Inter-processor communication reduction functions.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:296
Wave propagation of information through grid. Every iteration information goes through one layer of c...
Definition: FaceCellWave.H:79
UList< Type > & allCellInfo()
Access allCellInfo.
Definition: FaceCellWave.H:302
DynamicList< label > changedCells_
Definition: FaceCellWave.H:118
void handleExplicitConnections()
Merge data across explicitly provided local connections (usually.
Definition: FaceCellWave.C:539
DynamicList< label > changedFaces_
List of changed faces.
Definition: FaceCellWave.H:112
void handleProcPatches()
Merge data from across processor boundaries.
Definition: FaceCellWave.C:372
label getUnsetCells() const
Get number of unvisited cells, i.e. cells that were not (yet)
Definition: FaceCellWave.C:778
void setFaceInfo(const labelList &changedFaces, const List< Type > &changedFacesInfo)
Set initial changed faces.
Definition: FaceCellWave.C:260
bool hasPatch() const
Has cyclic patch?
Definition: FaceCellWave.C:245
void transform(const polyPatch &patch, const label nFaces, const labelList &faceLabels, const transformer &transform, List< Type > &faceInfo)
Transform across an interface. Implementation referred to Type.
Definition: FaceCellWave.C:356
void mergeFaceInfo(const polyPatch &patch, const label nFaces, const labelList &, const List< Type > &)
Merge received patch data into global data.
Definition: FaceCellWave.C:290
UList< Type > & allFaceInfo()
Access allFaceInfo.
Definition: FaceCellWave.H:296
label getUnsetFaces() const
Get number of unvisited faces.
Definition: FaceCellWave.C:785
void checkCyclic(const polyPatch &pPatch) const
Debugging: check info on both sides of cyclic.
Definition: FaceCellWave.C:199
virtual label iterate(const label maxIter)
Iterate until no changes or maxIter reached. Returns actual.
Definition: FaceCellWave.C:951
void handleCyclicPatches()
Merge data from across cyclics.
Definition: FaceCellWave.C:474
bool updateCell(const label celli, const label neighbourFacei, const Type &neighbourInfo, const scalar tol, Type &cellInfo)
Updates cellInfo with information from neighbour. Updates all.
Definition: FaceCellWave.C:54
label getChangedPatchFaces(const polyPatch &patch, const label startFacei, const label nFaces, labelList &changedPatchFaces, List< Type > &changedPatchFacesInfo) const
Extract info for single patch only.
Definition: FaceCellWave.C:324
virtual label faceToCell()
Propagate from face to cell. Returns total number of cells.
Definition: FaceCellWave.C:792
const polyMesh & mesh_
Reference to mesh.
Definition: FaceCellWave.H:94
bool updateFace(const label facei, const label neighbourCelli, const Type &neighbourInfo, const scalar tol, Type &faceInfo)
Updates faceInfo with information from neighbour. Updates all.
Definition: FaceCellWave.C:103
virtual label cellToFace()
Propagate from cell to face. Returns total number of faces.
Definition: FaceCellWave.C:873
FaceCellWave(const polyMesh &, UList< Type > &allFaceInfo, UList< Type > &allCellInfo, TrackingData &td=defaultTrackingData_)
Definition: FaceCellWave.C:615
void append(const T &)
Append an element at the end of the list.
Definition: ListI.H:178
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Buffers for inter-processor communications streams (UOPstream, UIPstream).
void finishedSends(const bool block=true)
Mark all sends as having been done. This will start receives.
A List obtained as a section of another List.
Definition: SubList.H:56
Input inter-processor communications stream operating on external buffer.
Definition: UIPstream.H:57
label size() const
Return the number of elements in the UList.
Definition: UListI.H:311
Output inter-processor communications stream operating on external buffer.
Definition: UOPstream.H:58
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
Holds information regarding type of cell. Used in inside/outside determination in cellClassification.
Definition: cellInfo.H:67
bool updateCell(const polyMesh &, const label thisCelli, const label neighbourFacei, const cellInfo &neighbourInfo, const scalar tol, TrackingData &td)
Influence of neighbouring face.
Definition: cellInfoI.H:143
bool valid(TrackingData &td) const
Check whether origin has been changed at all or.
Definition: cellInfoI.H:110
A topoSetSource to select a faceSet from cells.
Definition: cellToFace.H:56
Cyclic plane patch.
const cyclicPolyPatch & nbrPatch() const
virtual const transformer & transform() const
Return transformation between the coupled patches.
A topoSetSource to select cells based on usage in faces.
Definition: faceToCell.H:52
Various mesh related information for a parallel run. Upon construction, constructs all info using par...
const labelList & processorPatches() const
Return list of processor patch labels.
const word & name() const
Return name.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:70
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:280
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:270
label nCells() const
label nFaces() const
Neighbour processor patch.
int neighbProcNo() const
Return neighbour processor number.
virtual const transformer & transform() const
Return null transform between processor patches.
Vector-tensor class used to perform translations, rotations and scaling operations in 3D space.
Definition: transformer.H:84
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
label patchi
const cellShapeList & cells
bool valid(const PtrList< ModelType > &l)
const dimensionedScalar e
Elementary charge.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
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)
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:504
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
error FatalError
static const char nl
Definition: Ostream.H:266
Basic run-time type information using word as the type's name. Used to enhance the standard RTTI to c...