lduPrimitiveMesh.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) 2013-2014 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 "lduPrimitiveMesh.H"
27 #include "processorLduInterface.H"
28 #include "EdgeMap.H"
29 #include "labelPair.H"
30 #include "processorGAMGInterface.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36  defineTypeNameAndDebug(lduPrimitiveMesh, 0);
37 
38  //- Less operator for pairs of <processor><index>
39  class procLess
40  {
41  const labelPairList& lst_;
42 
43  public:
44 
45  procLess(const labelPairList& lst)
46  :
47  lst_(lst)
48  {}
49 
50  bool operator()(const label a, const label b)
51  {
52  return lst_[a].first() < lst_[b].first();
53  }
54  };
55 
56 }
57 
58 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
59 
60 void Foam::lduPrimitiveMesh::checkUpperTriangular
61 (
62  const label size,
63  const labelUList& l,
64  const labelUList& u
65 )
66 {
67  forAll(l, faceI)
68  {
69  if (u[faceI] < l[faceI])
70  {
72  (
73  "checkUpperTriangular"
74  "(const label, const labelUList&, const labelUList&)"
75  ) << "Reversed face. Problem at face " << faceI
76  << " l:" << l[faceI] << " u:" << u[faceI]
77  << abort(FatalError);
78  }
79  if (l[faceI] < 0 || u[faceI] < 0 || u[faceI] >= size)
80  {
82  (
83  "checkUpperTriangular"
84  "(const label, const labelUList&, const labelUList&)"
85  ) << "Illegal cell label. Problem at face " << faceI
86  << " l:" << l[faceI] << " u:" << u[faceI]
87  << abort(FatalError);
88  }
89  }
90 
91  for (label faceI=1; faceI < l.size(); faceI++)
92  {
93  if (l[faceI-1] > l[faceI])
94  {
96  (
97  "checkUpperTriangular"
98  "(const label, const labelUList&, const labelUList&)"
99  ) << "Lower not in incremental cell order."
100  << " Problem at face " << faceI
101  << " l:" << l[faceI] << " u:" << u[faceI]
102  << " previous l:" << l[faceI-1]
103  << abort(FatalError);
104  }
105  else if (l[faceI-1] == l[faceI])
106  {
107  // Same cell.
108  if (u[faceI-1] > u[faceI])
109  {
111  (
112  "checkUpperTriangular"
113  "(const label, const labelUList&, const labelUList&)"
114  ) << "Upper not in incremental cell order."
115  << " Problem at face " << faceI
116  << " l:" << l[faceI] << " u:" << u[faceI]
117  << " previous u:" << u[faceI-1]
118  << abort(FatalError);
119  }
120  }
121  }
122 }
123 
124 
125 Foam::label Foam::lduPrimitiveMesh::totalSize
126 (
127  const PtrList<lduPrimitiveMesh>& meshes
128 )
129 {
130  label size = 0;
131 
132  forAll(meshes, i)
133  {
134  size += meshes[i].lduAddr().size();
135  }
136  return size;
137 }
138 
139 
140 Foam::labelList Foam::lduPrimitiveMesh::upperTriOrder
141 (
142  const label nCells,
143  const labelUList& lower,
144  const labelUList& upper
145 )
146 {
147  labelList nNbrs(nCells, 0);
148 
149  // Count number of upper neighbours
150  forAll(lower, faceI)
151  {
152  if (upper[faceI] < lower[faceI])
153  {
154  FatalErrorIn("lduPrimitiveMesh::upperTriOrder(..)")
155  << "Problem at face:" << faceI
156  << " lower:" << lower[faceI]
157  << " upper:" << upper[faceI]
158  << exit(FatalError);
159  }
160  nNbrs[lower[faceI]]++;
161  }
162 
163  // Construct cell-upper cell addressing
164  labelList offsets(nCells+1);
165  offsets[0] = 0;
166  forAll(nNbrs, cellI)
167  {
168  offsets[cellI+1] = offsets[cellI]+nNbrs[cellI];
169  }
170 
171  nNbrs = offsets;
172 
173  labelList cellToFaces(offsets.last());
174  forAll(upper, faceI)
175  {
176  label cellI = lower[faceI];
177  cellToFaces[nNbrs[cellI]++] = faceI;
178  }
179 
180  // Sort
181 
182  labelList oldToNew(lower.size());
183 
184  labelList order;
185  labelList nbr;
186 
187  label newFaceI = 0;
188 
189  for (label cellI = 0; cellI < nCells; cellI++)
190  {
191  label startOfCell = offsets[cellI];
192  label nNbr = offsets[cellI+1] - startOfCell;
193 
194  nbr.setSize(nNbr);
195  order.setSize(nNbr);
196  forAll(order, i)
197  {
198  nbr[i] = upper[cellToFaces[offsets[cellI]+i]];
199  }
200  sortedOrder(nbr, order);
201 
202  forAll(order, i)
203  {
204  label index = order[i];
205  oldToNew[cellToFaces[startOfCell + index]] = newFaceI++;
206  }
207  }
208 
209  return oldToNew;
210 }
211 
212 
213 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
214 
215 Foam::lduPrimitiveMesh::lduPrimitiveMesh
216 (
217  const label nCells,
218  labelList& l,
219  labelList& u,
220  const label comm,
221  bool reUse
222 )
223 :
224  lduAddressing(nCells),
225  lowerAddr_(l, reUse),
226  upperAddr_(u, reUse),
227  comm_(comm)
228 {}
229 
230 
232 (
233  lduInterfacePtrsList& interfaces,
234  const lduSchedule& ps
235 )
236 {
237  interfaces_ = interfaces;
238  patchSchedule_ = ps;
239 
240  // Create interfaces
241  primitiveInterfaces_.setSize(interfaces_.size());
242  forAll(interfaces_, i)
243  {
244  if (interfaces_.set(i))
245  {
246  primitiveInterfaces_.set(i, &interfaces_[i]);
247  }
248  }
249 }
250 
251 
252 Foam::lduPrimitiveMesh::lduPrimitiveMesh
253 (
254  const label nCells,
255  labelList& l,
256  labelList& u,
257  PtrList<const lduInterface>& primitiveInterfaces,
258  const lduSchedule& ps,
259  const label comm
260 )
261 :
262  lduAddressing(nCells),
263  lowerAddr_(l, true),
264  upperAddr_(u, true),
265  primitiveInterfaces_(0),
266  patchSchedule_(ps),
267  comm_(comm)
268 {
269  primitiveInterfaces_.transfer(primitiveInterfaces);
270 
271  // Create interfaces
272  interfaces_.setSize(primitiveInterfaces_.size());
273  forAll(primitiveInterfaces_, i)
274  {
275  if (primitiveInterfaces_.set(i))
276  {
277  interfaces_.set(i, &primitiveInterfaces_[i]);
278  }
279  }
280 }
281 
282 
283 Foam::lduPrimitiveMesh::lduPrimitiveMesh
284 (
285  const label comm,
286  const labelList& procAgglomMap,
287 
288  const labelList& procIDs,
289  const lduMesh& myMesh,
290  const PtrList<lduPrimitiveMesh>& otherMeshes,
291 
292  labelList& cellOffsets,
293  labelList& faceOffsets,
294  labelListList& faceMap,
295  labelListList& boundaryMap,
296  labelListListList& boundaryFaceMap
297 )
298 :
299  lduAddressing(myMesh.lduAddr().size() + totalSize(otherMeshes)),
300  lowerAddr_(0),
301  upperAddr_(0),
302  interfaces_(0),
303  patchSchedule_(0),
304  comm_(comm)
305 {
306  const label currentComm = myMesh.comm();
307 
308  forAll(otherMeshes, i)
309  {
310  if (otherMeshes[i].comm() != currentComm)
311  {
312  WarningIn
313  (
314  "lduPrimitiveMesh::lduPrimitiveMesh(..)"
315  ) << "Communicator " << otherMeshes[i].comm()
316  << " at index " << i
317  << " differs from that of predecessor "
318  << currentComm
319  << endl; //exit(FatalError);
320  }
321  }
322 
323  const label nMeshes = otherMeshes.size()+1;
324 
325  const label myAgglom = procAgglomMap[UPstream::myProcNo(currentComm)];
326 
327  if (lduPrimitiveMesh::debug)
328  {
329  Pout<< "I am " << UPstream::myProcNo(currentComm)
330  << " agglomerating into " << myAgglom
331  << " as are " << findIndices(procAgglomMap, myAgglom)
332  << endl;
333  }
334 
335 
336  forAll(procIDs, i)
337  {
338  if (procAgglomMap[procIDs[i]] != procAgglomMap[procIDs[0]])
339  {
341  (
342  "lduPrimitiveMesh::lduPrimitiveMesh(..)"
343  ) << "Processor " << procIDs[i]
344  << " agglomerates to " << procAgglomMap[procIDs[i]]
345  << " whereas other processors " << procIDs
346  << " agglomerate to "
347  << UIndirectList<label>(procAgglomMap, procIDs)
348  << exit(FatalError);
349  }
350  }
351 
352 
353  // Cells get added in order.
354  cellOffsets.setSize(nMeshes+1);
355  cellOffsets[0] = 0;
356  for (label procMeshI = 0; procMeshI < nMeshes; procMeshI++)
357  {
358  const lduMesh& procMesh = mesh(myMesh, otherMeshes, procMeshI);
359 
360  cellOffsets[procMeshI+1] =
361  cellOffsets[procMeshI]
362  + procMesh.lduAddr().size();
363  }
364 
365 
366  // Faces initially get added in order, sorted later
367  labelList internalFaceOffsets(nMeshes+1);
368  internalFaceOffsets[0] = 0;
369  for (label procMeshI = 0; procMeshI < nMeshes; procMeshI++)
370  {
371  const lduMesh& procMesh = mesh(myMesh, otherMeshes, procMeshI);
372 
373  internalFaceOffsets[procMeshI+1] =
374  internalFaceOffsets[procMeshI]
375  + procMesh.lduAddr().lowerAddr().size();
376  }
377 
378  // Count how faces get added. Interfaces inbetween get merged.
379 
380  // Merged interfaces: map from two coarse processors back to
381  // - procMeshes
382  // - interface in procMesh
383  // (estimate size from number of patches of mesh0)
384  EdgeMap<labelPairList> mergedMap(2*myMesh.interfaces().size());
385 
386  // Unmerged interfaces: map from two coarse processors back to
387  // - procMeshes
388  // - interface in procMesh
389  EdgeMap<labelPairList> unmergedMap(mergedMap.size());
390 
391  boundaryMap.setSize(nMeshes);
392  boundaryFaceMap.setSize(nMeshes);
393 
394 
395  label nOtherInterfaces = 0;
396  labelList nCoupledFaces(nMeshes, 0);
397 
398  for (label procMeshI = 0; procMeshI < nMeshes; procMeshI++)
399  {
400  const lduInterfacePtrsList interfaces =
401  mesh(myMesh, otherMeshes, procMeshI).interfaces();
402 
403  // Inialise all boundaries as merged
404  boundaryMap[procMeshI].setSize(interfaces.size(), -1);
405  boundaryFaceMap[procMeshI].setSize(interfaces.size());
406 
407  // Get sorted order of processors
408  forAll(interfaces, intI)
409  {
410  if (interfaces.set(intI))
411  {
412  const lduInterface& ldui = interfaces[intI];
413 
414  if (isA<processorLduInterface>(ldui))
415  {
416  const processorLduInterface& pldui =
417  refCast<const processorLduInterface>(ldui);
418 
419  label agglom0 = procAgglomMap[pldui.myProcNo()];
420  label agglom1 = procAgglomMap[pldui.neighbProcNo()];
421 
422  const edge procEdge(agglom0, agglom1);
423 
424  if (agglom0 != myAgglom && agglom1 != myAgglom)
425  {
426  FatalErrorIn("lduPrimitiveMesh::lduPrimitiveMesh(..)")
427  << "At mesh from processor " << procIDs[procMeshI]
428  << " have interface " << intI
429  << " with myProcNo:" << pldui.myProcNo()
430  << " with neighbProcNo:" << pldui.neighbProcNo()
431  << exit(FatalError);
432  }
433  else if (agglom0 == myAgglom && agglom1 == myAgglom)
434  {
435  // Merged interface
436  if (debug)
437  {
438  Pout<< "merged interface: myProcNo:"
439  << pldui.myProcNo()
440  << " nbr:" << pldui.neighbProcNo()
441  << " size:" << ldui.faceCells().size()
442  << endl;
443  }
444 
445  label nbrProcMeshI = findIndex
446  (
447  procIDs,
448  pldui.neighbProcNo()
449  );
450 
451  if (procMeshI < nbrProcMeshI)
452  {
453  // I am 'master' since get lowest numbered cells
454  nCoupledFaces[procMeshI] += ldui.faceCells().size();
455  }
456 
458  mergedMap.find(procEdge);
459 
460  if (iter != mergedMap.end())
461  {
462  iter().append(labelPair(procMeshI, intI));
463  }
464  else
465  {
466  mergedMap.insert
467  (
468  procEdge,
469  labelPairList(1, labelPair(procMeshI, intI))
470  );
471  }
472  }
473  else
474  {
475  if (debug)
476  {
477  Pout<< "external interface: myProcNo:"
478  << pldui.myProcNo()
479  << " nbr:" << pldui.neighbProcNo()
480  << " size:" << ldui.faceCells().size()
481  << endl;
482  }
483 
485  unmergedMap.find(procEdge);
486 
487  if (iter != unmergedMap.end())
488  {
489  iter().append(labelPair(procMeshI, intI));
490  }
491  else
492  {
493  unmergedMap.insert
494  (
495  procEdge,
496  labelPairList(1, labelPair(procMeshI, intI))
497  );
498  }
499  }
500  }
501  else
502  {
503  // Still external (non proc) interface
504  FatalErrorIn("lduPrimitiveMesh::lduPrimitiveMesh(..)")
505  << "At mesh from processor " << procIDs[procMeshI]
506  << " have interface " << intI
507  << " of unhandled type " << interfaces[intI].type()
508  << exit(FatalError);
509 
510  nOtherInterfaces++;
511  }
512  }
513  }
514  }
515 
516 
517 
518  if (debug)
519  {
520  Pout<< "Remaining interfaces:" << endl;
521  forAllConstIter(EdgeMap<labelPairList>, unmergedMap, iter)
522  {
523  Pout<< " agglom procEdge:" << iter.key() << endl;
524  const labelPairList& elems = iter();
525  forAll(elems, i)
526  {
527  label procMeshI = elems[i][0];
528  label interfaceI = elems[i][1];
529  const lduInterfacePtrsList interfaces =
530  mesh(myMesh, otherMeshes, procMeshI).interfaces();
531 
532  const processorLduInterface& pldui =
533  refCast<const processorLduInterface>
534  (
535  interfaces[interfaceI]
536  );
537 
538  Pout<< " proc:" << procIDs[procMeshI]
539  << " interfaceI:" << interfaceI
540  << " between:" << pldui.myProcNo()
541  << " and:" << pldui.neighbProcNo()
542  << endl;
543  }
544  Pout<< endl;
545  }
546  }
547  if (debug)
548  {
549  Pout<< "Merged interfaces:" << endl;
550  forAllConstIter(EdgeMap<labelPairList>, mergedMap, iter)
551  {
552  Pout<< " agglom procEdge:" << iter.key() << endl;
553  const labelPairList& elems = iter();
554 
555  forAll(elems, i)
556  {
557  label procMeshI = elems[i][0];
558  label interfaceI = elems[i][1];
559  const lduInterfacePtrsList interfaces =
560  mesh(myMesh, otherMeshes, procMeshI).interfaces();
561  const processorLduInterface& pldui =
562  refCast<const processorLduInterface>
563  (
564  interfaces[interfaceI]
565  );
566 
567  Pout<< " proc:" << procIDs[procMeshI]
568  << " interfaceI:" << interfaceI
569  << " between:" << pldui.myProcNo()
570  << " and:" << pldui.neighbProcNo()
571  << endl;
572  }
573  Pout<< endl;
574  }
575  }
576 
577 
578  // Adapt faceOffsets for internal interfaces
579  faceOffsets.setSize(nMeshes+1);
580  faceOffsets[0] = 0;
581  faceMap.setSize(nMeshes);
582  for (label procMeshI = 0; procMeshI < nMeshes; procMeshI++)
583  {
584  const lduMesh& procMesh = mesh(myMesh, otherMeshes, procMeshI);
585  label nInternal = procMesh.lduAddr().lowerAddr().size();
586 
587  faceOffsets[procMeshI+1] =
588  faceOffsets[procMeshI]
589  + nInternal
590  + nCoupledFaces[procMeshI];
591 
592  labelList& map = faceMap[procMeshI];
593  map.setSize(nInternal);
594  forAll(map, i)
595  {
596  map[i] = faceOffsets[procMeshI] + i;
597  }
598  }
599 
600 
601  // Combine upper and lower
602  lowerAddr_.setSize(faceOffsets.last(), -1);
603  upperAddr_.setSize(lowerAddr_.size(), -1);
604 
605 
606  // Old internal faces and resolved coupled interfaces
607  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
608 
609  for (label procMeshI = 0; procMeshI < nMeshes; procMeshI++)
610  {
611  const lduMesh& procMesh = mesh(myMesh, otherMeshes, procMeshI);
612 
613  const labelUList& l = procMesh.lduAddr().lowerAddr();
614  const labelUList& u = procMesh.lduAddr().upperAddr();
615 
616  // Add internal faces
617  label allFaceI = faceOffsets[procMeshI];
618 
619  forAll(l, faceI)
620  {
621  lowerAddr_[allFaceI] = cellOffsets[procMeshI]+l[faceI];
622  upperAddr_[allFaceI] = cellOffsets[procMeshI]+u[faceI];
623  allFaceI++;
624  }
625 
626 
627  // Add merged interfaces
628  const lduInterfacePtrsList interfaces = procMesh.interfaces();
629 
630  forAll(interfaces, intI)
631  {
632  if (interfaces.set(intI))
633  {
634  const lduInterface& ldui = interfaces[intI];
635 
636  if (isA<processorLduInterface>(ldui))
637  {
638  const processorLduInterface& pldui =
639  refCast<const processorLduInterface>(ldui);
640 
641  // Look up corresponding interfaces
642  label myP = pldui.myProcNo();
643  label nbrP = pldui.neighbProcNo();
644  label nbrProcMeshI = findIndex(procIDs, nbrP);
645 
646  if (procMeshI < nbrProcMeshI)
647  {
648  // I am 'master' since my cell numbers will be lower
649  // since cells get added in procMeshI order.
650 
651  label agglom0 = procAgglomMap[myP];
652  label agglom1 = procAgglomMap[nbrP];
653 
655  mergedMap.find(edge(agglom0, agglom1));
656 
657  if (fnd != mergedMap.end())
658  {
659  const labelPairList& elems = fnd();
660 
661  // Find nbrP in elems
662  label nbrIntI = -1;
663  forAll(elems, i)
664  {
665  label procI = elems[i][0];
666  label interfaceI = elems[i][1];
667  const lduInterfacePtrsList interfaces =
668  mesh
669  (
670  myMesh,
671  otherMeshes,
672  procI
673  ).interfaces();
674  const processorLduInterface& pldui =
675  refCast<const processorLduInterface>
676  (
677  interfaces[interfaceI]
678  );
679 
680  if
681  (
682  elems[i][0] == nbrProcMeshI
683  && pldui.neighbProcNo() == procIDs[procMeshI]
684  )
685  {
686  nbrIntI = elems[i][1];
687  break;
688  }
689  }
690 
691 
692  if (nbrIntI == -1)
693  {
695  (
696  "lduPrimitiveMesh::lduPrimitiveMesh(..)"
697  ) << "elems:" << elems << abort(FatalError);
698  }
699 
700 
701  const lduInterfacePtrsList nbrInterfaces = mesh
702  (
703  myMesh,
704  otherMeshes,
705  nbrProcMeshI
706  ).interfaces();
707 
708 
709  const labelUList& faceCells =
710  ldui.faceCells();
711  const labelUList& nbrFaceCells =
712  nbrInterfaces[nbrIntI].faceCells();
713 
714  if (faceCells.size() != nbrFaceCells.size())
715  {
717  (
718  "lduPrimitiveMesh::lduPrimitiveMesh(..)"
719  ) << "faceCells:" << faceCells
720  << " nbrFaceCells:" << nbrFaceCells
721  << abort(FatalError);
722  }
723 
724 
725  labelList& bfMap =
726  boundaryFaceMap[procMeshI][intI];
727  labelList& nbrBfMap =
728  boundaryFaceMap[nbrProcMeshI][nbrIntI];
729 
730  bfMap.setSize(faceCells.size());
731  nbrBfMap.setSize(faceCells.size());
732 
733  forAll(faceCells, pfI)
734  {
735  lowerAddr_[allFaceI] =
736  cellOffsets[procMeshI]+faceCells[pfI];
737  bfMap[pfI] = allFaceI;
738  upperAddr_[allFaceI] =
739  cellOffsets[nbrProcMeshI]+nbrFaceCells[pfI];
740  nbrBfMap[pfI] = (-allFaceI-1);
741  allFaceI++;
742  }
743  }
744  }
745  }
746  }
747  }
748  }
749 
750 
751  // Sort upper-tri order
752  {
753  labelList oldToNew
754  (
755  upperTriOrder
756  (
757  cellOffsets.last(), //nCells
758  lowerAddr_,
759  upperAddr_
760  )
761  );
762 
763  forAll(faceMap, procMeshI)
764  {
765  labelList& map = faceMap[procMeshI];
766  forAll(map, i)
767  {
768  if (map[i] >= 0)
769  {
770  map[i] = oldToNew[map[i]];
771  }
772  else
773  {
774  label allFaceI = -map[i]-1;
775  map[i] = -oldToNew[allFaceI]-1;
776  }
777  }
778  }
779 
780 
781  inplaceReorder(oldToNew, lowerAddr_);
782  inplaceReorder(oldToNew, upperAddr_);
783 
784  forAll(boundaryFaceMap, procI)
785  {
786  const labelList& bMap = boundaryMap[procI];
787  forAll(bMap, intI)
788  {
789  if (bMap[intI] == -1)
790  {
791  // Merged interface
792  labelList& bfMap = boundaryFaceMap[procI][intI];
793 
794  forAll(bfMap, i)
795  {
796  if (bfMap[i] >= 0)
797  {
798  bfMap[i] = oldToNew[bfMap[i]];
799  }
800  else
801  {
802  label allFaceI = -bfMap[i]-1;
803  bfMap[i] = (-oldToNew[allFaceI]-1);
804  }
805  }
806  }
807  }
808  }
809  }
810 
811 
812  // Kept interfaces
813  // ~~~~~~~~~~~~~~~
814 
815  interfaces_.setSize(unmergedMap.size() + nOtherInterfaces);
816  primitiveInterfaces_.setSize(interfaces_.size());
817 
818  label allInterfaceI = 0;
819 
820  forAllConstIter(EdgeMap<labelPairList>, unmergedMap, iter)
821  {
822  const labelPairList& elems = iter();
823 
824  // Sort processors in increasing order so both sides walk through in
825  // same order.
826  labelPairList procPairs(elems.size());
827  forAll(elems, i)
828  {
829  const labelPair& elem = elems[i];
830  label procMeshI = elem[0];
831  label interfaceI = elem[1];
832  const lduInterfacePtrsList interfaces = mesh
833  (
834  myMesh,
835  otherMeshes,
836  procMeshI
837  ).interfaces();
838 
839  const processorLduInterface& pldui =
840  refCast<const processorLduInterface>
841  (
842  interfaces[interfaceI]
843  );
844  label myProcNo = pldui.myProcNo();
845  label nbrProcNo = pldui.neighbProcNo();
846  procPairs[i] = labelPair
847  (
848  min(myProcNo, nbrProcNo),
849  max(myProcNo, nbrProcNo)
850  );
851  }
852 
853  labelList order;
854  sortedOrder(procPairs, order);
855 
856  // Count
857  label n = 0;
858 
859  forAll(order, i)
860  {
861  const labelPair& elem = elems[order[i]];
862  label procMeshI = elem[0];
863  label interfaceI = elem[1];
864  const lduInterfacePtrsList interfaces = mesh
865  (
866  myMesh,
867  otherMeshes,
868  procMeshI
869  ).interfaces();
870 
871  n += interfaces[interfaceI].faceCells().size();
872  }
873 
874  // Size
875  labelField allFaceCells(n);
876  labelField allFaceRestrictAddressing(n);
877  n = 0;
878 
879  // Fill
880  forAll(order, i)
881  {
882  const labelPair& elem = elems[order[i]];
883  label procMeshI = elem[0];
884  label interfaceI = elem[1];
885  const lduInterfacePtrsList interfaces = mesh
886  (
887  myMesh,
888  otherMeshes,
889  procMeshI
890  ).interfaces();
891 
892  boundaryMap[procMeshI][interfaceI] = allInterfaceI;
893  labelList& bfMap = boundaryFaceMap[procMeshI][interfaceI];
894 
895  const labelUList& l = interfaces[interfaceI].faceCells();
896  bfMap.setSize(l.size());
897 
898  forAll(l, faceI)
899  {
900  allFaceCells[n] = cellOffsets[procMeshI]+l[faceI];
901  allFaceRestrictAddressing[n] = n;
902  bfMap[faceI] = n;
903  n++;
904  }
905  }
906 
907 
908  // Find out local and remote processor in new communicator
909 
910  label neighbProcNo = -1;
911 
912  // See what the two processors map onto
913 
914  if (iter.key()[0] == myAgglom)
915  {
916  if (iter.key()[1] == myAgglom)
917  {
919  (
920  "lduPrimitiveMesh::lduPrimitiveMesh(..)"
921  ) << "problem procEdge:" << iter.key()
922  << exit(FatalError);
923  }
924 
925  neighbProcNo = iter.key()[1];
926  }
927  else
928  {
929  if (iter.key()[1] != myAgglom)
930  {
932  (
933  "lduPrimitiveMesh::lduPrimitiveMesh(..)"
934  ) << "problem procEdge:" << iter.key()
935  << exit(FatalError);
936  }
937 
938  neighbProcNo = iter.key()[0];
939  }
940 
941  primitiveInterfaces_.set
942  (
943  allInterfaceI,
945  (
946  allInterfaceI,
947  interfaces_,
948  allFaceCells,
949  allFaceRestrictAddressing,
950  comm_,
951  myAgglom,
952  neighbProcNo,
953  tensorField(), // forwardT
954  Pstream::msgType() // tag
955  )
956  );
957  interfaces_.set(allInterfaceI, &primitiveInterfaces_[allInterfaceI]);
958 
959  if (debug)
960  {
961  Pout<< "Created " << interfaces_[allInterfaceI].type()
962  << " interface at " << allInterfaceI
963  << " comm:" << comm_
964  << " myProcNo:" << myAgglom
965  << " neighbProcNo:" << neighbProcNo
966  << " nFaces:" << allFaceCells.size()
967  << endl;
968  }
969 
970 
971  allInterfaceI++;
972  }
973 
974 
975  patchSchedule_ = nonBlockingSchedule<processorGAMGInterface>(interfaces_);
976 
977  if (debug)
978  {
979  checkUpperTriangular(cellOffsets.last(), lowerAddr_, upperAddr_);
980  }
981 }
982 
983 
984 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
985 
987 (
988  const lduMesh& myMesh,
989  const PtrList<lduPrimitiveMesh>& otherMeshes,
990  const label meshI
991 )
992 {
993  return (meshI == 0 ? myMesh : otherMeshes[meshI-1]);
994 }
995 
996 
998 (
999  const label comm,
1000  const lduMesh& mesh,
1001  const labelList& procIDs,
1002  PtrList<lduPrimitiveMesh>& otherMeshes
1003 )
1004 {
1005  // Force calculation of schedule (since does parallel comms)
1006  (void)mesh.lduAddr().patchSchedule();
1007 
1008  if (Pstream::myProcNo(comm) == procIDs[0])
1009  {
1010  otherMeshes.setSize(procIDs.size()-1);
1011 
1012  // Slave meshes
1013  for (label i = 1; i < procIDs.size(); i++)
1014  {
1015  //Pout<< "on master :"
1016  // << " receiving from slave " << procIDs[i] << endl;
1017 
1018  IPstream fromSlave
1019  (
1021  procIDs[i],
1022  0, // bufSize
1023  Pstream::msgType(),
1024  comm
1025  );
1026 
1027  label nCells = readLabel(fromSlave);
1028  labelList lowerAddr(fromSlave);
1029  labelList upperAddr(fromSlave);
1030  boolList validInterface(fromSlave);
1031 
1032 
1033  // Construct mesh without interfaces
1034  otherMeshes.set
1035  (
1036  i-1,
1037  new lduPrimitiveMesh
1038  (
1039  nCells,
1040  lowerAddr,
1041  upperAddr,
1042  comm,
1043  true // reuse
1044  )
1045  );
1046 
1047  // Construct GAMGInterfaces
1048  lduInterfacePtrsList newInterfaces(validInterface.size());
1049  forAll(validInterface, intI)
1050  {
1051  if (validInterface[intI])
1052  {
1053  word coupleType(fromSlave);
1054 
1055  newInterfaces.set
1056  (
1057  intI,
1059  (
1060  coupleType,
1061  intI,
1062  otherMeshes[i-1].rawInterfaces(),
1063  fromSlave
1064  ).ptr()
1065  );
1066  }
1067  }
1068 
1069  otherMeshes[i-1].addInterfaces
1070  (
1071  newInterfaces,
1072  nonBlockingSchedule<processorGAMGInterface>
1073  (
1074  newInterfaces
1075  )
1076  );
1077  }
1078  }
1079  else if (findIndex(procIDs, Pstream::myProcNo(comm)) != -1)
1080  {
1081  // Send to master
1082 
1083  const lduAddressing& addressing = mesh.lduAddr();
1084  lduInterfacePtrsList interfaces(mesh.interfaces());
1085  boolList validInterface(interfaces.size());
1086  forAll(interfaces, intI)
1087  {
1088  validInterface[intI] = interfaces.set(intI);
1089  }
1090 
1091  OPstream toMaster
1092  (
1094  procIDs[0],
1095  0,
1096  Pstream::msgType(),
1097  comm
1098  );
1099 
1100  toMaster
1101  << addressing.size()
1102  << addressing.lowerAddr()
1103  << addressing.upperAddr()
1104  << validInterface;
1105 
1106  forAll(interfaces, intI)
1107  {
1108  if (interfaces.set(intI))
1109  {
1110  const GAMGInterface& interface = refCast<const GAMGInterface>
1111  (
1112  interfaces[intI]
1113  );
1114 
1115  toMaster << interface.type();
1116  interface.write(toMaster);
1117  }
1118  }
1119  }
1120 }
1121 
1122 
1123 // ************************************************************************* //
label size() const
Return number of equations.
virtual const labelUList & upperAddr() const =0
Return upper addressing.
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Definition: PtrList.C:142
Map from edge (expressed as its endpoints) to value.
Definition: EdgeMap.H:47
label size() const
Return the number of elements in the PtrList.
Definition: PtrListI.H:32
bool set(const label) const
Is element set.
Definition: PtrListI.H:107
void addInterfaces(lduInterfacePtrsList &interfaces, const lduSchedule &ps)
Add interfaces to a mesh constructed without.
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurence of given element and return index,.
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: contiguous.H:49
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:28
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:61
virtual label comm() const =0
Return communicator used for parallel communication.
List< labelPair > labelPairList
List of labelPairs.
Definition: labelPair.H:49
virtual lduInterfacePtrsList interfaces() const =0
Return a list of pointers for each patch.
bool insert(const edge &, const T &newElmt)
Insert a new hashedEntry.
T & last()
Return the last element of the list.
Definition: UListI.H:131
void inplaceReorder(const labelUList &oldToNew, ListType &)
Inplace reorder the elements of a list.
void sortedOrder(const UList< T > &, labelList &order)
Generate the (stable) sort order for the list.
virtual const labelUList & faceCells() const =0
Return faceCell addressing.
A class for handling words, derived from string.
Definition: word.H:59
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 size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
Output inter-processor communications stream.
Definition: OPstream.H:50
Abstract base class for meshes which provide LDU addressing for the construction of lduMatrix and LDU...
Definition: lduMesh.H:51
dynamicFvMesh & mesh
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
T & first()
Return the first element of the list.
Definition: UListI.H:117
static const lduMesh & mesh(const lduMesh &mesh0, const PtrList< lduPrimitiveMesh > &otherMeshes, const label meshI)
Select either mesh0 (meshI is 0) or otherMeshes[meshI-1].
Namespace for OpenFOAM.
label readLabel(Istream &is)
Definition: label.H:64
Field< tensor > tensorField
Specialisation of Field<T> for tensor.
virtual const lduAddressing & lduAddr() const =0
Return ldu addressing.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
label n
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
static autoPtr< GAMGInterface > New(const label index, const lduInterfacePtrsList &coarseInterfaces, const lduInterface &fineInterface, const labelField &localRestrictAddressing, const labelField &neighbourRestrictAddressing, const label fineLevelIndex, const label coarseComm)
Return a pointer to a new interface created on freestore given.
void setSize(const label)
Reset size of List.
Definition: List.C:318
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
An abstract base class for implicitly-coupled interfaces e.g. processor and cyclic patches...
Definition: lduInterface.H:53
#define WarningIn(functionName)
Report a warning using Foam::Warning.
A List with indirect addressing.
Definition: fvMatrix.H:106
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual int myProcNo() const =0
Return processor number (rank in communicator)
#define forAll(list, i)
Definition: UList.H:421
GAMG agglomerated processor interface.
static int & msgType()
Message tag of standard messages.
Definition: UPstream.H:451
label size() const
Return the number of elements in the UList.
Definition: UListI.H:299
virtual const labelUList & lowerAddr() const =0
Return lower addressing.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:404
bool set(const label) const
Is element set.
Definition: UPtrListI.H:80
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
virtual int neighbProcNo() const =0
Return neigbour processor number (rank in communicator)
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:31
error FatalError
virtual const lduSchedule & patchSchedule() const =0
Input inter-processor communications stream.
Definition: IPstream.H:50
static void gather(const label comm, const lduMesh &mesh, const labelList &procIDs, PtrList< lduPrimitiveMesh > &otherMeshes)
Gather meshes from other processors onto procIDs[0].
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
Pair< label > labelPair
Label pair.
Definition: labelPair.H:48
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:106
Abstract base class for GAMG agglomerated interfaces.
Definition: GAMGInterface.H:51
The class contains the addressing required by the lduMatrix: upper, lower and losort.
To & refCast(From &r)
Reference type cast template function.
Definition: typeInfo.H:106
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
Definition: HashTable.C:139
bool operator()(const label a, const label b)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
procLess(const labelPairList &lst)
interfaceProperties interface(alpha1, U, mixture())
labelList findIndices(const ListType &, typename ListType::const_reference, const label start=0)
Find all occurences of given element. Linear search.
Less operator for pairs of <processor><index>
Simplest contrete lduMesh which stores the addressing needed by lduMatrix.
defineTypeNameAndDebug(combustionModel, 0)
An abstract base class for processor coupled interfaces.
prefixOSstream Pout(cout,"Pout")
Definition: IOstreams.H:53