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