syncToolsTemplates.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-2018 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 "syncTools.H"
27 #include "polyMesh.H"
28 #include "processorPolyPatch.H"
29 #include "cyclicPolyPatch.H"
30 #include "globalMeshData.H"
31 #include "contiguous.H"
32 #include "transform.H"
33 #include "transformList.H"
34 #include "SubField.H"
35 
36 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
37 
38 template<class T, class CombineOp>
39 void Foam::syncTools::combine
40 (
41  Map<T>& pointValues,
42  const CombineOp& cop,
43  const label index,
44  const T& val
45 )
46 {
47  typename Map<T>::iterator iter = pointValues.find(index);
48 
49  if (iter != pointValues.end())
50  {
51  cop(iter(), val);
52  }
53  else
54  {
55  pointValues.insert(index, val);
56  }
57 }
58 
59 
60 template<class T, class CombineOp>
61 void Foam::syncTools::combine
62 (
63  EdgeMap<T>& edgeValues,
64  const CombineOp& cop,
65  const edge& index,
66  const T& val
67 )
68 {
69  typename EdgeMap<T>::iterator iter = edgeValues.find(index);
70 
71  if (iter != edgeValues.end())
72  {
73  cop(iter(), val);
74  }
75  else
76  {
77  edgeValues.insert(index, val);
78  }
79 }
80 
81 
82 template<class T, class CombineOp, class TransformOp>
84 (
85  const polyMesh& mesh,
86  Map<T>& pointValues, // from mesh point label to value
87  const CombineOp& cop,
88  const TransformOp& top
89 )
90 {
91  const polyBoundaryMesh& patches = mesh.boundaryMesh();
92 
93  // Synchronize multiple shared points.
94  const globalMeshData& pd = mesh.globalData();
95 
96  // Values on shared points. Keyed on global shared index.
97  Map<T> sharedPointValues(0);
98 
99  if (pd.nGlobalPoints() > 0)
100  {
101  // meshPoint per local index
102  const labelList& sharedPtLabels = pd.sharedPointLabels();
103  // global shared index per local index
104  const labelList& sharedPtAddr = pd.sharedPointAddr();
105 
106  sharedPointValues.resize(sharedPtAddr.size());
107 
108  // Fill my entries in the shared points
109  forAll(sharedPtLabels, i)
110  {
111  label meshPointi = sharedPtLabels[i];
112 
113  typename Map<T>::const_iterator fnd =
114  pointValues.find(meshPointi);
115 
116  if (fnd != pointValues.end())
117  {
118  combine
119  (
120  sharedPointValues,
121  cop,
122  sharedPtAddr[i], // index
123  fnd() // value
124  );
125  }
126  }
127  }
128 
129 
130  if (Pstream::parRun())
131  {
133 
134  // Send
135 
136  forAll(patches, patchi)
137  {
138  if
139  (
140  isA<processorPolyPatch>(patches[patchi])
141  && patches[patchi].nPoints() > 0
142  )
143  {
144  const processorPolyPatch& procPatch =
145  refCast<const processorPolyPatch>(patches[patchi]);
146 
147  // Get data per patchPoint in neighbouring point numbers.
148 
149  const labelList& meshPts = procPatch.meshPoints();
150  const labelList& nbrPts = procPatch.neighbPoints();
151 
152  // Extract local values. Create map from nbrPoint to value.
153  // Note: how small initial size?
154  Map<T> patchInfo(meshPts.size() / 20);
155 
156  forAll(meshPts, i)
157  {
158  typename Map<T>::const_iterator iter =
159  pointValues.find(meshPts[i]);
160 
161  if (iter != pointValues.end())
162  {
163  patchInfo.insert(nbrPts[i], iter());
164  }
165  }
166 
167  UOPstream toNeighb(procPatch.neighbProcNo(), pBufs);
168  toNeighb << patchInfo;
169  }
170  }
171 
172  pBufs.finishedSends();
173 
174  // Receive and combine.
175 
176  forAll(patches, patchi)
177  {
178  if
179  (
180  isA<processorPolyPatch>(patches[patchi])
181  && patches[patchi].nPoints() > 0
182  )
183  {
184  const processorPolyPatch& procPatch =
185  refCast<const processorPolyPatch>(patches[patchi]);
186 
187  UIPstream fromNb(procPatch.neighbProcNo(), pBufs);
188  Map<T> nbrPatchInfo(fromNb);
189 
190  // Transform
191  top(procPatch, nbrPatchInfo);
192 
193  const labelList& meshPts = procPatch.meshPoints();
194 
195  // Only update those values which come from neighbour
196  forAllConstIter(typename Map<T>, nbrPatchInfo, nbrIter)
197  {
198  combine
199  (
200  pointValues,
201  cop,
202  meshPts[nbrIter.key()],
203  nbrIter()
204  );
205  }
206  }
207  }
208  }
209 
210  // Do the cyclics.
211  forAll(patches, patchi)
212  {
213  if (isA<cyclicPolyPatch>(patches[patchi]))
214  {
215  const cyclicPolyPatch& cycPatch =
216  refCast<const cyclicPolyPatch>(patches[patchi]);
217 
218  if (cycPatch.owner())
219  {
220  // Owner does all.
221 
222  const cyclicPolyPatch& nbrPatch = cycPatch.neighbPatch();
223  const edgeList& coupledPoints = cycPatch.coupledPoints();
224  const labelList& meshPtsA = cycPatch.meshPoints();
225  const labelList& meshPtsB = nbrPatch.meshPoints();
226 
227  // Extract local values. Create map from coupled-edge to value.
228  Map<T> half0Values(meshPtsA.size() / 20);
229  Map<T> half1Values(half0Values.size());
230 
231  forAll(coupledPoints, i)
232  {
233  const edge& e = coupledPoints[i];
234 
235  typename Map<T>::const_iterator point0Fnd =
236  pointValues.find(meshPtsA[e[0]]);
237 
238  if (point0Fnd != pointValues.end())
239  {
240  half0Values.insert(i, point0Fnd());
241  }
242 
243  typename Map<T>::const_iterator point1Fnd =
244  pointValues.find(meshPtsB[e[1]]);
245 
246  if (point1Fnd != pointValues.end())
247  {
248  half1Values.insert(i, point1Fnd());
249  }
250  }
251 
252  // Transform to receiving side
253  top(cycPatch, half1Values);
254  top(nbrPatch, half0Values);
255 
256  forAll(coupledPoints, i)
257  {
258  const edge& e = coupledPoints[i];
259 
260  typename Map<T>::const_iterator half0Fnd =
261  half0Values.find(i);
262 
263  if (half0Fnd != half0Values.end())
264  {
265  combine
266  (
267  pointValues,
268  cop,
269  meshPtsB[e[1]],
270  half0Fnd()
271  );
272  }
273 
274  typename Map<T>::const_iterator half1Fnd =
275  half1Values.find(i);
276 
277  if (half1Fnd != half1Values.end())
278  {
279  combine
280  (
281  pointValues,
282  cop,
283  meshPtsA[e[0]],
284  half1Fnd()
285  );
286  }
287  }
288  }
289  }
290  }
291 
292  // Synchronize multiple shared points.
293  if (pd.nGlobalPoints() > 0)
294  {
295  // meshPoint per local index
296  const labelList& sharedPtLabels = pd.sharedPointLabels();
297  // global shared index per local index
298  const labelList& sharedPtAddr = pd.sharedPointAddr();
299 
300  // Reduce on master.
301 
302  if (Pstream::parRun())
303  {
304  if (Pstream::master())
305  {
306  // Receive the edges using shared points from the slave.
307  for
308  (
309  int slave=Pstream::firstSlave();
310  slave<=Pstream::lastSlave();
311  slave++
312  )
313  {
314  IPstream fromSlave(Pstream::commsTypes::scheduled, slave);
315  Map<T> nbrValues(fromSlave);
316 
317  // Merge neighbouring values with my values
318  forAllConstIter(typename Map<T>, nbrValues, iter)
319  {
320  combine
321  (
322  sharedPointValues,
323  cop,
324  iter.key(), // edge
325  iter() // value
326  );
327  }
328  }
329 
330  // Send back
331  for
332  (
333  int slave=Pstream::firstSlave();
334  slave<=Pstream::lastSlave();
335  slave++
336  )
337  {
339  toSlave << sharedPointValues;
340  }
341  }
342  else
343  {
344  // Slave: send to master
345  {
346  OPstream toMaster
347  (
350  );
351  toMaster << sharedPointValues;
352  }
353  // Receive merged values
354  {
355  IPstream fromMaster
356  (
359  );
360  fromMaster >> sharedPointValues;
361  }
362  }
363  }
364 
365 
366  // Merge sharedPointValues (keyed on sharedPointAddr) into
367  // pointValues (keyed on mesh points).
368 
369  // Map from global shared index to meshpoint
370  Map<label> sharedToMeshPoint(2*sharedPtAddr.size());
371  forAll(sharedPtAddr, i)
372  {
373  sharedToMeshPoint.insert(sharedPtAddr[i], sharedPtLabels[i]);
374  }
375 
376  forAllConstIter(Map<label>, sharedToMeshPoint, iter)
377  {
378  // Do I have a value for my shared point
379  typename Map<T>::const_iterator sharedFnd =
380  sharedPointValues.find(iter.key());
381 
382  if (sharedFnd != sharedPointValues.end())
383  {
384  pointValues.set(iter(), sharedFnd());
385  }
386  }
387  }
388 }
389 
390 
391 template<class T, class CombineOp, class TransformOp>
393 (
394  const polyMesh& mesh,
395  EdgeMap<T>& edgeValues,
396  const CombineOp& cop,
397  const TransformOp& top
398 )
399 {
400  const polyBoundaryMesh& patches = mesh.boundaryMesh();
401 
402 
403  // Do synchronisation without constructing globalEdge addressing
404  // (since this constructs mesh edge addressing)
405 
406 
407  // Swap proc patch info
408  // ~~~~~~~~~~~~~~~~~~~~
409 
410  if (Pstream::parRun())
411  {
413 
414  // Send
415 
416  forAll(patches, patchi)
417  {
418  if
419  (
420  isA<processorPolyPatch>(patches[patchi])
421  && patches[patchi].nEdges() > 0
422  )
423  {
424  const processorPolyPatch& procPatch =
425  refCast<const processorPolyPatch>(patches[patchi]);
426 
427 
428  // Get data per patch edge in neighbouring edge.
429 
430  const edgeList& edges = procPatch.edges();
431  const labelList& meshPts = procPatch.meshPoints();
432  const labelList& nbrPts = procPatch.neighbPoints();
433 
434  EdgeMap<T> patchInfo(edges.size() / 20);
435 
436  forAll(edges, i)
437  {
438  const edge& e = edges[i];
439  const edge meshEdge(meshPts[e[0]], meshPts[e[1]]);
440 
441  typename EdgeMap<T>::const_iterator iter =
442  edgeValues.find(meshEdge);
443 
444  if (iter != edgeValues.end())
445  {
446  const edge nbrEdge(nbrPts[e[0]], nbrPts[e[1]]);
447  patchInfo.insert(nbrEdge, iter());
448  }
449  }
450 
451  UOPstream toNeighb(procPatch.neighbProcNo(), pBufs);
452  toNeighb << patchInfo;
453  }
454  }
455 
456  pBufs.finishedSends();
457 
458  // Receive and combine.
459 
460  forAll(patches, patchi)
461  {
462  if
463  (
464  isA<processorPolyPatch>(patches[patchi])
465  && patches[patchi].nEdges() > 0
466  )
467  {
468  const processorPolyPatch& procPatch =
469  refCast<const processorPolyPatch>(patches[patchi]);
470 
471  EdgeMap<T> nbrPatchInfo;
472  {
473  UIPstream fromNbr(procPatch.neighbProcNo(), pBufs);
474  fromNbr >> nbrPatchInfo;
475  }
476 
477  // Apply transform to convert to this side properties.
478  top(procPatch, nbrPatchInfo);
479 
480 
481  // Only update those values which come from neighbour
482  const labelList& meshPts = procPatch.meshPoints();
483 
484  forAllConstIter(typename EdgeMap<T>, nbrPatchInfo, nbrIter)
485  {
486  const edge& e = nbrIter.key();
487  const edge meshEdge(meshPts[e[0]], meshPts[e[1]]);
488 
489  combine
490  (
491  edgeValues,
492  cop,
493  meshEdge, // edge
494  nbrIter() // value
495  );
496  }
497  }
498  }
499  }
500 
501 
502  // Swap cyclic info
503  // ~~~~~~~~~~~~~~~~
504 
505  forAll(patches, patchi)
506  {
507  if (isA<cyclicPolyPatch>(patches[patchi]))
508  {
509  const cyclicPolyPatch& cycPatch =
510  refCast<const cyclicPolyPatch>(patches[patchi]);
511 
512  if (cycPatch.owner())
513  {
514  // Owner does all.
515 
516  const edgeList& coupledEdges = cycPatch.coupledEdges();
517  const labelList& meshPtsA = cycPatch.meshPoints();
518  const edgeList& edgesA = cycPatch.edges();
519  const cyclicPolyPatch& nbrPatch = cycPatch.neighbPatch();
520  const labelList& meshPtsB = nbrPatch.meshPoints();
521  const edgeList& edgesB = nbrPatch.edges();
522 
523  // Extract local values. Create map from edge to value.
524  Map<T> half0Values(edgesA.size() / 20);
525  Map<T> half1Values(half0Values.size());
526 
527  forAll(coupledEdges, i)
528  {
529  const edge& twoEdges = coupledEdges[i];
530 
531  {
532  const edge& e0 = edgesA[twoEdges[0]];
533  const edge meshEdge0(meshPtsA[e0[0]], meshPtsA[e0[1]]);
534 
535  typename EdgeMap<T>::const_iterator iter =
536  edgeValues.find(meshEdge0);
537 
538  if (iter != edgeValues.end())
539  {
540  half0Values.insert(i, iter());
541  }
542  }
543  {
544  const edge& e1 = edgesB[twoEdges[1]];
545  const edge meshEdge1(meshPtsB[e1[0]], meshPtsB[e1[1]]);
546 
547  typename EdgeMap<T>::const_iterator iter =
548  edgeValues.find(meshEdge1);
549 
550  if (iter != edgeValues.end())
551  {
552  half1Values.insert(i, iter());
553  }
554  }
555  }
556 
557  // Transform to this side
558  top(cycPatch, half1Values);
559  top(nbrPatch, half0Values);
560 
561 
562  // Extract and combine information
563 
564  forAll(coupledEdges, i)
565  {
566  const edge& twoEdges = coupledEdges[i];
567 
568  typename Map<T>::const_iterator half1Fnd =
569  half1Values.find(i);
570 
571  if (half1Fnd != half1Values.end())
572  {
573  const edge& e0 = edgesA[twoEdges[0]];
574  const edge meshEdge0(meshPtsA[e0[0]], meshPtsA[e0[1]]);
575 
576  combine
577  (
578  edgeValues,
579  cop,
580  meshEdge0, // edge
581  half1Fnd() // value
582  );
583  }
584 
585  typename Map<T>::const_iterator half0Fnd =
586  half0Values.find(i);
587  if (half0Fnd != half0Values.end())
588  {
589  const edge& e1 = edgesB[twoEdges[1]];
590  const edge meshEdge1(meshPtsB[e1[0]], meshPtsB[e1[1]]);
591 
592  combine
593  (
594  edgeValues,
595  cop,
596  meshEdge1, // edge
597  half0Fnd() // value
598  );
599  }
600  }
601  }
602  }
603  }
604 
605  // Synchronize multiple shared points.
606  // Problem is that we don't want to construct shared edges so basically
607  // we do here like globalMeshData but then using sparse edge representation
608  // (EdgeMap instead of mesh.edges())
609 
610  const globalMeshData& pd = mesh.globalData();
611  const labelList& sharedPtAddr = pd.sharedPointAddr();
612  const labelList& sharedPtLabels = pd.sharedPointLabels();
613 
614  // 1. Create map from meshPoint to globalShared index.
615  Map<label> meshToShared(2*sharedPtLabels.size());
616  forAll(sharedPtLabels, i)
617  {
618  meshToShared.insert(sharedPtLabels[i], sharedPtAddr[i]);
619  }
620 
621  // Values on shared points. From two sharedPtAddr indices to a value.
622  EdgeMap<T> sharedEdgeValues(meshToShared.size());
623 
624  // From shared edge to mesh edge. Used for merging later on.
625  EdgeMap<edge> potentialSharedEdge(meshToShared.size());
626 
627  // 2. Find any edges using two global shared points. These will always be
628  // on the outside of the mesh. (though might not be on coupled patch
629  // if is single edge and not on coupled face)
630  // Store value (if any) on sharedEdgeValues
631  for (label facei = mesh.nInternalFaces(); facei < mesh.nFaces(); facei++)
632  {
633  const face& f = mesh.faces()[facei];
634 
635  forAll(f, fp)
636  {
637  label v0 = f[fp];
638  label v1 = f[f.fcIndex(fp)];
639 
640  Map<label>::const_iterator v0Fnd = meshToShared.find(v0);
641 
642  if (v0Fnd != meshToShared.end())
643  {
644  Map<label>::const_iterator v1Fnd = meshToShared.find(v1);
645 
646  if (v1Fnd != meshToShared.end())
647  {
648  const edge meshEdge(v0, v1);
649 
650  // edge in shared point labels
651  const edge sharedEdge(v0Fnd(), v1Fnd());
652 
653  // Store mesh edge as a potential shared edge.
654  potentialSharedEdge.insert(sharedEdge, meshEdge);
655 
656  typename EdgeMap<T>::const_iterator edgeFnd =
657  edgeValues.find(meshEdge);
658 
659  if (edgeFnd != edgeValues.end())
660  {
661  // edge exists in edgeValues. See if already in map
662  // (since on same processor, e.g. cyclic)
663  combine
664  (
665  sharedEdgeValues,
666  cop,
667  sharedEdge, // edge
668  edgeFnd() // value
669  );
670  }
671  }
672  }
673  }
674  }
675 
676 
677  // Now sharedEdgeValues will contain per potential sharedEdge the value.
678  // (potential since an edge having two shared points is not necessary a
679  // shared edge).
680  // Reduce this on the master.
681 
682  if (Pstream::parRun())
683  {
684  if (Pstream::master())
685  {
686  // Receive the edges using shared points from the slave.
687  for
688  (
689  int slave=Pstream::firstSlave();
690  slave<=Pstream::lastSlave();
691  slave++
692  )
693  {
694  IPstream fromSlave(Pstream::commsTypes::scheduled, slave);
695  EdgeMap<T> nbrValues(fromSlave);
696 
697  // Merge neighbouring values with my values
698  forAllConstIter(typename EdgeMap<T>, nbrValues, iter)
699  {
700  combine
701  (
702  sharedEdgeValues,
703  cop,
704  iter.key(), // edge
705  iter() // value
706  );
707  }
708  }
709 
710  // Send back
711  for
712  (
713  int slave=Pstream::firstSlave();
714  slave<=Pstream::lastSlave();
715  slave++
716  )
717  {
718 
720  toSlave << sharedEdgeValues;
721  }
722  }
723  else
724  {
725  // Send to master
726  {
727  OPstream toMaster
728  (
731  );
732  toMaster << sharedEdgeValues;
733  }
734  // Receive merged values
735  {
736  IPstream fromMaster
737  (
740  );
741  fromMaster >> sharedEdgeValues;
742  }
743  }
744  }
745 
746 
747  // Merge sharedEdgeValues (keyed on sharedPointAddr) into edgeValues
748  // (keyed on mesh points).
749 
750  // Loop over all my shared edges.
751  forAllConstIter(typename EdgeMap<edge>, potentialSharedEdge, iter)
752  {
753  const edge& sharedEdge = iter.key();
754  const edge& meshEdge = iter();
755 
756  // Do I have a value for the shared edge?
757  typename EdgeMap<T>::const_iterator sharedFnd =
758  sharedEdgeValues.find(sharedEdge);
759 
760  if (sharedFnd != sharedEdgeValues.end())
761  {
762  combine
763  (
764  edgeValues,
765  cop,
766  meshEdge, // edge
767  sharedFnd() // value
768  );
769  }
770  }
771 }
772 
773 
774 //template<class T, class CombineOp, class TransformOp>
775 //void Foam::syncTools::syncPointList
776 //(
777 // const polyMesh& mesh,
778 // List<T>& pointValues,
779 // const CombineOp& cop,
780 // const T& nullValue,
781 // const TransformOp& top
782 //)
783 //{
784 // if (pointValues.size() != mesh.nPoints())
785 // {
786 // FatalErrorInFunction
787 // << "Number of values " << pointValues.size()
788 // << " is not equal to the number of points in the mesh "
789 // << mesh.nPoints() << abort(FatalError);
790 // }
791 //
792 // const polyBoundaryMesh& patches = mesh.boundaryMesh();
793 //
794 // // Synchronize multiple shared points.
795 // const globalMeshData& pd = mesh.globalData();
796 //
797 // // Values on shared points.
798 // Field<T> sharedPts(0);
799 // if (pd.nGlobalPoints() > 0)
800 // {
801 // // Values on shared points.
802 // sharedPts.setSize(pd.nGlobalPoints(), nullValue);
803 //
804 // forAll(pd.sharedPointLabels(), i)
805 // {
806 // label meshPointi = pd.sharedPointLabels()[i];
807 // // Fill my entries in the shared points
808 // sharedPts[pd.sharedPointAddr()[i]] = pointValues[meshPointi];
809 // }
810 // }
811 //
812 // if (Pstream::parRun())
813 // {
814 // PstreamBuffers pBufs(Pstream::commsTypes::nonBlocking);
815 //
816 // // Send
817 //
818 // forAll(patches, patchi)
819 // {
820 // if
821 // (
822 // isA<processorPolyPatch>(patches[patchi])
823 // && patches[patchi].nPoints() > 0
824 // )
825 // {
826 // const processorPolyPatch& procPatch =
827 // refCast<const processorPolyPatch>(patches[patchi]);
828 //
829 // // Get data per patchPoint in neighbouring point numbers.
830 // Field<T> patchInfo(procPatch.nPoints());
831 //
832 // const labelList& meshPts = procPatch.meshPoints();
833 // const labelList& nbrPts = procPatch.neighbPoints();
834 //
835 // forAll(nbrPts, pointi)
836 // {
837 // label nbrPointi = nbrPts[pointi];
838 // patchInfo[nbrPointi] = pointValues[meshPts[pointi]];
839 // }
840 //
841 // UOPstream toNbr(procPatch.neighbProcNo(), pBufs);
842 // toNbr << patchInfo;
843 // }
844 // }
845 //
846 // pBufs.finishedSends();
847 //
848 // // Receive and combine.
849 //
850 // forAll(patches, patchi)
851 // {
852 // if
853 // (
854 // isA<processorPolyPatch>(patches[patchi])
855 // && patches[patchi].nPoints() > 0
856 // )
857 // {
858 // const processorPolyPatch& procPatch =
859 // refCast<const processorPolyPatch>(patches[patchi]);
860 //
861 // Field<T> nbrPatchInfo(procPatch.nPoints());
862 // {
863 // UIPstream fromNbr(procPatch.neighbProcNo(), pBufs);
864 // fromNbr >> nbrPatchInfo;
865 // }
866 //
867 // // Transform to this side
868 // top(procPatch, nbrPatchInfo);
869 //
870 // const labelList& meshPts = procPatch.meshPoints();
871 //
872 // forAll(meshPts, pointi)
873 // {
874 // label meshPointi = meshPts[pointi];
875 // cop(pointValues[meshPointi], nbrPatchInfo[pointi]);
876 // }
877 // }
878 // }
879 // }
880 //
881 // // Do the cyclics.
882 // forAll(patches, patchi)
883 // {
884 // if (isA<cyclicPolyPatch>(patches[patchi]))
885 // {
886 // const cyclicPolyPatch& cycPatch =
887 // refCast<const cyclicPolyPatch>(patches[patchi]);
888 //
889 // if (cycPatch.owner())
890 // {
891 // // Owner does all.
892 //
893 // const edgeList& coupledPoints = cycPatch.coupledPoints();
894 // const labelList& meshPts = cycPatch.meshPoints();
895 // const cyclicPolyPatch& nbrPatch = cycPatch.neighbPatch();
896 // const labelList& nbrMeshPoints = nbrPatch.meshPoints();
897 //
898 // Field<T> half0Values(coupledPoints.size());
899 // Field<T> half1Values(coupledPoints.size());
900 //
901 // forAll(coupledPoints, i)
902 // {
903 // const edge& e = coupledPoints[i];
904 // half0Values[i] = pointValues[meshPts[e[0]]];
905 // half1Values[i] = pointValues[nbrMeshPoints[e[1]]];
906 // }
907 //
908 // // SubField<T> slice0(half0Values, half0Values.size());
909 // // SubField<T> slice1(half1Values, half1Values.size());
910 // // top(cycPatch, reinterpret_cast<Field<T>&>(slice1));
911 // // top(nbrPatch, reinterpret_cast<Field<T>&>(slice0));
912 //
913 // top(cycPatch, half1Values);
914 // top(nbrPatch, half0Values);
915 //
916 // forAll(coupledPoints, i)
917 // {
918 // const edge& e = coupledPoints[i];
919 // cop(pointValues[meshPts[e[0]]], half1Values[i]);
920 // cop(pointValues[nbrMeshPoints[e[1]]], half0Values[i]);
921 // }
922 // }
923 // }
924 // }
925 //
926 // // Synchronize multiple shared points.
927 // const globalMeshData& pd = mesh.globalData();
928 //
929 // if (pd.nGlobalPoints() > 0)
930 // {
931 // // Combine on master.
932 // Pstream::listCombineGather(sharedPts, cop);
933 // Pstream::listCombineScatter(sharedPts);
934 //
935 // // Now we will all have the same information. Merge it back with
936 // // my local information.
937 // forAll(pd.sharedPointLabels(), i)
938 // {
939 // label meshPointi = pd.sharedPointLabels()[i];
940 // pointValues[meshPointi] = sharedPts[pd.sharedPointAddr()[i]];
941 // }
942 // }
943 //}
944 
945 
946 //template<class T, class CombineOp, class TransformOp>
947 //void Foam::syncTools::syncPointList
948 //(
949 // const polyMesh& mesh,
950 // const labelList& meshPoints,
951 // List<T>& pointValues,
952 // const CombineOp& cop,
953 // const T& nullValue,
954 // const TransformOp& top
955 //)
956 //{
957 // if (pointValues.size() != meshPoints.size())
958 // {
959 // FatalErrorInFunction
960 // << "Number of values " << pointValues.size()
961 // << " is not equal to the number of points "
962 // << meshPoints.size() << abort(FatalError);
963 // }
964 //
965 // if (!hasCouples(mesh.boundaryMesh()))
966 // {
967 // return;
968 // }
969 //
970 // Field<T> meshValues(mesh.nPoints(), nullValue);
971 //
972 // forAll(meshPoints, i)
973 // {
974 // meshValues[meshPoints[i]] = pointValues[i];
975 // }
976 //
977 // syncTools::syncPointList
978 // (
979 // mesh,
980 // meshValues,
981 // cop, // combine op
982 // nullValue, // null value
983 // top // position or field
984 // );
985 //
986 // forAll(meshPoints, i)
987 // {
988 // pointValues[i] = meshValues[meshPoints[i]];
989 // }
990 //}
991 
992 template<class T, class CombineOp, class TransformOp>
994 (
995  const polyMesh& mesh,
996  List<T>& pointValues,
997  const CombineOp& cop,
998  const T& nullValue,
999  const TransformOp& top
1000 )
1001 {
1002  if (pointValues.size() != mesh.nPoints())
1003  {
1005  << "Number of values " << pointValues.size()
1006  << " is not equal to the number of points in the mesh "
1007  << mesh.nPoints() << abort(FatalError);
1008  }
1009 
1010  mesh.globalData().syncPointData(pointValues, cop, top);
1011 }
1012 
1013 
1014 //template<class CombineOp>
1015 //void Foam::syncTools::syncPointPositions
1016 //(
1017 // const polyMesh& mesh,
1018 // List<point>& pointValues,
1019 // const CombineOp& cop,
1020 // const point& nullValue
1021 //)
1022 //{
1023 // if (pointValues.size() != mesh.nPoints())
1024 // {
1025 // FatalErrorInFunction
1026 // << "Number of values " << pointValues.size()
1027 // << " is not equal to the number of points in the mesh "
1028 // << mesh.nPoints() << abort(FatalError);
1029 // }
1030 //
1031 // mesh.globalData().syncPointData(pointValues, cop, true);
1032 //}
1033 
1034 
1035 template<class T, class CombineOp, class TransformOp>
1038  const polyMesh& mesh,
1039  const labelList& meshPoints,
1040  List<T>& pointValues,
1041  const CombineOp& cop,
1042  const T& nullValue,
1043  const TransformOp& top
1044 )
1045 {
1046  if (pointValues.size() != meshPoints.size())
1047  {
1049  << "Number of values " << pointValues.size()
1050  << " is not equal to the number of meshPoints "
1051  << meshPoints.size() << abort(FatalError);
1052  }
1053  const globalMeshData& gd = mesh.globalData();
1054  const indirectPrimitivePatch& cpp = gd.coupledPatch();
1055  const Map<label>& mpm = cpp.meshPointMap();
1056 
1057  List<T> cppFld(cpp.nPoints(), nullValue);
1058 
1059  forAll(meshPoints, i)
1060  {
1061  label pointi = meshPoints[i];
1062  Map<label>::const_iterator iter = mpm.find(pointi);
1063  if (iter != mpm.end())
1064  {
1065  cppFld[iter()] = pointValues[i];
1066  }
1067  }
1068 
1070  (
1071  cppFld,
1072  gd.globalPointSlaves(),
1074  gd.globalPointSlavesMap(),
1075  gd.globalTransforms(),
1076  cop,
1077  top
1078  );
1079 
1080  forAll(meshPoints, i)
1081  {
1082  label pointi = meshPoints[i];
1083  Map<label>::const_iterator iter = mpm.find(pointi);
1084  if (iter != mpm.end())
1085  {
1086  pointValues[i] = cppFld[iter()];
1087  }
1088  }
1089 }
1090 
1091 
1092 //template<class CombineOp>
1093 //void Foam::syncTools::syncPointPositions
1094 //(
1095 // const polyMesh& mesh,
1096 // const labelList& meshPoints,
1097 // List<point>& pointValues,
1098 // const CombineOp& cop,
1099 // const point& nullValue
1100 //)
1101 //{
1102 // if (pointValues.size() != meshPoints.size())
1103 // {
1104 // FatalErrorInFunction
1105 // << "Number of values " << pointValues.size()
1106 // << " is not equal to the number of meshPoints "
1107 // << meshPoints.size() << abort(FatalError);
1108 // }
1109 // const globalMeshData& gd = mesh.globalData();
1110 // const indirectPrimitivePatch& cpp = gd.coupledPatch();
1111 // const Map<label>& mpm = cpp.meshPointMap();
1112 //
1113 // List<point> cppFld(cpp.nPoints(), nullValue);
1114 //
1115 // forAll(meshPoints, i)
1116 // {
1117 // label pointi = meshPoints[i];
1118 // Map<label>::const_iterator iter = mpm.find(pointi);
1119 // if (iter != mpm.end())
1120 // {
1121 // cppFld[iter()] = pointValues[i];
1122 // }
1123 // }
1124 //
1125 // globalMeshData::syncData
1126 // (
1127 // cppFld,
1128 // gd.globalPointSlaves(),
1129 // gd.globalPointTransformedSlaves(),
1130 // gd.globalPointSlavesMap(),
1131 // gd.globalTransforms(),
1132 // cop,
1133 // true, // position?
1134 // mapDistribute::transform() // not used
1135 // );
1136 //
1137 // forAll(meshPoints, i)
1138 // {
1139 // label pointi = meshPoints[i];
1140 // Map<label>::const_iterator iter = mpm.find(pointi);
1141 // if (iter != mpm.end())
1142 // {
1143 // pointValues[i] = cppFld[iter()];
1144 // }
1145 // }
1146 //}
1147 
1148 
1149 template<class T, class CombineOp, class TransformOp>
1152  const polyMesh& mesh,
1153  List<T>& edgeValues,
1154  const CombineOp& cop,
1155  const T& nullValue,
1156  const TransformOp& top
1157 )
1158 {
1159  if (edgeValues.size() != mesh.nEdges())
1160  {
1162  << "Number of values " << edgeValues.size()
1163  << " is not equal to the number of edges in the mesh "
1164  << mesh.nEdges() << abort(FatalError);
1165  }
1166 
1167  const globalMeshData& gd = mesh.globalData();
1168  const labelList& meshEdges = gd.coupledPatchMeshEdges();
1169  const globalIndexAndTransform& git = gd.globalTransforms();
1170  const mapDistribute& edgeMap = gd.globalEdgeSlavesMap();
1171 
1172  List<T> cppFld(UIndirectList<T>(edgeValues, meshEdges));
1173 
1175  (
1176  cppFld,
1177  gd.globalEdgeSlaves(),
1179  edgeMap,
1180  git,
1181  cop,
1182  top
1183  );
1184 
1185  // Extract back onto mesh
1186  forAll(meshEdges, i)
1187  {
1188  edgeValues[meshEdges[i]] = cppFld[i];
1189  }
1190 }
1191 
1192 
1193 //template<class CombineOp>
1194 //void Foam::syncTools::syncEdgePositions
1195 //(
1196 // const polyMesh& mesh,
1197 // List<point>& edgeValues,
1198 // const CombineOp& cop,
1199 // const point& nullValue
1200 //)
1201 //{
1202 // if (edgeValues.size() != mesh.nEdges())
1203 // {
1204 // FatalErrorInFunction
1205 // << "Number of values " << edgeValues.size()
1206 // << " is not equal to the number of edges in the mesh "
1207 // << mesh.nEdges() << abort(FatalError);
1208 // }
1209 //
1210 // const globalMeshData& gd = mesh.globalData();
1211 // const labelList& meshEdges = gd.coupledPatchMeshEdges();
1212 // const globalIndexAndTransform& git = gd.globalTransforms();
1213 // const mapDistribute& map = gd.globalEdgeSlavesMap();
1214 //
1215 // List<point> cppFld(UIndirectList<point>(edgeValues, meshEdges));
1216 //
1217 // globalMeshData::syncData
1218 // (
1219 // cppFld,
1220 // gd.globalEdgeSlaves(),
1221 // gd.globalEdgeTransformedSlaves(),
1222 // map,
1223 // git,
1224 // cop,
1225 // true, // position?
1226 // mapDistribute::transform() // not used
1227 // );
1228 //
1229 // // Extract back onto mesh
1230 // forAll(meshEdges, i)
1231 // {
1232 // edgeValues[meshEdges[i]] = cppFld[i];
1233 // }
1234 //}
1235 
1236 
1237 template<class T, class CombineOp, class TransformOp>
1240  const polyMesh& mesh,
1241  const labelList& meshEdges,
1242  List<T>& edgeValues,
1243  const CombineOp& cop,
1244  const T& nullValue,
1245  const TransformOp& top
1246 )
1247 {
1248  if (edgeValues.size() != meshEdges.size())
1249  {
1251  << "Number of values " << edgeValues.size()
1252  << " is not equal to the number of meshEdges "
1253  << meshEdges.size() << abort(FatalError);
1254  }
1255  const globalMeshData& gd = mesh.globalData();
1256  const indirectPrimitivePatch& cpp = gd.coupledPatch();
1257  const Map<label>& mpm = gd.coupledPatchMeshEdgeMap();
1258 
1259  List<T> cppFld(cpp.nEdges(), nullValue);
1260 
1261  forAll(meshEdges, i)
1262  {
1263  label edgeI = meshEdges[i];
1264  Map<label>::const_iterator iter = mpm.find(edgeI);
1265  if (iter != mpm.end())
1266  {
1267  cppFld[iter()] = edgeValues[i];
1268  }
1269  }
1270 
1272  (
1273  cppFld,
1274  gd.globalEdgeSlaves(),
1276  gd.globalEdgeSlavesMap(),
1277  gd.globalTransforms(),
1278  cop,
1279  top
1280  );
1281 
1282  forAll(meshEdges, i)
1283  {
1284  label edgeI = meshEdges[i];
1285  Map<label>::const_iterator iter = mpm.find(edgeI);
1286  if (iter != mpm.end())
1287  {
1288  edgeValues[i] = cppFld[iter()];
1289  }
1290  }
1291 }
1292 
1293 template<class T, class CombineOp, class TransformOp>
1296  const polyMesh& mesh,
1297  UList<T>& faceValues,
1298  const CombineOp& cop,
1299  const TransformOp& top,
1300  const bool parRun
1301 )
1302 {
1303  const label nBFaces = mesh.nFaces() - mesh.nInternalFaces();
1304 
1305  if (faceValues.size() != nBFaces)
1306  {
1308  << "Number of values " << faceValues.size()
1309  << " is not equal to the number of boundary faces in the mesh "
1310  << nBFaces << abort(FatalError);
1311  }
1312 
1313  const polyBoundaryMesh& patches = mesh.boundaryMesh();
1314 
1315  if (parRun)
1316  {
1318 
1319  // Send
1320 
1321  forAll(patches, patchi)
1322  {
1323  if
1324  (
1325  isA<processorPolyPatch>(patches[patchi])
1326  && patches[patchi].size() > 0
1327  )
1328  {
1329  const processorPolyPatch& procPatch =
1330  refCast<const processorPolyPatch>(patches[patchi]);
1331 
1332  label patchStart = procPatch.start()-mesh.nInternalFaces();
1333 
1334  UOPstream toNbr(procPatch.neighbProcNo(), pBufs);
1335  toNbr << SubField<T>(faceValues, procPatch.size(), patchStart);
1336  }
1337  }
1338 
1339 
1340  pBufs.finishedSends();
1341 
1342 
1343  // Receive and combine.
1344 
1345  forAll(patches, patchi)
1346  {
1347  if
1348  (
1349  isA<processorPolyPatch>(patches[patchi])
1350  && patches[patchi].size() > 0
1351  )
1352  {
1353  const processorPolyPatch& procPatch =
1354  refCast<const processorPolyPatch>(patches[patchi]);
1355 
1356  Field<T> nbrPatchInfo(procPatch.size());
1357 
1358  UIPstream fromNeighb(procPatch.neighbProcNo(), pBufs);
1359  fromNeighb >> nbrPatchInfo;
1360 
1361  top(procPatch, nbrPatchInfo);
1362 
1363  label bFacei = procPatch.start()-mesh.nInternalFaces();
1364 
1365  forAll(nbrPatchInfo, i)
1366  {
1367  cop(faceValues[bFacei++], nbrPatchInfo[i]);
1368  }
1369  }
1370  }
1371  }
1372 
1373  // Do the cyclics.
1374  forAll(patches, patchi)
1375  {
1376  if (isA<cyclicPolyPatch>(patches[patchi]))
1377  {
1378  const cyclicPolyPatch& cycPatch =
1379  refCast<const cyclicPolyPatch>(patches[patchi]);
1380 
1381  if (cycPatch.owner())
1382  {
1383  // Owner does all.
1384  const cyclicPolyPatch& nbrPatch = cycPatch.neighbPatch();
1385  label ownStart = cycPatch.start()-mesh.nInternalFaces();
1386  label nbrStart = nbrPatch.start()-mesh.nInternalFaces();
1387 
1388  label sz = cycPatch.size();
1389 
1390  // Transform (copy of) data on both sides
1391  Field<T> ownVals(SubField<T>(faceValues, sz, ownStart));
1392  top(nbrPatch, ownVals);
1393 
1394  Field<T> nbrVals(SubField<T>(faceValues, sz, nbrStart));
1395  top(cycPatch, nbrVals);
1396 
1397  label i0 = ownStart;
1398  forAll(nbrVals, i)
1399  {
1400  cop(faceValues[i0++], nbrVals[i]);
1401  }
1402 
1403  label i1 = nbrStart;
1404  forAll(ownVals, i)
1405  {
1406  cop(faceValues[i1++], ownVals[i]);
1407  }
1408  }
1409  }
1410  }
1411 }
1412 
1413 
1414 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
1415 
1416 template<unsigned nBits, class CombineOp>
1419  const polyMesh& mesh,
1420  PackedList<nBits>& faceValues,
1421  const CombineOp& cop,
1422  const bool parRun
1423 )
1424 {
1425  if (faceValues.size() != mesh.nFaces())
1426  {
1428  << "Number of values " << faceValues.size()
1429  << " is not equal to the number of faces in the mesh "
1430  << mesh.nFaces() << abort(FatalError);
1431  }
1432 
1433  const polyBoundaryMesh& patches = mesh.boundaryMesh();
1434 
1435  if (parRun)
1436  {
1438 
1439  // Send
1440 
1441  forAll(patches, patchi)
1442  {
1443  if
1444  (
1445  isA<processorPolyPatch>(patches[patchi])
1446  && patches[patchi].size() > 0
1447  )
1448  {
1449  const processorPolyPatch& procPatch =
1450  refCast<const processorPolyPatch>(patches[patchi]);
1451 
1452  List<unsigned int> patchInfo(procPatch.size());
1453  forAll(procPatch, i)
1454  {
1455  patchInfo[i] = faceValues[procPatch.start()+i];
1456  }
1457 
1458  UOPstream toNbr(procPatch.neighbProcNo(), pBufs);
1459  toNbr << patchInfo;
1460  }
1461  }
1462 
1463 
1464  pBufs.finishedSends();
1465 
1466  // Receive and combine.
1467 
1468  forAll(patches, patchi)
1469  {
1470  if
1471  (
1472  isA<processorPolyPatch>(patches[patchi])
1473  && patches[patchi].size() > 0
1474  )
1475  {
1476  const processorPolyPatch& procPatch =
1477  refCast<const processorPolyPatch>(patches[patchi]);
1478 
1479  List<unsigned int> patchInfo(procPatch.size());
1480  {
1481  UIPstream fromNbr(procPatch.neighbProcNo(), pBufs);
1482  fromNbr >> patchInfo;
1483  }
1484 
1485  // Combine (bitwise)
1486  forAll(procPatch, i)
1487  {
1488  unsigned int patchVal = patchInfo[i];
1489  label meshFacei = procPatch.start()+i;
1490  unsigned int faceVal = faceValues[meshFacei];
1491  cop(faceVal, patchVal);
1492  faceValues[meshFacei] = faceVal;
1493  }
1494  }
1495  }
1496  }
1497 
1498  // Do the cyclics.
1499  forAll(patches, patchi)
1500  {
1501  if (isA<cyclicPolyPatch>(patches[patchi]))
1502  {
1503  const cyclicPolyPatch& cycPatch =
1504  refCast<const cyclicPolyPatch>(patches[patchi]);
1505 
1506  if (cycPatch.owner())
1507  {
1508  // Owner does all.
1509  const cyclicPolyPatch& nbrPatch = cycPatch.neighbPatch();
1510 
1511  for (label i = 0; i < cycPatch.size(); i++)
1512  {
1513  label meshFace0 = cycPatch.start()+i;
1514  unsigned int val0 = faceValues[meshFace0];
1515  label meshFace1 = nbrPatch.start()+i;
1516  unsigned int val1 = faceValues[meshFace1];
1517 
1518  unsigned int t = val0;
1519  cop(t, val1);
1520  faceValues[meshFace0] = t;
1521 
1522  cop(val1, val0);
1523  faceValues[meshFace1] = val1;
1524  }
1525  }
1526  }
1527  }
1528 }
1529 
1530 
1531 template<class T>
1534  const polyMesh& mesh,
1535  const UList<T>& cellData,
1536  List<T>& neighbourCellData
1537 )
1538 {
1539  if (cellData.size() != mesh.nCells())
1540  {
1542  << "Number of cell values " << cellData.size()
1543  << " is not equal to the number of cells in the mesh "
1544  << mesh.nCells() << abort(FatalError);
1545  }
1546 
1547  const polyBoundaryMesh& patches = mesh.boundaryMesh();
1548 
1549  label nBnd = mesh.nFaces()-mesh.nInternalFaces();
1550 
1551  neighbourCellData.setSize(nBnd);
1552 
1553  forAll(patches, patchi)
1554  {
1555  const polyPatch& pp = patches[patchi];
1556  const labelUList& faceCells = pp.faceCells();
1557  forAll(faceCells, i)
1558  {
1559  label bFacei = pp.start()+i-mesh.nInternalFaces();
1560  neighbourCellData[bFacei] = cellData[faceCells[i]];
1561  }
1562  }
1563  syncTools::swapBoundaryFaceList(mesh, neighbourCellData);
1564 }
1565 
1566 
1567 template<unsigned nBits>
1570  const polyMesh& mesh,
1571  PackedList<nBits>& faceValues
1572 )
1573 {
1574  syncFaceList(mesh, faceValues, eqOp<unsigned int>());
1575 }
1576 
1577 
1578 template<unsigned nBits, class CombineOp>
1581  const polyMesh& mesh,
1582  PackedList<nBits>& pointValues,
1583  const CombineOp& cop,
1584  const unsigned int nullValue
1585 )
1586 {
1587  if (pointValues.size() != mesh.nPoints())
1588  {
1590  << "Number of values " << pointValues.size()
1591  << " is not equal to the number of points in the mesh "
1592  << mesh.nPoints() << abort(FatalError);
1593  }
1594 
1595  const globalMeshData& gd = mesh.globalData();
1596  const labelList& meshPoints = gd.coupledPatch().meshPoints();
1597 
1599  forAll(meshPoints, i)
1600  {
1601  cppFld[i] = pointValues[meshPoints[i]];
1602  }
1603 
1605  (
1606  cppFld,
1607  gd.globalPointSlaves(),
1609  gd.globalPointSlavesMap(),
1610  cop
1611  );
1612 
1613  // Extract back to mesh
1614  forAll(meshPoints, i)
1615  {
1616  pointValues[meshPoints[i]] = cppFld[i];
1617  }
1618 }
1619 
1620 
1621 template<unsigned nBits, class CombineOp>
1624  const polyMesh& mesh,
1625  PackedList<nBits>& edgeValues,
1626  const CombineOp& cop,
1627  const unsigned int nullValue
1628 )
1629 {
1630  if (edgeValues.size() != mesh.nEdges())
1631  {
1633  << "Number of values " << edgeValues.size()
1634  << " is not equal to the number of edges in the mesh "
1635  << mesh.nEdges() << abort(FatalError);
1636  }
1637 
1638  const globalMeshData& gd = mesh.globalData();
1639  const labelList& meshEdges = gd.coupledPatchMeshEdges();
1640 
1642  forAll(meshEdges, i)
1643  {
1644  cppFld[i] = edgeValues[meshEdges[i]];
1645  }
1646 
1648  (
1649  cppFld,
1650  gd.globalEdgeSlaves(),
1652  gd.globalEdgeSlavesMap(),
1653  cop
1654  );
1655 
1656  // Extract back to mesh
1657  forAll(meshEdges, i)
1658  {
1659  edgeValues[meshEdges[i]] = cppFld[i];
1660  }
1661 }
1662 
1663 
1664 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
static void syncBoundaryFaceList(const polyMesh &, UList< T > &, const CombineOp &cop, const TransformOp &top, const bool parRun=Pstream::parRun())
Synchronize values on boundary faces only.
const labelList & sharedPointLabels() const
Return indices of local points that are globally shared.
const Map< label > & coupledPatchMeshEdgeMap() const
Return map from mesh edges to coupledPatch edges.
label nPoints() const
Return number of points supporting patch faces.
const labelListList & globalPointSlaves() const
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
const labelList & neighbPoints() const
Return neighbour point labels. WIP.
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
static int masterNo()
Process index of the master.
Definition: UPstream.H:417
const labelListList & globalPointTransformedSlaves() const
void finishedSends(const bool block=true)
Mark all sends as having been done. This will start receives.
An STL-conforming const_iterator.
Definition: HashTable.H:481
const edgeList & coupledPoints() const
Return connected points (from patch local to neighbour patch local)
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:112
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
label nInternalFaces() const
static int firstSlave()
Process index of first slave.
Definition: UPstream.H:446
const mapDistribute & globalEdgeSlavesMap() const
HashTable< T, label, Hash< label > >::iterator iterator
Definition: Map.H:56
label nFaces() const
Various mesh related information for a parallel run. Upon construction, constructs all info using par...
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
label fcIndex(const label i) const
Return the forward circular index, i.e. the next index.
Definition: UListI.H:58
label nCells() const
Template function to specify if the data of a type are contiguous.
const mapDistribute & globalPointSlavesMap() const
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
patches[0]
const labelListList & globalEdgeSlaves() const
const labelList & coupledPatchMeshEdges() const
Return map from coupledPatch edges to mesh edges.
Pre-declare related SubField type.
Definition: Field.H:60
const edgeList & coupledEdges() const
Return connected edges (from patch local to neighbour patch local).
const labelList & meshPoints() const
Return labelList of mesh points in patch. They are constructed.
Neighbour processor patch.
Input inter-processor communications stream.
Definition: IPstream.H:50
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
Definition: HashTable.C:142
Input inter-processor communications stream operating on external buffer.
Definition: UIPstream.H:53
A list of faces which address into the list of points.
virtual bool owner() const
Does this side own the patch ?
label nGlobalPoints() const
Return number of globally shared points.
const Map< label > & meshPointMap() const
Mesh point map. Given the global point index find its.
const globalIndexAndTransform & globalTransforms() const
Global transforms numbering.
Definition: ops.H:70
A dynamically allocatable list of packed unsigned integers.
Definition: PackedList.H:117
const labelUList & faceCells() const
Return face-cell addressing.
Definition: polyPatch.C:315
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
static void syncEdgeList(const polyMesh &, List< T > &, const CombineOp &cop, const T &nullValue, const TransformOp &top)
Synchronize values on all mesh edges.
label nPoints
Cyclic plane patch.
3D tensor transformation operations.
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1394
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
static void swapFaceList(const polyMesh &mesh, UList< T > &l)
Swap coupled face values.
Definition: syncTools.H:463
static void syncPointMap(const polyMesh &, Map< T > &pointValues, const CombineOp &cop, const TransformOp &top)
Synchronize values on selected points.
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1156
int neighbProcNo() const
Return neighbour processor number.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
errorManip< error > abort(error &err)
Definition: errorManip.H:131
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
Foam::polyBoundaryMesh.
Output inter-processor communications stream operating on external buffer.
Definition: UOPstream.H:54
label nEdges() const
Return number of edges in patch.
Map from edge (expressed as its endpoints) to value.
Definition: EdgeMap.H:47
labelList f(nPoints)
Output inter-processor communications stream.
Definition: OPstream.H:50
Buffers for inter-processor communications streams (UOPstream, UIPstream).
static void syncEdgeMap(const polyMesh &, EdgeMap< T > &edgeValues, const CombineOp &cop, const TransformOp &top)
Synchronize values on selected edges.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const volScalarField & T
Spatial transformation functions for primitive fields.
label nEdges() const
static void swapBoundaryCellList(const polyMesh &mesh, const UList< T > &cellData, List< T > &neighbourCellData)
Swap to obtain neighbour cell values for all boundary faces.
const indirectPrimitivePatch & coupledPatch() const
Return patch of all coupled faces.
static void syncData(List< Type > &pointData, const labelListList &slaves, const labelListList &transformedSlaves, const mapDistribute &slavesMap, const globalIndexAndTransform &, const CombineOp &cop, const TransformOp &top)
Helper: synchronise data with transforms.
label constructSize() const
Constructed data size.
static void syncFaceList(const polyMesh &mesh, UList< T > &l, const CombineOp &cop)
Synchronize values on all mesh faces.
Definition: syncTools.H:381
void setSize(const label)
Reset size of List.
Definition: List.C:281
bool set(const label &, const T &newElmt)
Assign a new hashedEntry, overwriting existing entries.
Definition: HashTableI.H:91
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
label patchi
Class containing processor-to-processor mapping information.
A List with indirect addressing.
Definition: fvMatrix.H:106
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:303
static void syncPointList(const polyMesh &, List< T > &, const CombineOp &cop, const T &nullValue, const TransformOp &top)
Synchronize values on all mesh points.
label nPoints() const
void resize(const label newSize)
Resize the hash table for efficiency.
Definition: HashTable.C:432
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:98
label size() const
Number of entries.
Definition: PackedListI.H:711
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
static void swapBoundaryFaceList(const polyMesh &mesh, UList< T > &l)
Swap coupled boundary face values.
Definition: syncTools.H:430
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
const cyclicPolyPatch & neighbPatch() const
label size() const
Return the number of elements in the UList.
Definition: UListI.H:299
friend class iterator
Declare friendship with the iterator.
Definition: HashTable.H:194
const labelList & sharedPointAddr() const
Return addressing into the complete globally shared points.
Determination and storage of the possible independent transforms introduced by coupledPolyPatches, as well as all of the possible permutations of these transforms generated by the presence of multiple coupledPolyPatches, i.e. more than one cyclic boundary. Note that any given point can be on maximum 3 transforms only (and these transforms have to be perpendicular)
static int lastSlave(const label communicator=0)
Process index of last slave.
Definition: UPstream.H:452
A HashTable to objects of type <T> with a label key.
Definition: Map.H:49
const labelListList & globalEdgeTransformedSlaves() const
void syncPointData(List< Type > &pointData, const CombineOp &cop, const TransformOp &top) const
Helper to synchronise coupled patch point data.