meshToMeshParallelOps.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) 2012-2017 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 "meshToMesh.H"
27 #include "OFstream.H"
28 #include "Time.H"
29 #include "globalIndex.H"
30 #include "mergePoints.H"
31 #include "processorPolyPatch.H"
32 #include "SubField.H"
33 
34 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
35 
36 Foam::label Foam::meshToMesh::calcDistribution
37 (
38  const polyMesh& src,
39  const polyMesh& tgt
40 ) const
41 {
42  label proci = 0;
43 
44  if (Pstream::parRun())
45  {
46  List<label> cellsPresentOnProc(Pstream::nProcs(), 0);
47  if ((src.nCells() > 0) || (tgt.nCells() > 0))
48  {
49  cellsPresentOnProc[Pstream::myProcNo()] = 1;
50  }
51  else
52  {
53  cellsPresentOnProc[Pstream::myProcNo()] = 0;
54  }
55 
56  Pstream::gatherList(cellsPresentOnProc);
57  Pstream::scatterList(cellsPresentOnProc);
58 
59  label nHaveCells = sum(cellsPresentOnProc);
60 
61  if (nHaveCells > 1)
62  {
63  proci = -1;
64  if (debug)
65  {
66  InfoInFunction
67  << "Meshes split across multiple processors" << endl;
68  }
69  }
70  else if (nHaveCells == 1)
71  {
72  proci = findIndex(cellsPresentOnProc, 1);
73  if (debug)
74  {
75  InfoInFunction
76  << "Meshes local to processor" << proci << endl;
77  }
78  }
79  }
80 
81  return proci;
82 }
83 
84 
85 Foam::label Foam::meshToMesh::calcOverlappingProcs
86 (
87  const List<boundBox>& procBb,
88  const boundBox& bb,
89  boolList& overlaps
90 ) const
91 {
92  overlaps = false;
93 
94  label nOverlaps = 0;
95 
96  forAll(procBb, proci)
97  {
98  const boundBox& bbp = procBb[proci];
99 
100  if (bbp.overlaps(bb))
101  {
102  overlaps[proci] = true;
103  nOverlaps++;
104  }
105  }
106 
107  return nOverlaps;
108 }
109 
110 
111 Foam::autoPtr<Foam::mapDistribute> Foam::meshToMesh::calcProcMap
112 (
113  const polyMesh& src,
114  const polyMesh& tgt
115 ) const
116 {
117  // get decomposition of cells on src mesh
118  List<boundBox> procBb(Pstream::nProcs());
119 
120  if (src.nCells() > 0)
121  {
122  // bounding box for my mesh - do not parallel reduce
123  procBb[Pstream::myProcNo()] = boundBox(src.points(), false);
124 
125  // slightly increase size of bounding boxes to allow for cases where
126  // bounding boxes are perfectly alligned
127  procBb[Pstream::myProcNo()].inflate(0.01);
128  }
129  else
130  {
131  procBb[Pstream::myProcNo()] = boundBox();
132  }
133 
134 
135  Pstream::gatherList(procBb);
136  Pstream::scatterList(procBb);
137 
138 
139  if (debug)
140  {
141  InfoInFunction
142  << "Determining extent of src mesh per processor:" << nl
143  << "\tproc\tbb" << endl;
144  forAll(procBb, proci)
145  {
146  Info<< '\t' << proci << '\t' << procBb[proci] << endl;
147  }
148  }
149 
150 
151  // determine which cells of tgt mesh overlaps src mesh per proc
152  const cellList& cells = tgt.cells();
153  const faceList& faces = tgt.faces();
154  const pointField& points = tgt.points();
155 
156  labelListList sendMap;
157 
158  {
159  // per processor indices into all segments to send
160  List<DynamicList<label>> dynSendMap(Pstream::nProcs());
161  label iniSize = floor(tgt.nCells()/Pstream::nProcs());
162 
163  forAll(dynSendMap, proci)
164  {
165  dynSendMap[proci].setCapacity(iniSize);
166  }
167 
168  // work array - whether src processor bb overlaps the tgt cell bounds
169  boolList procBbOverlaps(Pstream::nProcs());
170  forAll(cells, celli)
171  {
172  const cell& c = cells[celli];
173 
174  // determine bounding box of tgt cell
175  boundBox cellBb(point::max, point::min);
176  forAll(c, facei)
177  {
178  const face& f = faces[c[facei]];
179  forAll(f, fp)
180  {
181  cellBb.min() = min(cellBb.min(), points[f[fp]]);
182  cellBb.max() = max(cellBb.max(), points[f[fp]]);
183  }
184  }
185 
186  // find the overlapping tgt cells on each src processor
187  (void)calcOverlappingProcs(procBb, cellBb, procBbOverlaps);
188 
189  forAll(procBbOverlaps, proci)
190  {
191  if (procBbOverlaps[proci])
192  {
193  dynSendMap[proci].append(celli);
194  }
195  }
196  }
197 
198  // convert dynamicList to labelList
199  sendMap.setSize(Pstream::nProcs());
200  forAll(sendMap, proci)
201  {
202  sendMap[proci].transfer(dynSendMap[proci]);
203  }
204  }
205 
206  // debug printing
207  if (debug)
208  {
209  Pout<< "Of my " << cells.size() << " target cells I need to send to:"
210  << nl << "\tproc\tcells" << endl;
211  forAll(sendMap, proci)
212  {
213  Pout<< '\t' << proci << '\t' << sendMap[proci].size() << endl;
214  }
215  }
216 
217 
218  // send over how many tgt cells I need to receive from each processor
219  labelListList sendSizes(Pstream::nProcs());
220  sendSizes[Pstream::myProcNo()].setSize(Pstream::nProcs());
221  forAll(sendMap, proci)
222  {
223  sendSizes[Pstream::myProcNo()][proci] = sendMap[proci].size();
224  }
225  Pstream::gatherList(sendSizes);
226  Pstream::scatterList(sendSizes);
227 
228 
229  // determine order of receiving
230  labelListList constructMap(Pstream::nProcs());
231 
232  label segmentI = 0;
233  forAll(constructMap, proci)
234  {
235  // what I need to receive is what other processor is sending to me
236  label nRecv = sendSizes[proci][Pstream::myProcNo()];
237  constructMap[proci].setSize(nRecv);
238 
239  for (label i = 0; i < nRecv; i++)
240  {
241  constructMap[proci][i] = segmentI++;
242  }
243  }
244 
245  autoPtr<mapDistribute> mapPtr
246  (
247  new mapDistribute
248  (
249  segmentI, // size after construction
250  sendMap.xfer(),
251  constructMap.xfer()
252  )
253  );
254 
255  return mapPtr;
256 }
257 
258 
259 void Foam::meshToMesh::distributeCells
260 (
261  const mapDistribute& map,
262  const polyMesh& tgtMesh,
263  const globalIndex& globalI,
264  List<pointField>& points,
265  List<label>& nInternalFaces,
266  List<faceList>& faces,
267  List<labelList>& faceOwner,
268  List<labelList>& faceNeighbour,
269  List<labelList>& cellIDs,
270  List<labelList>& nbrProcIDs,
271  List<labelList>& procLocalFaceIDs
272 ) const
273 {
274  PstreamBuffers pBufs(Pstream::commsTypes::nonBlocking);
275 
276  points.setSize(Pstream::nProcs());
277  nInternalFaces.setSize(Pstream::nProcs(), 0);
278  faces.setSize(Pstream::nProcs());
279  faceOwner.setSize(Pstream::nProcs());
280  faceNeighbour.setSize(Pstream::nProcs());
281  cellIDs.setSize(Pstream::nProcs());
282 
283  nbrProcIDs.setSize(Pstream::nProcs());;
284  procLocalFaceIDs.setSize(Pstream::nProcs());;
285 
286 
287  for (label domain = 0; domain < Pstream::nProcs(); domain++)
288  {
289  const labelList& sendElems = map.subMap()[domain];
290 
291  if (sendElems.size())
292  {
293  // reverse cell map
294  labelList reverseCellMap(tgtMesh.nCells(), -1);
295  forAll(sendElems, subCelli)
296  {
297  reverseCellMap[sendElems[subCelli]] = subCelli;
298  }
299 
300  DynamicList<face> subFaces(tgtMesh.nFaces());
301  DynamicList<label> subFaceOwner(tgtMesh.nFaces());
302  DynamicList<label> subFaceNeighbour(tgtMesh.nFaces());
303 
304  DynamicList<label> subNbrProcIDs(tgtMesh.nFaces());
305  DynamicList<label> subProcLocalFaceIDs(tgtMesh.nFaces());
306 
307  label nInternal = 0;
308 
309  // internal faces
310  forAll(tgtMesh.faceNeighbour(), facei)
311  {
312  label own = tgtMesh.faceOwner()[facei];
313  label nbr = tgtMesh.faceNeighbour()[facei];
314  label subOwn = reverseCellMap[own];
315  label subNbr = reverseCellMap[nbr];
316 
317  if (subOwn != -1 && subNbr != -1)
318  {
319  nInternal++;
320 
321  if (subOwn < subNbr)
322  {
323  subFaces.append(tgtMesh.faces()[facei]);
324  subFaceOwner.append(subOwn);
325  subFaceNeighbour.append(subNbr);
326  subNbrProcIDs.append(-1);
327  subProcLocalFaceIDs.append(-1);
328  }
329  else
330  {
331  subFaces.append(tgtMesh.faces()[facei].reverseFace());
332  subFaceOwner.append(subNbr);
333  subFaceNeighbour.append(subOwn);
334  subNbrProcIDs.append(-1);
335  subProcLocalFaceIDs.append(-1);
336  }
337  }
338  }
339 
340  // boundary faces for new region
341  forAll(tgtMesh.faceNeighbour(), facei)
342  {
343  label own = tgtMesh.faceOwner()[facei];
344  label nbr = tgtMesh.faceNeighbour()[facei];
345  label subOwn = reverseCellMap[own];
346  label subNbr = reverseCellMap[nbr];
347 
348  if (subOwn != -1 && subNbr == -1)
349  {
350  subFaces.append(tgtMesh.faces()[facei]);
351  subFaceOwner.append(subOwn);
352  subFaceNeighbour.append(subNbr);
353  subNbrProcIDs.append(-1);
354  subProcLocalFaceIDs.append(-1);
355  }
356  else if (subOwn == -1 && subNbr != -1)
357  {
358  subFaces.append(tgtMesh.faces()[facei].reverseFace());
359  subFaceOwner.append(subNbr);
360  subFaceNeighbour.append(subOwn);
361  subNbrProcIDs.append(-1);
362  subProcLocalFaceIDs.append(-1);
363  }
364  }
365 
366  // boundary faces of existing region
367  forAll(tgtMesh.boundaryMesh(), patchi)
368  {
369  const polyPatch& pp = tgtMesh.boundaryMesh()[patchi];
370 
371  label nbrProci = -1;
372 
373  // store info for faces on processor patches
374  if (isA<processorPolyPatch>(pp))
375  {
376  const processorPolyPatch& ppp =
377  dynamic_cast<const processorPolyPatch&>(pp);
378 
379  nbrProci = ppp.neighbProcNo();
380  }
381 
382  forAll(pp, i)
383  {
384  label facei = pp.start() + i;
385  label own = tgtMesh.faceOwner()[facei];
386 
387  if (reverseCellMap[own] != -1)
388  {
389  subFaces.append(tgtMesh.faces()[facei]);
390  subFaceOwner.append(reverseCellMap[own]);
391  subFaceNeighbour.append(-1);
392  subNbrProcIDs.append(nbrProci);
393  subProcLocalFaceIDs.append(i);
394  }
395  }
396  }
397 
398  // reverse point map
399  labelList reversePointMap(tgtMesh.nPoints(), -1);
400  DynamicList<point> subPoints(tgtMesh.nPoints());
401  forAll(subFaces, subFacei)
402  {
403  face& f = subFaces[subFacei];
404  forAll(f, fp)
405  {
406  label pointi = f[fp];
407  if (reversePointMap[pointi] == -1)
408  {
409  reversePointMap[pointi] = subPoints.size();
410  subPoints.append(tgtMesh.points()[pointi]);
411  }
412 
413  f[fp] = reversePointMap[pointi];
414  }
415  }
416 
417  // tgt cells into global numbering
418  labelList globalElems(sendElems.size());
419  forAll(sendElems, i)
420  {
421  if (debug)
422  {
423  Pout<< "tgtProc:" << Pstream::myProcNo()
424  << " sending tgt cell " << sendElems[i]
425  << "[" << globalI.toGlobal(sendElems[i]) << "]"
426  << " to srcProc " << domain << endl;
427  }
428 
429  globalElems[i] = globalI.toGlobal(sendElems[i]);
430  }
431 
432  // pass data
433  if (domain == Pstream::myProcNo())
434  {
435  // allocate my own data
436  points[Pstream::myProcNo()] = subPoints;
437  nInternalFaces[Pstream::myProcNo()] = nInternal;
438  faces[Pstream::myProcNo()] = subFaces;
439  faceOwner[Pstream::myProcNo()] = subFaceOwner;
440  faceNeighbour[Pstream::myProcNo()] = subFaceNeighbour;
441  cellIDs[Pstream::myProcNo()] = globalElems;
442  nbrProcIDs[Pstream::myProcNo()] = subNbrProcIDs;
443  procLocalFaceIDs[Pstream::myProcNo()] = subProcLocalFaceIDs;
444  }
445  else
446  {
447  // send data to other processor domains
448  UOPstream toDomain(domain, pBufs);
449 
450  toDomain
451  << subPoints
452  << nInternal
453  << subFaces
454  << subFaceOwner
455  << subFaceNeighbour
456  << globalElems
457  << subNbrProcIDs
458  << subProcLocalFaceIDs;
459  }
460  }
461  }
462 
463  // Start receiving
464  pBufs.finishedSends();
465 
466  // Consume
467  for (label domain = 0; domain < Pstream::nProcs(); domain++)
468  {
469  const labelList& recvElems = map.constructMap()[domain];
470 
471  if (domain != Pstream::myProcNo() && recvElems.size())
472  {
473  UIPstream str(domain, pBufs);
474 
475  str >> points[domain]
476  >> nInternalFaces[domain]
477  >> faces[domain]
478  >> faceOwner[domain]
479  >> faceNeighbour[domain]
480  >> cellIDs[domain]
481  >> nbrProcIDs[domain]
482  >> procLocalFaceIDs[domain];
483  }
484 
485  if (debug)
486  {
487  Pout<< "Target mesh send sizes[" << domain << "]"
488  << ": points="<< points[domain].size()
489  << ", faces=" << faces[domain].size()
490  << ", nInternalFaces=" << nInternalFaces[domain]
491  << ", faceOwn=" << faceOwner[domain].size()
492  << ", faceNbr=" << faceNeighbour[domain].size()
493  << ", cellIDs=" << cellIDs[domain].size() << endl;
494  }
495  }
496 }
497 
498 
499 void Foam::meshToMesh::distributeAndMergeCells
500 (
501  const mapDistribute& map,
502  const polyMesh& tgt,
503  const globalIndex& globalI,
504  pointField& tgtPoints,
505  faceList& tgtFaces,
506  labelList& tgtFaceOwners,
507  labelList& tgtFaceNeighbours,
508  labelList& tgtCellIDs
509 ) const
510 {
511  // Exchange per-processor data
512  List<pointField> allPoints;
513  List<label> allNInternalFaces;
514  List<faceList> allFaces;
515  List<labelList> allFaceOwners;
516  List<labelList> allFaceNeighbours;
517  List<labelList> allTgtCellIDs;
518 
519  // Per target mesh face the neighbouring proc and index in
520  // processor patch (all -1 for normal boundary face)
521  List<labelList> allNbrProcIDs;
522  List<labelList> allProcLocalFaceIDs;
523 
524  distributeCells
525  (
526  map,
527  tgt,
528  globalI,
529  allPoints,
530  allNInternalFaces,
531  allFaces,
532  allFaceOwners,
533  allFaceNeighbours,
534  allTgtCellIDs,
535  allNbrProcIDs,
536  allProcLocalFaceIDs
537  );
538 
539  // Convert lists into format that can be used to generate a valid polyMesh
540  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
541  //
542  // Points and cells are collected into single flat lists:
543  // - i.e. proc0, proc1 ... procN
544  //
545  // Faces need to be sorted after collection to that internal faces are
546  // contiguous, followed by all boundary faces
547  //
548  // Processor patch faces between included cells on neighbouring processors
549  // are converted into internal faces
550  //
551  // Face list structure:
552  // - Per processor:
553  // - internal faces
554  // - processor faces that have been converted into internal faces
555  // - Followed by all boundary faces
556  // - from 'normal' boundary faces
557  // - from singularly-sided processor patch faces
558 
559 
560  // Number of internal+coupled faces
561  labelList allNIntCoupledFaces(allNInternalFaces);
562 
563  // Starting offset for points
564  label nPoints = 0;
565  labelList pointOffset(Pstream::nProcs(), 0);
566  forAll(allPoints, proci)
567  {
568  pointOffset[proci] = nPoints;
569  nPoints += allPoints[proci].size();
570  }
571 
572  // Starting offset for cells
573  label nCells = 0;
574  labelList cellOffset(Pstream::nProcs(), 0);
575  forAll(allTgtCellIDs, proci)
576  {
577  cellOffset[proci] = nCells;
578  nCells += allTgtCellIDs[proci].size();
579  }
580 
581  // Count any coupled faces
582  typedef FixedList<label, 3> label3;
583  typedef HashTable<label, label3, label3::Hash<>> procCoupleInfo;
584  procCoupleInfo procFaceToGlobalCell;
585 
586  forAll(allNbrProcIDs, proci)
587  {
588  const labelList& nbrProci = allNbrProcIDs[proci];
589  const labelList& localFacei = allProcLocalFaceIDs[proci];
590 
591  forAll(nbrProci, i)
592  {
593  if (nbrProci[i] != -1 && localFacei[i] != -1)
594  {
595  label3 key;
596  key[0] = min(proci, nbrProci[i]);
597  key[1] = max(proci, nbrProci[i]);
598  key[2] = localFacei[i];
599 
600  procCoupleInfo::const_iterator fnd =
601  procFaceToGlobalCell.find(key);
602 
603  if (fnd == procFaceToGlobalCell.end())
604  {
605  procFaceToGlobalCell.insert(key, -1);
606  }
607  else
608  {
609  if (debug)
610  {
611  Pout<< "Additional internal face between procs:"
612  << key[0] << " and " << key[1]
613  << " across local face " << key[2] << endl;
614  }
615 
616  allNIntCoupledFaces[proci]++;
617  }
618  }
619  }
620  }
621 
622 
623  // Starting offset for internal faces
624  label nIntFaces = 0;
625  label nFacesTotal = 0;
626  labelList internalFaceOffset(Pstream::nProcs(), 0);
627  forAll(allNIntCoupledFaces, proci)
628  {
629  label nCoupledFaces =
630  allNIntCoupledFaces[proci] - allNInternalFaces[proci];
631 
632  internalFaceOffset[proci] = nIntFaces;
633  nIntFaces += allNIntCoupledFaces[proci];
634  nFacesTotal += allFaceOwners[proci].size() - nCoupledFaces;
635  }
636 
637  tgtPoints.setSize(nPoints);
638  tgtFaces.setSize(nFacesTotal);
639  tgtFaceOwners.setSize(nFacesTotal);
640  tgtFaceNeighbours.setSize(nFacesTotal);
641  tgtCellIDs.setSize(nCells);
642 
643  // Insert points
644  forAll(allPoints, proci)
645  {
646  const pointField& pts = allPoints[proci];
647  SubList<point>(tgtPoints, pts.size(), pointOffset[proci]) = pts;
648  }
649 
650  // Insert cellIDs
651  forAll(allTgtCellIDs, proci)
652  {
653  const labelList& cellIDs = allTgtCellIDs[proci];
654  SubList<label>(tgtCellIDs, cellIDs.size(), cellOffset[proci]) = cellIDs;
655  }
656 
657 
658  // Insert internal faces (from internal faces)
659  forAll(allFaces, proci)
660  {
661  const faceList& fcs = allFaces[proci];
662  const labelList& faceOs = allFaceOwners[proci];
663  const labelList& faceNs = allFaceNeighbours[proci];
664 
665  SubList<face> slice
666  (
667  tgtFaces,
668  allNInternalFaces[proci],
669  internalFaceOffset[proci]
670  );
671  slice = SubList<face>(fcs, allNInternalFaces[proci]);
672  forAll(slice, i)
673  {
674  add(slice[i], pointOffset[proci]);
675  }
676 
677  SubField<label> ownSlice
678  (
679  tgtFaceOwners,
680  allNInternalFaces[proci],
681  internalFaceOffset[proci]
682  );
683  ownSlice = SubField<label>(faceOs, allNInternalFaces[proci]);
684  add(ownSlice, cellOffset[proci]);
685 
686  SubField<label> nbrSlice
687  (
688  tgtFaceNeighbours,
689  allNInternalFaces[proci],
690  internalFaceOffset[proci]
691  );
692  nbrSlice = SubField<label>(faceNs, allNInternalFaces[proci]);
693  add(nbrSlice, cellOffset[proci]);
694 
695  internalFaceOffset[proci] += allNInternalFaces[proci];
696  }
697 
698 
699  // Insert internal faces (from coupled face-pairs)
700  forAll(allNbrProcIDs, proci)
701  {
702  const labelList& nbrProci = allNbrProcIDs[proci];
703  const labelList& localFacei = allProcLocalFaceIDs[proci];
704  const labelList& faceOs = allFaceOwners[proci];
705  const faceList& fcs = allFaces[proci];
706 
707  forAll(nbrProci, i)
708  {
709  if (nbrProci[i] != -1 && localFacei[i] != -1)
710  {
711  label3 key;
712  key[0] = min(proci, nbrProci[i]);
713  key[1] = max(proci, nbrProci[i]);
714  key[2] = localFacei[i];
715 
716  procCoupleInfo::iterator fnd = procFaceToGlobalCell.find(key);
717 
718  if (fnd != procFaceToGlobalCell.end())
719  {
720  label tgtFacei = fnd();
721  if (tgtFacei == -1)
722  {
723  // on first visit store the new cell on this side
724  fnd() = cellOffset[proci] + faceOs[i];
725  }
726  else
727  {
728  // get owner and neighbour in new cell numbering
729  label newOwn = cellOffset[proci] + faceOs[i];
730  label newNbr = fnd();
731  label tgtFacei = internalFaceOffset[proci]++;
732 
733  if (debug)
734  {
735  Pout<< " proc " << proci
736  << "\tinserting face:" << tgtFacei
737  << " connection between owner " << newOwn
738  << " and neighbour " << newNbr
739  << endl;
740  }
741 
742  if (newOwn < newNbr)
743  {
744  // we have correct orientation
745  tgtFaces[tgtFacei] = fcs[i];
746  tgtFaceOwners[tgtFacei] = newOwn;
747  tgtFaceNeighbours[tgtFacei] = newNbr;
748  }
749  else
750  {
751  // reverse orientation
752  tgtFaces[tgtFacei] = fcs[i].reverseFace();
753  tgtFaceOwners[tgtFacei] = newNbr;
754  tgtFaceNeighbours[tgtFacei] = newOwn;
755  }
756 
757  add(tgtFaces[tgtFacei], pointOffset[proci]);
758 
759  // mark with unique value
760  fnd() = -2;
761  }
762  }
763  }
764  }
765  }
766 
767 
768  forAll(allNbrProcIDs, proci)
769  {
770  const labelList& nbrProci = allNbrProcIDs[proci];
771  const labelList& localFacei = allProcLocalFaceIDs[proci];
772  const labelList& faceOs = allFaceOwners[proci];
773  const labelList& faceNs = allFaceNeighbours[proci];
774  const faceList& fcs = allFaces[proci];
775 
776  forAll(nbrProci, i)
777  {
778  // coupled boundary face
779  if (nbrProci[i] != -1 && localFacei[i] != -1)
780  {
781  label3 key;
782  key[0] = min(proci, nbrProci[i]);
783  key[1] = max(proci, nbrProci[i]);
784  key[2] = localFacei[i];
785 
786  label tgtFacei = procFaceToGlobalCell[key];
787 
788  if (tgtFacei == -1)
789  {
790  FatalErrorInFunction
791  << "Unvisited " << key
792  << abort(FatalError);
793  }
794  else if (tgtFacei != -2)
795  {
796  label newOwn = cellOffset[proci] + faceOs[i];
797  label tgtFacei = nIntFaces++;
798 
799  if (debug)
800  {
801  Pout<< " proc " << proci
802  << "\tinserting boundary face:" << tgtFacei
803  << " from coupled face " << key
804  << endl;
805  }
806 
807  tgtFaces[tgtFacei] = fcs[i];
808  add(tgtFaces[tgtFacei], pointOffset[proci]);
809 
810  tgtFaceOwners[tgtFacei] = newOwn;
811  tgtFaceNeighbours[tgtFacei] = -1;
812  }
813  }
814  // normal boundary face
815  else
816  {
817  label own = faceOs[i];
818  label nbr = faceNs[i];
819  if ((own != -1) && (nbr == -1))
820  {
821  label newOwn = cellOffset[proci] + faceOs[i];
822  label tgtFacei = nIntFaces++;
823 
824  tgtFaces[tgtFacei] = fcs[i];
825  add(tgtFaces[tgtFacei], pointOffset[proci]);
826 
827  tgtFaceOwners[tgtFacei] = newOwn;
828  tgtFaceNeighbours[tgtFacei] = -1;
829  }
830  }
831  }
832  }
833 
834 
835  if (debug)
836  {
837  // only merging points in debug mode
838 
839  labelList oldToNew;
840  pointField newTgtPoints;
841  bool hasMerged = mergePoints
842  (
843  tgtPoints,
844  SMALL,
845  false,
846  oldToNew,
847  newTgtPoints
848  );
849 
850  if (hasMerged)
851  {
852  if (debug)
853  {
854  Pout<< "Merged from " << tgtPoints.size()
855  << " down to " << newTgtPoints.size() << " points" << endl;
856  }
857 
858  tgtPoints.transfer(newTgtPoints);
859  forAll(tgtFaces, i)
860  {
861  inplaceRenumber(oldToNew, tgtFaces[i]);
862  }
863  }
864  }
865 }
866 
867 
868 // ************************************************************************* //
Foam::labelListList
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:57
Foam::Pstream::scatterList
static void scatterList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Scatter data. Reverse of gatherList.
Definition: gatherScatterList.C:212
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
Foam::label
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
Foam::inplaceRenumber
void inplaceRenumber(const labelUList &oldToNew, ListType &)
Inplace renumber the values of a list.
Definition: ListOpsTemplates.C:57
Foam::FatalError
error FatalError
Foam::max
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Foam::UPstream::myProcNo
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:418
Foam::faceList
List< face > faceList
Definition: faceListFwd.H:43
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
Foam::sum
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
Definition: DimensionedFieldFunctions.C:333
Foam::boolList
List< bool > boolList
Bool container classes.
Definition: boolList.H:50
Foam::pointField
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
nPoints
label nPoints
Definition: gmvOutputHeader.H:2
Foam::labelList
List< label > labelList
A List of labels.
Definition: labelList.H:56
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Foam::nl
static const char nl
Definition: Ostream.H:262
mergePoints.H
Merge points. See below.
Foam::findIndex
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurence of given element and return index,.
Foam::min
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
Foam::UPstream::parRun
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:394
Foam::UPstream::nProcs
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:400
patchi
label patchi
Definition: getPatchFieldScalar.H:1
Foam::Pout
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
Foam::mergePoints
label mergePoints(const UList< Type > &points, const scalar mergeTol, const bool verbose, labelList &pointMap, const Type &origin=Type::zero)
Sorts and merges points. All points closer than/equal mergeTol get merged.
Foam::Info
messageStream Info
Foam::DelaunayMeshTools::allPoints
tmp< pointField > allPoints(const Triangulation &t)
Extract all points in vertex-index order.
Foam::Pstream::gatherList
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
Definition: gatherScatterList.C:49
Foam::cellList
List< cell > cellList
list of cells
Definition: cellList.H:42
InfoInFunction
#define InfoInFunction
Report an information message using Foam::Info.
Definition: messageStream.H:281
Generated by   doxygen 1.8.13