1 /*=========================================================================
2 
3   Program:   Visualization Toolkit
4   Module:    vtkGenerateGlobalIds.cxx
5 
6   Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
7   All rights reserved.
8   See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
9 
10      This software is distributed WITHOUT ANY WARRANTY; without even
11      the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
12      PURPOSE.  See the above copyright notice for more information.
13 
14 =========================================================================*/
15 
16 // Hide VTK_DEPRECATED_IN_9_1_0() warning for this class
17 #define VTK_DEPRECATION_LEVEL 0
18 
19 #include "vtkGenerateGlobalIds.h"
20 
21 #include "vtkBoundingBox.h"
22 #include "vtkCellData.h"
23 #include "vtkCompositeDataIterator.h"
24 #include "vtkCompositeDataSet.h"
25 #include "vtkDIYExplicitAssigner.h"
26 #include "vtkDIYUtilities.h"
27 #include "vtkDataSet.h"
28 #include "vtkIdList.h"
29 #include "vtkIdTypeArray.h"
30 #include "vtkLogger.h"
31 #include "vtkMath.h"
32 #include "vtkMultiProcessController.h"
33 #include "vtkNew.h"
34 #include "vtkObjectFactory.h"
35 #include "vtkPointData.h"
36 #include "vtkPoints.h"
37 #include "vtkSMPThreadLocalObject.h"
38 #include "vtkSMPTools.h"
39 #include "vtkStaticPointLocator.h"
40 #include "vtkTuple.h"
41 #include "vtkUnsignedCharArray.h"
42 #include "vtkUnstructuredGrid.h"
43 
44 #include <algorithm>
45 #include <climits>
46 #include <tuple>
47 #include <type_traits>
48 
49 // clang-format off
50 #include "vtk_diy2.h"
51 // #define DIY_USE_SPDLOG
52 #include VTK_DIY2(diy/mpi.hpp)
53 #include VTK_DIY2(diy/master.hpp)
54 #include VTK_DIY2(diy/link.hpp)
55 #include VTK_DIY2(diy/reduce.hpp)
56 #include VTK_DIY2(diy/reduce-operations.hpp)
57 #include VTK_DIY2(diy/partners/swap.hpp)
58 #include VTK_DIY2(diy/assigner.hpp)
59 #include VTK_DIY2(diy/algorithms.hpp)
60 // clang-format on
61 
62 namespace impl
63 {
64 
AllReduceBounds(diy::mpi::communicator & comm,std::vector<vtkSmartPointer<vtkPoints>> points)65 static vtkBoundingBox AllReduceBounds(
66   diy::mpi::communicator& comm, std::vector<vtkSmartPointer<vtkPoints>> points)
67 {
68   vtkBoundingBox bbox;
69   for (auto& pts : points)
70   {
71     if (pts)
72     {
73       double bds[6];
74       pts->GetBounds(bds);
75       bbox.AddBounds(bds);
76     }
77   }
78   vtkDIYUtilities::AllReduce(comm, bbox);
79   return bbox;
80 }
81 
82 class ExplicitAssigner : public ::diy::StaticAssigner
83 {
84   std::vector<int> GIDs;
85 
86 public:
ExplicitAssigner(const std::vector<int> & counts)87   ExplicitAssigner(const std::vector<int>& counts)
88     : diy::StaticAssigner(
89         static_cast<int>(counts.size()), std::accumulate(counts.begin(), counts.end(), 0))
90     , GIDs(counts)
91   {
92     for (size_t cc = 1; cc < this->GIDs.size(); ++cc)
93     {
94       this->GIDs[cc] += this->GIDs[cc - 1];
95     }
96   }
97 
98   //! returns the process rank of the block with global id gid (need not be local)
rank(int gid) const99   int rank(int gid) const override
100   {
101     for (size_t cc = 0; cc < this->GIDs.size(); ++cc)
102     {
103       if (gid < this->GIDs[cc])
104       {
105         return static_cast<int>(cc);
106       }
107     }
108     abort();
109   }
110 
111   //! gets the local gids for a given process rank
local_gids(int rank,std::vector<int> & gids) const112   void local_gids(int rank, std::vector<int>& gids) const override
113   {
114     const auto min = rank == 0 ? 0 : this->GIDs[rank - 1];
115     const auto max = this->GIDs[rank];
116     gids.resize(max - min);
117     std::iota(gids.begin(), gids.end(), min);
118   }
119 };
120 
121 /**
122  * This is the main implementation of the global id generation algorithm.
123  * The code is similar for both point and cell ids generation except small
124  * differences that are implemented using `BlockTraits` and `BlockT`.
125  *
126  * The general algorithm can be described as:
127  * - sort points (or cells) globally so that all "coincident" points (or cells)
128  *   are within the same block;
129  * - merge coincident points (or cells) per block and assign unique ids for
130  *   unique points (or cells) -- note this is local to each block since we know
131  *   all coincident points are same block after earlier step;
132  * - uniquify the generated ids globally by exchanging information of local
133  *   unique id counts;
134  * - communicate back the assigned unique id to the source block where the point
135  *   (or cell) came from.
136  */
137 template <typename ElementBlockT>
GenerateIds(vtkDataObject * dobj,vtkGenerateGlobalIds * self,bool cell_centers)138 static bool GenerateIds(vtkDataObject* dobj, vtkGenerateGlobalIds* self, bool cell_centers)
139 {
140   const double tolerance = self->GetTolerance();
141 
142   self->UpdateProgress(0.0);
143   diy::mpi::communicator comm = vtkDIYUtilities::GetCommunicator(self->GetController());
144 
145   vtkLogStartScope(TRACE, "extract points");
146   auto datasets = vtkCompositeDataSet::GetDataSets(dobj);
147   datasets.erase(std::remove_if(datasets.begin(), datasets.end(),
148                    [&cell_centers](vtkDataSet* ds) {
149                      return ds == nullptr || ds->GetNumberOfPoints() == 0 ||
150                        (cell_centers && ds->GetNumberOfCells() == 0);
151                    }),
152     datasets.end());
153   const auto points = vtkDIYUtilities::ExtractPoints(datasets, cell_centers);
154   vtkLogEndScope("extract points");
155 
156   // get the bounds for the domain globally.
157   const diy::ContinuousBounds gdomain =
158     vtkDIYUtilities::Convert(impl::AllReduceBounds(comm, points));
159 
160   const int local_num_blocks = static_cast<int>(points.size());
161   vtkDIYExplicitAssigner assigner(comm, local_num_blocks, /*pow-of-2*/ true);
162 
163   diy::Master master(
164     comm, 1, -1, []() { return static_cast<void*>(new ElementBlockT); },
165     [](void* b) -> void { delete static_cast<ElementBlockT*>(b); });
166 
167   vtkLogStartScope(TRACE, "populate master");
168   std::vector<int> gids;
169   assigner.local_gids(comm.rank(), gids);
170   for (size_t lid = 0; lid < gids.size(); ++lid)
171   {
172     auto block = new ElementBlockT();
173     if (lid < points.size() && points[lid] != nullptr)
174     {
175       assert(datasets[lid] != nullptr);
176       block->Initialize(gids[lid], points[lid], datasets[lid]);
177     }
178 
179     auto link = new diy::RegularContinuousLink(3, gdomain, gdomain);
180     master.add(gids[lid], block, link);
181   }
182   vtkLogEndScope("populate master");
183   self->UpdateProgress(0.25);
184 
185   if (assigner.nblocks() > 1)
186   {
187     vtkLogStartScope(TRACE, "kdtree");
188     // use diy::kdtree to shuffle points around so that all spatially co-located
189     // points are within a block.
190     diy::kdtree(master, assigner, 3, gdomain, &ElementBlockT::Elements, /*hist_bins=*/512);
191     vtkLogEndScope("kdtree");
192   }
193   self->UpdateProgress(0.50);
194 
195   vtkLogStartScope(TRACE, "merge-points");
196   // iterate over all local blocks to give them unique ids.
197   master.foreach ([&tolerance](ElementBlockT* b,               // local block
198                     const diy::Master::ProxyWithLink&) -> void // communication proxy
199     { b->MergeElements(tolerance); });
200   vtkLogEndScope("merge-points");
201   self->UpdateProgress(0.75);
202 
203   // now communicate point ownership information and assign ids to locally owned
204   // points.
205   vtkLogStartScope(TRACE, "exchange-ownership-ids");
206   diy::all_to_all(master, assigner, [](ElementBlockT* b, const diy::ReduceProxy& rp) {
207     if (rp.round() == 0)
208     {
209       // now enqueue ownership information.
210       b->EnqueueOwnershipInformation(rp);
211     }
212     else
213     {
214       // now dequeue owership information and process locally to assign ids
215       // to locally owned points and flag ghost points.
216       b->DequeueOwnershipInformation(rp);
217     }
218   });
219   vtkLogEndScope("exchange-ownership-ids");
220 
221   // exchange unique ids count so that we can determine global id offsets
222   vtkLogStartScope(TRACE, "exchange-unique-ids");
223   diy::all_to_all(master, assigner, [](ElementBlockT* b, const diy::ReduceProxy& rp) {
224     if (rp.round() == 0)
225     {
226       for (int i = rp.gid() + 1; i < rp.nblocks(); ++i)
227       {
228         rp.enqueue(rp.out_link().target(i), b->UniqueElementsCount);
229       }
230     }
231     else
232     {
233       vtkIdType offset = 0;
234       for (int src_gid = 0; src_gid < rp.gid(); ++src_gid)
235       {
236         vtkIdType msg;
237         rp.dequeue(src_gid, msg);
238         offset += msg;
239       }
240       b->AddOffset(offset);
241     }
242   });
243   vtkLogEndScope("exchange-unique-ids");
244 
245   // exchange assigned ids.
246   vtkLogStartScope(TRACE, "exchange-assigned-ids");
247   diy::all_to_all(master, assigner, [](ElementBlockT* b, const diy::ReduceProxy& rp) {
248     if (rp.round() == 0)
249     {
250       b->EnqueueReplies(rp);
251     }
252     else
253     {
254       b->DequeueReplies(rp);
255     }
256   });
257   vtkLogEndScope("exchange-assigned-ids");
258 
259   // final back communication to assign ids to ghosted points.
260   vtkLogStartScope(TRACE, "exchange-ghosted-ids");
261   diy::all_to_all(master, assigner, [](ElementBlockT* b, const diy::ReduceProxy& rp) {
262     if (rp.round() == 0)
263     {
264       b->EnqueueGhostedIds(rp);
265     }
266     else
267     {
268       b->DequeueGhostedIds(rp);
269     }
270   });
271   vtkLogEndScope("exchange-ghosted-ids");
272   self->UpdateProgress(1.0);
273   return true;
274 }
275 }
276 
277 namespace
278 {
279 
280 /**
281  * This is the point type that keeps the coordinates for each point in the
282  * dataset as well as enough information to track where that point came from so
283  * that we can communicate back to the source once a unique global id has been
284  * assigned.
285  */
286 struct PointTT
287 {
288   static const int attr_type = vtkDataObject::POINT;
289 
290   vtkTuple<double, 3> coords;
291   int source_gid;
292   vtkIdType source_id;
293 
294   // note: there's loss of precision here, but that's okay. this is only used by
295   // DIY when building the kdtree
operator []__anon440fe3b50911::PointTT296   float operator[](unsigned int index) const { return static_cast<float>(this->coords[index]); }
297 
GetElements__anon440fe3b50911::PointTT298   static std::vector<PointTT> GetElements(int gid, vtkPoints* pts, vtkDataSet*)
299   {
300     std::vector<PointTT> elems;
301     elems.resize(pts->GetNumberOfPoints());
302     vtkSMPTools::For(
303       0, pts->GetNumberOfPoints(), [&elems, pts, gid](vtkIdType start, vtkIdType end) {
304         for (vtkIdType cc = start; cc < end; ++cc)
305         {
306           auto& pt = elems[cc];
307           pts->GetPoint(cc, pt.coords.GetData());
308           pt.source_gid = gid;
309           pt.source_id = cc;
310         }
311       });
312     return elems;
313   }
314 
Sort__anon440fe3b50911::PointTT315   static void Sort(std::vector<PointTT>& points)
316   {
317     // let's sort the points by source-id. This ensures that when a point is
318     // duplicated among multiple blocks, the block with lower block-id owns the
319     // point. Thus, keeping the numbering consistent.
320     std::sort(points.begin(), points.end(), [](const PointTT& a, const PointTT& b) {
321       return (a.source_gid == b.source_gid) ? (a.source_id < b.source_id)
322                                             : (a.source_gid < b.source_gid);
323     });
324   }
325 
GenerateMergeMap__anon440fe3b50911::PointTT326   static std::vector<vtkIdType> GenerateMergeMap(
327     const std::vector<PointTT>& points, double tolerance)
328   {
329     std::vector<vtkIdType> mergemap(points.size(), -1);
330     if (points.empty())
331     {
332       return mergemap;
333     }
334 
335     auto numPts = static_cast<vtkIdType>(points.size());
336 
337     vtkNew<vtkUnstructuredGrid> grid;
338     vtkNew<vtkPoints> pts;
339     pts->SetDataTypeToDouble();
340     pts->SetNumberOfPoints(numPts);
341     vtkSMPTools::For(0, numPts, [&](vtkIdType start, vtkIdType end) {
342       for (vtkIdType cc = start; cc < end; ++cc)
343       {
344         pts->SetPoint(cc, points[cc].coords.GetData());
345       }
346     });
347     grid->SetPoints(pts);
348 
349     vtkNew<vtkStaticPointLocator> locator;
350     locator->SetDataSet(grid);
351     locator->SetTolerance(tolerance);
352     locator->BuildLocator();
353     locator->MergePoints(tolerance, &mergemap[0]);
354     return mergemap;
355   }
356 };
357 
358 struct CellTT
359 {
360   static const int attr_type = vtkDataObject::CELL;
361   vtkTuple<double, 3> center;
362   int source_gid;
363   vtkIdType source_id;
364   std::vector<vtkIdType> point_ids;
365 
366   // note: there's loss of precision here, but that's okay. this is only used by
367   // DIY when building the kdtree
operator []__anon440fe3b50911::CellTT368   float operator[](unsigned int index) const { return static_cast<float>(this->center[index]); }
369 
GetElements__anon440fe3b50911::CellTT370   static std::vector<CellTT> GetElements(int gid, vtkPoints* centers, vtkDataSet* ds)
371   {
372     const vtkIdType ncells = ds->GetNumberOfCells();
373     assert(centers->GetNumberOfPoints() == ncells);
374 
375     std::vector<CellTT> elems(ncells);
376 
377     vtkSMPThreadLocalObject<vtkIdList> tlIdList;
378     if (ncells > 0)
379     {
380       // so that we can call GetCellPoints in vtkSMPTools::For
381       ds->GetCellPoints(0, tlIdList.Local());
382     }
383 
384     auto pt_gids = vtkIdTypeArray::SafeDownCast(ds->GetPointData()->GetGlobalIds());
385     assert(ncells == 0 || pt_gids != nullptr);
386     vtkSMPTools::For(0, ncells, [&](vtkIdType start, vtkIdType end) {
387       auto ids = tlIdList.Local();
388       for (vtkIdType cc = start; cc < end; ++cc)
389       {
390         auto& cell = elems[cc];
391         centers->GetPoint(cc, cell.center.GetData());
392         cell.source_gid = gid;
393         cell.source_id = cc;
394 
395         ds->GetCellPoints(cc, ids);
396         cell.point_ids.resize(ids->GetNumberOfIds());
397         for (vtkIdType kk = 0, max = ids->GetNumberOfIds(); kk < max; ++kk)
398         {
399           cell.point_ids[kk] = pt_gids->GetTypedComponent(ids->GetId(kk), 0);
400         }
401       }
402     });
403 
404     return elems;
405   }
406 
Sort__anon440fe3b50911::CellTT407   static void Sort(std::vector<CellTT>& cells)
408   {
409     // here, we are sorting such that for duplicated cells, we always order the
410     // cell on the lower block before the one on the higher block. This is
411     // essential to keep the cell numbering consistent.
412     std::sort(cells.begin(), cells.end(), [](const CellTT& lhs, const CellTT& rhs) {
413       return (lhs.point_ids == rhs.point_ids ? lhs.source_gid < rhs.source_gid
414                                              : lhs.point_ids < rhs.point_ids);
415     });
416   }
417 
GenerateMergeMap__anon440fe3b50911::CellTT418   static std::vector<vtkIdType> GenerateMergeMap(
419     const std::vector<CellTT>& cells, double vtkNotUsed(tolerance))
420   {
421     std::vector<vtkIdType> mergemap(cells.size(), -1);
422     if (cells.empty())
423     {
424       return mergemap;
425     }
426 
427     mergemap[0] = 0;
428     for (size_t cc = 1, max = cells.size(); cc < max; ++cc)
429     {
430       if (cells[cc - 1].point_ids == cells[cc].point_ids)
431       {
432         mergemap[cc] = mergemap[cc - 1];
433       }
434       else
435       {
436         mergemap[cc] = static_cast<vtkIdType>(cc);
437       }
438     }
439     return mergemap;
440   }
441 };
442 
443 struct MessageItemTT
444 {
445   vtkIdType elem_id;
446   vtkIdType index;
447 };
448 
449 template <typename ElementT>
450 struct BlockT
451 {
452 private:
Enqueue__anon440fe3b50911::BlockT453   void Enqueue(const diy::ReduceProxy& rp)
454   {
455     for (const auto& pair : this->OutMessage)
456     {
457       rp.enqueue(rp.out_link().target(pair.first), pair.second);
458     }
459     this->OutMessage.clear();
460   }
461 
462 public:
463   vtkDataSet* Dataset{ nullptr };
464   std::vector<ElementT> Elements;
465   std::vector<vtkIdType> MergeMap;
466   vtkIdType UniqueElementsCount{ 0 };
467   std::map<int, std::vector<MessageItemTT>> OutMessage;
468 
469   vtkSmartPointer<vtkIdTypeArray> GlobalIds;
470   vtkSmartPointer<vtkUnsignedCharArray> GhostArray;
471 
Initialize__anon440fe3b50911::BlockT472   void Initialize(int self_gid, vtkPoints* points, vtkDataSet* dataset)
473   {
474     this->Dataset = dataset;
475     this->Elements = ElementT::GetElements(self_gid, points, dataset);
476 
477     if (dataset)
478     {
479       this->GlobalIds = vtkSmartPointer<vtkIdTypeArray>::New();
480       this->GlobalIds->SetName(
481         ElementT::attr_type == vtkDataObject::POINT ? "GlobalPointIds" : "GlobalCellIds");
482       this->GlobalIds->SetNumberOfTuples(points->GetNumberOfPoints());
483       this->GlobalIds->FillValue(-1);
484       dataset->GetAttributes(ElementT::attr_type)->SetGlobalIds(this->GlobalIds);
485 
486       this->GhostArray = vtkSmartPointer<vtkUnsignedCharArray>::New();
487       this->GhostArray->SetName(vtkDataSetAttributes::GhostArrayName());
488       this->GhostArray->SetNumberOfTuples(points->GetNumberOfPoints());
489       this->GhostArray->FillValue(vtkDataSetAttributes::DUPLICATEPOINT);
490 
491       // we're only adding ghost points, not cells.
492       if (ElementT::attr_type == vtkDataObject::POINT)
493       {
494         dataset->GetAttributes(ElementT::attr_type)->AddArray(this->GhostArray);
495       }
496     }
497   }
498 
MergeElements__anon440fe3b50911::BlockT499   void MergeElements(double tolerance)
500   {
501     // sort to make elements on lower gid's the primary elements
502     ElementT::Sort(this->Elements);
503     this->MergeMap = ElementT::GenerateMergeMap(this->Elements, tolerance);
504 
505     std::vector<char> needs_replies(this->MergeMap.size());
506     for (size_t cc = 0, max = this->MergeMap.size(); cc < max; ++cc)
507     {
508       if (this->MergeMap[cc] != static_cast<vtkIdType>(cc))
509       {
510         needs_replies[this->MergeMap[cc]] = 1;
511       }
512     }
513 
514     // populate out-message.
515     for (size_t cc = 0, max = this->MergeMap.size(); cc < max; ++cc)
516     {
517       if (this->MergeMap[cc] == static_cast<vtkIdType>(cc))
518       {
519         const auto& elem = this->Elements[cc];
520 
521         MessageItemTT datum;
522         datum.elem_id = elem.source_id;
523         datum.index = needs_replies[cc] == 1 ? static_cast<vtkIdType>(cc) : vtkIdType(-1);
524         this->OutMessage[elem.source_gid].push_back(datum);
525       }
526     }
527   }
528 
EnqueueOwnershipInformation__anon440fe3b50911::BlockT529   void EnqueueOwnershipInformation(const diy::ReduceProxy& rp) { this->Enqueue(rp); }
DequeueOwnershipInformation__anon440fe3b50911::BlockT530   void DequeueOwnershipInformation(const diy::ReduceProxy& rp)
531   {
532     std::map<int, std::vector<MessageItemTT>> inmessage;
533     for (int i = 0; i < rp.in_link().size(); ++i)
534     {
535       const int in_gid = rp.in_link().target(i).gid;
536       while (rp.incoming(in_gid))
537       {
538         auto& ownerships = inmessage[in_gid];
539         rp.dequeue(in_gid, ownerships);
540       }
541     }
542 
543     // we should not have received any message if we don't have a dataset.
544     assert(this->Dataset != nullptr || inmessage.size() == 0);
545     if (!this->Dataset)
546     {
547       return;
548     }
549 
550     for (const auto& pair : inmessage)
551     {
552       for (const auto& data : pair.second)
553       {
554         this->GhostArray->SetTypedComponent(data.elem_id, 0, 0);
555       }
556     }
557 
558     // Assign global ids starting with 0 for locally owned elems.
559     this->UniqueElementsCount = 0;
560     for (vtkIdType cc = 0, max = this->GhostArray->GetNumberOfTuples(); cc < max; ++cc)
561     {
562       if (this->GhostArray->GetTypedComponent(cc, 0) == 0)
563       {
564         this->GlobalIds->SetTypedComponent(cc, 0, this->UniqueElementsCount);
565         this->UniqueElementsCount++;
566       }
567     }
568 
569     // Generate message send back assigned global ids to requesting blocks.
570     for (const auto& pair : inmessage)
571     {
572       for (const auto& data : pair.second)
573       {
574         if (data.index != -1)
575         {
576           MessageItemTT reply;
577           reply.index = data.index;
578           reply.elem_id = this->GlobalIds->GetTypedComponent(data.elem_id, 0);
579           this->OutMessage[pair.first].push_back(reply);
580         }
581       }
582     }
583   }
584 
AddOffset__anon440fe3b50911::BlockT585   void AddOffset(const vtkIdType offset)
586   {
587     if (this->GlobalIds == nullptr || offset == 0)
588     {
589       return;
590     }
591     vtkSMPTools::For(
592       0, this->GlobalIds->GetNumberOfTuples(), [&offset, this](vtkIdType start, vtkIdType end) {
593         for (vtkIdType cc = start; cc < end; ++cc)
594         {
595           const auto id = this->GlobalIds->GetTypedComponent(cc, 0);
596           if (id != -1)
597           {
598             this->GlobalIds->SetTypedComponent(cc, 0, id + offset);
599           }
600         }
601       });
602 
603     // offset replies too.
604     for (auto& pair : this->OutMessage)
605     {
606       for (auto& data : pair.second)
607       {
608         data.elem_id += offset;
609       }
610     }
611   }
612 
EnqueueReplies__anon440fe3b50911::BlockT613   void EnqueueReplies(const diy::ReduceProxy& rp) { this->Enqueue(rp); }
DequeueReplies__anon440fe3b50911::BlockT614   void DequeueReplies(const diy::ReduceProxy& rp)
615   {
616     for (int i = 0; i < rp.in_link().size(); ++i)
617     {
618       const int in_gid = rp.in_link().target(i).gid;
619       while (rp.incoming(in_gid))
620       {
621         std::vector<MessageItemTT> ownerships;
622         rp.dequeue(in_gid, ownerships);
623         for (const auto& data : ownerships)
624         {
625           // we're changing the id in our local storage to now be the global
626           // id.
627           this->Elements[data.index].source_id = data.elem_id;
628         }
629       }
630     }
631 
632     const auto& mergemap = this->MergeMap;
633     for (size_t cc = 0, max = mergemap.size(); cc < max; ++cc)
634     {
635       if (mergemap[cc] != static_cast<vtkIdType>(cc))
636       {
637         auto& original_elem = this->Elements[mergemap[cc]];
638         auto& duplicate_elem = this->Elements[cc];
639 
640         MessageItemTT data;
641         data.elem_id = original_elem.source_id;
642         data.index = duplicate_elem.source_id;
643         this->OutMessage[duplicate_elem.source_gid].push_back(data);
644       }
645     }
646   }
647 
EnqueueGhostedIds__anon440fe3b50911::BlockT648   void EnqueueGhostedIds(const diy::ReduceProxy& rp) { this->Enqueue(rp); }
DequeueGhostedIds__anon440fe3b50911::BlockT649   void DequeueGhostedIds(const diy::ReduceProxy& rp)
650   {
651     for (int i = 0; i < rp.in_link().size(); ++i)
652     {
653       const int in_gid = rp.in_link().target(i).gid;
654       while (rp.incoming(in_gid))
655       {
656         std::vector<MessageItemTT> ownerships;
657         rp.dequeue(in_gid, ownerships);
658         assert(this->Dataset != nullptr || ownerships.size() == 0);
659         for (const auto& data : ownerships)
660         {
661           this->GlobalIds->SetTypedComponent(data.index, 0, data.elem_id);
662         }
663       }
664     }
665   }
666 };
667 }
668 
669 namespace diy
670 {
671 template <>
672 struct Serialization<::CellTT>
673 {
savediy::Serialization674   static void save(BinaryBuffer& bb, const CellTT& c)
675   {
676     diy::save(bb, c.center);
677     diy::save(bb, c.source_gid);
678     diy::save(bb, c.source_id);
679     diy::save(bb, c.point_ids);
680   }
loaddiy::Serialization681   static void load(BinaryBuffer& bb, CellTT& c)
682   {
683     c.point_ids.clear();
684 
685     diy::load(bb, c.center);
686     diy::load(bb, c.source_gid);
687     diy::load(bb, c.source_id);
688     diy::load(bb, c.point_ids);
689   }
690 };
691 }
692 
693 vtkStandardNewMacro(vtkGenerateGlobalIds);
694 vtkCxxSetObjectMacro(vtkGenerateGlobalIds, Controller, vtkMultiProcessController);
695 //------------------------------------------------------------------------------
vtkGenerateGlobalIds()696 vtkGenerateGlobalIds::vtkGenerateGlobalIds()
697   : Controller(nullptr)
698   , Tolerance(0)
699 {
700   this->SetNumberOfInputPorts(1);
701   this->SetNumberOfOutputPorts(1);
702   this->SetController(vtkMultiProcessController::GetGlobalController());
703 }
704 
705 //------------------------------------------------------------------------------
~vtkGenerateGlobalIds()706 vtkGenerateGlobalIds::~vtkGenerateGlobalIds()
707 {
708   this->SetController(nullptr);
709 }
710 
711 //------------------------------------------------------------------------------
RequestData(vtkInformation *,vtkInformationVector ** inputVector,vtkInformationVector * outputVector)712 int vtkGenerateGlobalIds::RequestData(
713   vtkInformation*, vtkInformationVector** inputVector, vtkInformationVector* outputVector)
714 {
715   auto inputDO = vtkDataObject::GetData(inputVector[0], 0);
716   auto outputDO = vtkDataObject::GetData(outputVector, 0);
717   outputDO->ShallowCopy(inputDO);
718 
719   // generate point ids first.
720   {
721     this->SetProgressShiftScale(0, 0.5);
722     vtkLogScopeF(TRACE, "generate global point ids");
723     if (!impl::GenerateIds<BlockT<PointTT>>(outputDO, this, false))
724     {
725       this->SetProgressShiftScale(0, 1.0);
726       return 0;
727     }
728   }
729 
730   // generate cell ids next (this needs global point ids)
731   {
732     this->SetProgressShiftScale(0.5, 0.5);
733     vtkLogScopeF(TRACE, "generate global cell ids");
734     if (!impl::GenerateIds<BlockT<CellTT>>(outputDO, this, true))
735     {
736       this->SetProgressShiftScale(0, 1.0);
737       return 0;
738     }
739   }
740 
741   this->SetProgressShiftScale(0, 1.0);
742   return 1;
743 }
744 
745 //------------------------------------------------------------------------------
PrintSelf(ostream & os,vtkIndent indent)746 void vtkGenerateGlobalIds::PrintSelf(ostream& os, vtkIndent indent)
747 {
748   this->Superclass::PrintSelf(os, indent);
749   os << indent << "Controller: " << this->Controller << endl;
750   os << indent << "Tolerance: " << this->Tolerance << endl;
751 }
752