Zoltan2
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Zoltan2_GraphModel.hpp
Go to the documentation of this file.
1 // @HEADER
2 //
3 // ***********************************************************************
4 //
5 // Zoltan2: A package of combinatorial algorithms for scientific computing
6 // Copyright 2012 Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Karen Devine (kddevin@sandia.gov)
39 // Erik Boman (egboman@sandia.gov)
40 // Siva Rajamanickam (srajama@sandia.gov)
41 //
42 // ***********************************************************************
43 //
44 // @HEADER
45 
50 #ifndef _ZOLTAN2_GRAPHMODEL_HPP_
51 #define _ZOLTAN2_GRAPHMODEL_HPP_
52 
53 #include <Zoltan2_Model.hpp>
54 #include <Zoltan2_ModelHelpers.hpp>
55 #include <Zoltan2_InputTraits.hpp>
57 #include <Zoltan2_GraphAdapter.hpp>
60 #include <Zoltan2_MeshAdapter.hpp>
61 #include <Zoltan2_StridedData.hpp>
62 #include <unordered_map>
63 
64 namespace Zoltan2 {
65 
67 
79 template <typename Adapter>
80 class GraphModel : public Model<Adapter>
81 {
82 public:
83 
84 #ifndef DOXYGEN_SHOULD_SKIP_THIS
85  typedef typename Adapter::scalar_t scalar_t;
86  typedef typename Adapter::gno_t gno_t;
87  typedef typename Adapter::lno_t lno_t;
88  typedef typename Adapter::node_t node_t;
89  typedef typename Adapter::user_t user_t;
90  typedef typename Adapter::userCoord_t userCoord_t;
91  typedef StridedData<lno_t, scalar_t> input_t;
92  typedef typename Adapter::offset_t offset_t;
93 #endif
94 
97 
109  GraphModel(const RCP<const MatrixAdapter<user_t,userCoord_t> > &ia,
110  const RCP<const Environment> &env, const RCP<const Comm<int> > &comm,
111  modelFlag_t &modelFlags);
112 
113  GraphModel(const RCP<const GraphAdapter<user_t,userCoord_t> > &ia,
114  const RCP<const Environment> &env, const RCP<const Comm<int> > &comm,
115  modelFlag_t &modelFlags);
116 
117  GraphModel(const RCP<const MeshAdapter<user_t> > &ia,
118  const RCP<const Environment> &env, const RCP<const Comm<int> > &comm,
119  modelFlag_t &modelflags);
120 
121  GraphModel(const RCP<const VectorAdapter<userCoord_t> > &/* ia */,
122  const RCP<const Environment> &/* env */, const RCP<const Comm<int> > &/* comm */,
123  modelFlag_t &/* flags */)
124  {
125  throw std::runtime_error("cannot build GraphModel from VectorAdapter");
126  }
127 
128  GraphModel(const RCP<const IdentifierAdapter<user_t> > &ia,
129  const RCP<const Environment> &env, const RCP<const Comm<int> > &comm,
130  modelFlag_t &flags)
131  {
132  throw std::runtime_error("cannot build GraphModel from IdentifierAdapter");
133  }
134 
137  const RCP<const Comm<int> > getComm() { return comm_; }
138 
141  size_t getLocalNumVertices() const { return nLocalVertices_; }
142 
145  size_t getGlobalNumVertices() const { return nGlobalVertices_; }
146 
150  size_t getLocalNumEdges() const { return nLocalEdges_; }
151 
155  size_t getGlobalNumEdges() const { return nGlobalEdges_; }
156 
159  int getNumWeightsPerVertex() const { return nWeightsPerVertex_; }
160 
163  int getNumWeightsPerEdge() const { return nWeightsPerEdge_; }
164 
167  int getCoordinateDim() const { return vCoordDim_; }
168 
178  ArrayView<const gno_t> &Ids,
179  ArrayView<input_t> &wgts) const
180  {
181  Ids = vGids_.view(0, nLocalVertices_);
182  wgts = vWeights_.view(0, nWeightsPerVertex_);
183  return nLocalVertices_;
184  }
185 
192  size_t getVertexCoords(ArrayView<input_t> &xyz) const
193  {
194  xyz = vCoords_.view(0, vCoordDim_);
195  return nLocalVertices_;
196  }
197 
209  // Implied Vertex LNOs from getVertexList are used as indices to offsets
210  // array.
211  // Vertex GNOs are returned as neighbors in edgeIds.
212 
213  size_t getEdgeList( ArrayView<const gno_t> &edgeIds,
214  ArrayView<const offset_t> &offsets,
215  ArrayView<input_t> &wgts) const
216  {
217  edgeIds = eGids_.view(0, nLocalEdges_);
218  offsets = eOffsets_.view(0, nLocalVertices_+1);
219  wgts = eWeights_.view(0, nWeightsPerEdge_);
220  return nLocalEdges_;
221  }
222 
227  inline void getVertexDist(ArrayView<size_t> &vtxdist) const
228  {
229  vtxdist = vtxDist_();
230  if (vtxDist_.size() == 0) {
231  throw std::runtime_error("getVertexDist is available only "
232  "when consecutiveIdsRequired");
233  }
234  }
235 
237  // The Model interface.
239 
240  size_t getLocalNumObjects() const { return nLocalVertices_; }
241 
242  size_t getGlobalNumObjects() const { return nGlobalVertices_; }
243 
244 private:
245 
246  void shared_constructor(const RCP<const Adapter>&ia, modelFlag_t &modelFlags);
247 
248  template <typename AdapterWithCoords>
249  void shared_GetVertexCoords(const AdapterWithCoords *ia);
250 
251  void print(); // For debugging
252 
253  const RCP<const Environment > env_;
254  const RCP<const Comm<int> > comm_;
255 
256  bool localGraph_; // Flag indicating whether this graph is
257  // LOCAL with respect to the process;
258  // if !localGraph_, graph is GLOBAL with respect to
259  // the communicator.
260 
261 
262  size_t nLocalVertices_; // # local vertices in built graph
263  size_t nGlobalVertices_; // # global vertices in built graph
264  ArrayRCP<gno_t> vGids_; // vertices of graph built in model;
265  // may be same as adapter's input
266  // or may be renumbered 0 to (N-1).
267 
268  int nWeightsPerVertex_;
269  ArrayRCP<input_t> vWeights_;
270 
271  int vCoordDim_;
272  ArrayRCP<input_t> vCoords_;
273 
274  // Note: in some cases, size of these arrays
275  // may be larger than nLocalEdges_. So do not use .size().
276  // Use nLocalEdges_, nGlobalEdges_
277 
278  size_t nLocalEdges_; // # local edges in built graph
279  size_t nGlobalEdges_; // # global edges in built graph
280  ArrayRCP<gno_t> eGids_; // edges of graph built in model
281  ArrayRCP<offset_t> eOffsets_; // edge offsets build in model
282  // May be same as adapter's input
283  // or may differ
284  // due to renumbering, self-edge
285  // removal, or local graph.
286 
287  int nWeightsPerEdge_;
288  ArrayRCP<input_t> eWeights_; // edge weights in built graph
289  // May be same as adapter's input
290  // or may differ due to self-edge
291  // removal, or local graph.
292 
293  ArrayRCP<size_t> vtxDist_; // If consecutiveIdsRequired,
294  // vtxDist (as needed by ParMETIS
295  // and Scotch) is also created.
296  // Otherwise, it is Teuchos::null.
297 };
298 
299 
301 // GraphModel from MatrixAdapter
302 template <typename Adapter>
304  const RCP<const MatrixAdapter<user_t,userCoord_t> > &ia,
305  const RCP<const Environment> &env,
306  const RCP<const Comm<int> > &comm,
307  modelFlag_t &modelFlags):
308  env_(env),
309  comm_(comm),
310  localGraph_(false),
311  nLocalVertices_(0),
312  nGlobalVertices_(0),
313  vGids_(),
314  nWeightsPerVertex_(0),
315  vWeights_(),
316  vCoordDim_(0),
317  vCoords_(),
318  nLocalEdges_(0),
319  nGlobalEdges_(0),
320  eGids_(),
321  eOffsets_(),
322  nWeightsPerEdge_(0),
323  eWeights_(),
324  vtxDist_()
325 {
326  // Model creation flags
327  localGraph_ = modelFlags.test(BUILD_LOCAL_GRAPH);
328 
329  bool symTranspose = modelFlags.test(SYMMETRIZE_INPUT_TRANSPOSE);
330  bool symBipartite = modelFlags.test(SYMMETRIZE_INPUT_BIPARTITE);
331  bool vertexCols = modelFlags.test(VERTICES_ARE_MATRIX_COLUMNS);
332  bool vertexNz = modelFlags.test(VERTICES_ARE_MATRIX_NONZEROS);
333 
334  if (symTranspose || symBipartite || vertexCols || vertexNz){
335  throw std::runtime_error("graph build option not yet implemented");
336  }
337 
338  // Get the matrix from the input adapter
339  gno_t const *vtxIds=NULL, *nborIds=NULL;
340  offset_t const *offsets=NULL;
341  try{
342  nLocalVertices_ = ia->getLocalNumIDs();
343  ia->getIDsView(vtxIds);
344  }
346  try{
347  if (ia->CRSViewAvailable()) {
348  ia->getCRSView(offsets, nborIds);
349  }
350  else {
351  // TODO: Add support for CCS matrix layout
352  throw std::runtime_error("Only MatrixAdapter::getCRSView is supported "
353  "in graph model");
354  }
355  }
357 
358  // Save the pointers from the input adapter
359  nLocalEdges_ = offsets[nLocalVertices_];
360  vGids_ = arcp_const_cast<gno_t>(
361  arcp<const gno_t>(vtxIds, 0, nLocalVertices_, false));
362  eGids_ = arcp_const_cast<gno_t>(
363  arcp<const gno_t>(nborIds, 0, nLocalEdges_, false));
364  eOffsets_ = arcp_const_cast<offset_t>(
365  arcp<const offset_t>(offsets, 0, nLocalVertices_+1, false));
366 
367  // Edge weights
368  nWeightsPerEdge_ = 0; // no edge weights from a matrix yet.
369  // TODO: use matrix values as edge weights
370 
371  // Do constructor common to all adapters
372  shared_constructor(ia, modelFlags);
373 
374  // Get vertex coordinates, if available
375  if (ia->coordinatesAvailable()) {
376  typedef VectorAdapter<userCoord_t> adapterWithCoords_t;
377  shared_GetVertexCoords<adapterWithCoords_t>(ia->getCoordinateInput());
378  }
379  // print();
380 }
381 
382 
384 // GraphModel from GraphAdapter
385 template <typename Adapter>
387  const RCP<const GraphAdapter<user_t,userCoord_t> > &ia,
388  const RCP<const Environment> &env,
389  const RCP<const Comm<int> > &comm,
390  modelFlag_t &modelFlags):
391  env_(env),
392  comm_(comm),
393  localGraph_(false),
394  nLocalVertices_(0),
395  nGlobalVertices_(0),
396  vGids_(),
397  nWeightsPerVertex_(0),
398  vWeights_(),
399  vCoordDim_(0),
400  vCoords_(),
401  nLocalEdges_(0),
402  nGlobalEdges_(0),
403  eGids_(),
404  eOffsets_(),
405  nWeightsPerEdge_(0),
406  eWeights_(),
407  vtxDist_()
408 {
409  // Model creation flags
410  localGraph_ = modelFlags.test(BUILD_LOCAL_GRAPH);
411 
412  // This GraphModel is built with vertices == GRAPH_VERTEX from GraphAdapter.
413  // It is not ready to use vertices == GRAPH_EDGE from GraphAdapter.
414  env_->localInputAssertion(__FILE__, __LINE__,
415  "GraphModel from GraphAdapter is implemented only for "
416  "Graph Vertices as primary object, not for Graph Edges",
417  ia->getPrimaryEntityType() == Zoltan2::GRAPH_VERTEX, BASIC_ASSERTION);
418 
419  // Get the graph from the input adapter
420 
421  gno_t const *vtxIds=NULL, *nborIds=NULL;
422  offset_t const *offsets=NULL;
423  try{
424  nLocalVertices_ = ia->getLocalNumVertices();
425  ia->getVertexIDsView(vtxIds);
426  ia->getEdgesView(offsets, nborIds);
427  }
429 
430  // Save the pointers from the input adapter
431  nLocalEdges_ = offsets[nLocalVertices_];
432  vGids_ = arcp_const_cast<gno_t>(
433  arcp<const gno_t>(vtxIds, 0, nLocalVertices_, false));
434  eGids_ = arcp_const_cast<gno_t>(
435  arcp<const gno_t>(nborIds, 0, nLocalEdges_, false));
436  eOffsets_ = arcp_const_cast<offset_t>(
437  arcp<const offset_t>(offsets, 0, nLocalVertices_+1, false));
438 
439  // Edge weights
440  nWeightsPerEdge_ = ia->getNumWeightsPerEdge();
441  if (nWeightsPerEdge_ > 0){
442  input_t *wgts = new input_t [nWeightsPerEdge_];
443  eWeights_ = arcp(wgts, 0, nWeightsPerEdge_, true);
444  }
445 
446  for (int w=0; w < nWeightsPerEdge_; w++){
447  const scalar_t *ewgts=NULL;
448  int stride=0;
449 
450  ia->getEdgeWeightsView(ewgts, stride, w);
451 
452  ArrayRCP<const scalar_t> wgtArray(ewgts, 0, nLocalEdges_, false);
453  eWeights_[w] = input_t(wgtArray, stride);
454  }
455 
456  // Do constructor common to all adapters
457  shared_constructor(ia, modelFlags);
458 
459  // Get vertex coordinates, if available
460  if (ia->coordinatesAvailable()) {
461  typedef VectorAdapter<userCoord_t> adapterWithCoords_t;
462  shared_GetVertexCoords<adapterWithCoords_t>(ia->getCoordinateInput());
463  }
464  // print();
465 }
466 
468 // GraphModel from MeshAdapter
469 template <typename Adapter>
471  const RCP<const MeshAdapter<user_t> > &ia,
472  const RCP<const Environment> &env,
473  const RCP<const Comm<int> > &comm,
474  modelFlag_t &modelFlags):
475  env_(env),
476  comm_(comm),
477  localGraph_(false),
478  nLocalVertices_(0),
479  nGlobalVertices_(0),
480  vGids_(),
481  nWeightsPerVertex_(0),
482  vWeights_(),
483  vCoordDim_(0),
484  vCoords_(),
485  nLocalEdges_(0),
486  nGlobalEdges_(0),
487  eGids_(),
488  eOffsets_(),
489  nWeightsPerEdge_(0),
490  eWeights_(),
491  vtxDist_()
492 {
493  env_->timerStart(MACRO_TIMERS, "GraphModel constructed from MeshAdapter");
494 
495  // Model creation flags
496  localGraph_ = modelFlags.test(BUILD_LOCAL_GRAPH);
497 
498  // This GraphModel is built with vertices == ia->getPrimaryEntityType()
499  // from MeshAdapter.
500 
501  // Get the graph from the input adapter
502 
503  Zoltan2::MeshEntityType primaryEType = ia->getPrimaryEntityType();
504  Zoltan2::MeshEntityType secondAdjEType = ia->getSecondAdjacencyEntityType();
505 
506  // Get the IDs of the primary entity type; these are graph vertices
507 
508  gno_t const *vtxIds=NULL;
509  try {
510  nLocalVertices_ = ia->getLocalNumOf(primaryEType);
511  ia->getIDsViewOf(primaryEType, vtxIds);
512  }
514 
515  vGids_ = arcp_const_cast<gno_t>(
516  arcp<const gno_t>(vtxIds, 0, nLocalVertices_, false));
517 
518  // Get the second adjacencies to construct edges of the dual graph.
519 
520  if (!ia->avail2ndAdjs(primaryEType, secondAdjEType)) {
521  // KDDKDD TODO Want to do this differently for local and global graphs?
522  // KDDKDD TODO Currently getting global 2nd Adjs and filtering them for
523  // KDDKDD TODO local graphs. That approach is consistent with other
524  // KDDKDD TODO adapters, but is more expensive -- why build them just to
525  // KDDKDD TODO throw them away? Instead, perhaps should build
526  // KDDKDD TODO only local adjacencies.
527  // KDDKDD TODO Does it suffice to pass a serial comm for local graph?
528  try {
529  get2ndAdjsViewFromAdjs(ia, comm_, primaryEType, secondAdjEType, eOffsets_,
530  eGids_);
531  nLocalEdges_ = eOffsets_[nLocalVertices_];
532  }
534  }
535  else { // avail2ndAdjs
536  // Get the edges
537  try {
538  gno_t const *nborIds=NULL;
539  offset_t const *offsets=NULL;
540 
541  ia->get2ndAdjsView(primaryEType, secondAdjEType, offsets, nborIds);
542  // Save the pointers from the input adapter; we do not control the
543  // offsets and nborIds memory
544  nLocalEdges_ = offsets[nLocalVertices_];
545  eGids_ = arcp_const_cast<gno_t>(
546  arcp<const gno_t>(nborIds, 0, nLocalEdges_, false));
547  eOffsets_ = arcp_const_cast<offset_t>(
548  arcp<const offset_t>(offsets, 0, nLocalVertices_+1, false));
549  }
551  }
552 
553 
554  // Edge weights
555  // Cannot specify edge weights if Zoltan2 computes the second adjacencies;
556  // there's no way to know the correct order for the adjacencies and weights.
557  // InputAdapter must provide 2nd adjs in order for edge weights to be used.
558  if (ia->avail2ndAdjs(primaryEType, secondAdjEType)) {
559  nWeightsPerEdge_ = ia->getNumWeightsPer2ndAdj(primaryEType, secondAdjEType);
560  if (nWeightsPerEdge_ > 0){
561  input_t *wgts = new input_t [nWeightsPerEdge_];
562  eWeights_ = arcp(wgts, 0, nWeightsPerEdge_, true);
563  }
564 
565  for (int w=0; w < nWeightsPerEdge_; w++){
566  const scalar_t *ewgts=NULL;
567  int stride=0;
568 
569  ia->get2ndAdjWeightsView(primaryEType, secondAdjEType,
570  ewgts, stride, w);
571 
572  ArrayRCP<const scalar_t> wgtArray(ewgts, 0,
573  nLocalEdges_*stride, false);
574  eWeights_[w] = input_t(wgtArray, stride);
575  }
576  }
577 
578  // Do constructor common to all adapters
579  shared_constructor(ia, modelFlags);
580 
581  // Get vertex coordinates
582  typedef MeshAdapter<user_t> adapterWithCoords_t;
583  shared_GetVertexCoords<adapterWithCoords_t>(&(*ia));
584 
585  env_->timerStop(MACRO_TIMERS, "GraphModel constructed from MeshAdapter");
586  // print();
587 }
588 
589 
591 template <typename Adapter>
593  const RCP<const Adapter> &ia,
594  modelFlag_t &modelFlags)
595 {
596  // Model creation flags
597  bool consecutiveIdsRequired = modelFlags.test(GENERATE_CONSECUTIVE_IDS);
598  bool removeSelfEdges = modelFlags.test(REMOVE_SELF_EDGES);
599  bool subsetGraph = modelFlags.test(BUILD_SUBSET_GRAPH);
600 
601  // May modify the graph provided from input adapter; save pointers to
602  // the input adapter's data.
603  size_t adapterNLocalEdges = nLocalEdges_;
604  ArrayRCP<gno_t> adapterVGids = vGids_; // vertices of graph from adapter
605  ArrayRCP<gno_t> adapterEGids = eGids_; // edges of graph from adapter
606  ArrayRCP<offset_t> adapterEOffsets = eOffsets_; // edge offsets from adapter
607  ArrayRCP<input_t> adapterEWeights = eWeights_; // edge weights from adapter
608 
609  if (localGraph_) {
610  // Local graph is requested.
611  // Renumber vertices 0 to nLocalVertices_-1
612  // Filter out off-process edges
613  // Renumber edge neighbors to be in range [0,nLocalVertices_-1]
614 
615  // Allocate new space for local graph info
616  // Note that eGids_ and eWeights_[w] may be larger than needed;
617  // we would have to pre-count the number of local edges to avoid overalloc
618  vGids_ = arcp(new gno_t[nLocalVertices_],
619  0, nLocalVertices_, true);
620  eGids_ = arcp(new gno_t[adapterNLocalEdges],
621  0, adapterNLocalEdges, true);
622  eOffsets_ = arcp(new offset_t[nLocalVertices_+1],
623  0, nLocalVertices_+1, true);
624 
625  scalar_t **tmpEWeights = NULL;
626  if (nWeightsPerEdge_ > 0){
627  eWeights_ = arcp(new input_t[nWeightsPerEdge_], 0,
628  nWeightsPerEdge_, true);
629  // Need to use temporary array because StridedData has const data
630  // so we can't write to it.
631  tmpEWeights = new scalar_t*[nWeightsPerEdge_];
632  for (int w = 0; w < nWeightsPerEdge_; w++)
633  tmpEWeights[w] = new scalar_t[adapterNLocalEdges];
634  }
635 
636  // Build map between global and local vertex numbers
637  std::unordered_map<gno_t, lno_t> globalToLocal(nLocalVertices_);
638  for (size_t i = 0; i < nLocalVertices_; i++)
639  globalToLocal[adapterVGids[i]] = i;
640 
641  // Loop over edges; keep only those that are local (i.e., on-rank)
642  eOffsets_[0] = 0;
643  lno_t ecnt = 0;
644  for (size_t i = 0; i < nLocalVertices_; i++) {
645  vGids_[i] = gno_t(i);
646  for (offset_t j = adapterEOffsets[i]; j < adapterEOffsets[i+1]; j++) {
647 
648  if (removeSelfEdges && (adapterEGids[j] == adapterVGids[i]))
649  continue; // Skipping self edge
650 
651  // Determine whether neighbor vertex is local
652  typename std::unordered_map<gno_t, lno_t>::iterator localidx;
653  if ((localidx = globalToLocal.find(adapterEGids[j])) !=
654  globalToLocal.end()) {
655  // neighbor vertex is local
656  // Keep the edge and its weights
657  eGids_[ecnt] = localidx->second;
658  for (int w = 0; w < nWeightsPerEdge_; w++)
659  tmpEWeights[w][ecnt] = adapterEWeights[w][j];
660 
661  ecnt++;
662  }
663  }
664  eOffsets_[i+1] = ecnt;
665  }
666  nLocalEdges_ = eOffsets_[nLocalVertices_];
667  if (nWeightsPerEdge_) {
668  for (int w = 0; w < nWeightsPerEdge_; w++) {
669  ArrayRCP<const scalar_t> wgtArray(tmpEWeights[w],
670  0, adapterNLocalEdges, true);
671  eWeights_[w] = input_t(wgtArray, 0);
672  }
673  delete [] tmpEWeights;
674  }
675  } // localGraph_
676 
677  else if (consecutiveIdsRequired || removeSelfEdges || subsetGraph) {
678  // Build a Global graph
679  // If we are here, we expect SOMETHING in the graph to change from input:
680  // removing self edges, or converting to consecutive IDs, or subsetting
681  // the graph.
682 
683 
684  // Determine vertex GIDs for the global GraphModel
685  if (consecutiveIdsRequired) {
686  // Allocate new memory for vertices for consecutiveIds
687  vGids_ = arcp(new gno_t[nLocalVertices_], 0, nLocalVertices_, true);
688 
689  // Build vtxDist_ array with starting vGid on each rank
690  int np = comm_->getSize();
691  vtxDist_ = arcp(new size_t[np+1], 0, np+1, true);
692  vtxDist_[0] = 0;
693  Teuchos::gatherAll(*comm_, 1, &nLocalVertices_, np, &vtxDist_[1]);
694  for (int i = 0; i < np; i++)
695  vtxDist_[i+1] += vtxDist_[i];
696  }
697 
698  // Allocate new memory for edges and offsets, as needed
699  // Note that eGids_ may or may not be larger than needed;
700  // would have to pre-count number of edges kept otherwise
701  eGids_ = arcp(new gno_t[adapterNLocalEdges],
702  0, adapterNLocalEdges, true);
703 
704  scalar_t **tmpEWeights = NULL;
705  if (subsetGraph || removeSelfEdges) {
706  // May change number of edges and, thus, the offsets
707  eOffsets_ = arcp(new offset_t[nLocalVertices_+1],
708  0, nLocalVertices_+1, true);
709  eOffsets_[0] = 0;
710 
711  // Need to copy weights if remove edges
712  if (nWeightsPerEdge_ > 0){
713  eWeights_ = arcp(new input_t[nWeightsPerEdge_], 0,
714  nWeightsPerEdge_, true);
715  // Need to use temporary array because StridedData has const data
716  // so we can't write to it.
717  tmpEWeights = new scalar_t*[nWeightsPerEdge_];
718  for (int w = 0; w < nWeightsPerEdge_; w++)
719  tmpEWeights[w] = new scalar_t[adapterNLocalEdges];
720  }
721  }
722 
723  // If needed, determine the owning ranks and its local index off-proc
724  Teuchos::ArrayRCP<int> edgeRemoteRanks;
725  Teuchos::ArrayRCP<lno_t> edgeRemoteLids;
726  std::unordered_map<gno_t, size_t> edgeRemoteUniqueMap;
727 
728  if (subsetGraph || consecutiveIdsRequired) {
729 
730  // Find global minGID for map construction
731  gno_t myMinGID = std::numeric_limits<gno_t>::max();
732  size_t nVtx = adapterVGids.size();
733  for (size_t i = 0; i < nVtx; i++)
734  if (adapterVGids[i] < myMinGID) myMinGID = adapterVGids[i];
735 
736  gno_t minGID;
737  reduceAll<int, gno_t>(*comm_, Teuchos::REDUCE_MIN, 1,
738  &myMinGID, &minGID);
739 
740  gno_t dummy = Teuchos::OrdinalTraits<gno_t>::invalid();
741  Tpetra::Map<lno_t,gno_t> vtxMap(dummy, adapterVGids(), minGID, comm_);
742 
743  // Need to filter requested edges to make a unique list,
744  // as Tpetra::Map does not return correct info for duplicated entries
745  // (See bug 6412)
746  // The local filter may be more efficient anyway -- fewer communicated
747  // values in the Tpetra directory
748  Teuchos::ArrayRCP<gno_t> edgeRemoteUniqueGids =
749  arcp(new gno_t[adapterNLocalEdges], 0, adapterNLocalEdges, true);
750 
751  size_t nEdgeUnique = 0;
752  for (size_t i = 0; i < adapterNLocalEdges; i++) {
753  if (edgeRemoteUniqueMap.find(adapterEGids[i]) ==
754  edgeRemoteUniqueMap.end()) {
755  edgeRemoteUniqueGids[nEdgeUnique] = adapterEGids[i];
756  edgeRemoteUniqueMap[adapterEGids[i]] = nEdgeUnique;
757  nEdgeUnique++;
758  }
759  }
760 
761  edgeRemoteRanks = arcp(new int[nEdgeUnique], 0, nEdgeUnique, true);
762  edgeRemoteLids = arcp(new lno_t[nEdgeUnique], 0, nEdgeUnique, true);
763  vtxMap.getRemoteIndexList(edgeRemoteUniqueGids(0, nEdgeUnique),
764  edgeRemoteRanks(), edgeRemoteLids());
765  }
766 
767  // Renumber and/or filter the edges and vertices
768  lno_t ecnt = 0;
769  int me = comm_->getRank();
770  for (size_t i = 0; i < nLocalVertices_; i++) {
771 
772  if (consecutiveIdsRequired)
773  vGids_[i] = vtxDist_[me] + i;
774 
775  for (offset_t j = adapterEOffsets[i]; j < adapterEOffsets[i+1]; j++) {
776 
777  if (removeSelfEdges && (adapterVGids[i] == adapterEGids[j]))
778  continue; // Skipping self edge
779 
780  size_t remoteIdx = edgeRemoteUniqueMap[adapterEGids[j]];
781 
782  if (subsetGraph && (edgeRemoteRanks[remoteIdx] == -1))
783  continue; // Skipping edge with neighbor vertex that was not found
784  // in communicator
785 
786  if (consecutiveIdsRequired)
787  // Renumber edge using local number on remote rank plus first
788  // vtx number for remote rank
789  eGids_[ecnt] = edgeRemoteLids[remoteIdx]
790  + vtxDist_[edgeRemoteRanks[remoteIdx]];
791  else
792  eGids_[ecnt] = adapterEGids[j];
793 
794  if (subsetGraph || removeSelfEdges) {
795  // Need to copy weights only if number of edges might change
796  for (int w = 0; w < nWeightsPerEdge_; w++)
797  tmpEWeights[w][ecnt] = adapterEWeights[w][j];
798  }
799 
800  ecnt++;
801  }
802  if (subsetGraph || removeSelfEdges)
803  eOffsets_[i+1] = ecnt;
804  }
805  nLocalEdges_ = ecnt;
806  if (nWeightsPerEdge_ && (subsetGraph || removeSelfEdges)) {
807  for (int w = 0; w < nWeightsPerEdge_; w++) {
808  ArrayRCP<const scalar_t> wgtArray(tmpEWeights[w],
809  0, nLocalEdges_, true);
810  eWeights_[w] = input_t(wgtArray, 1);
811  }
812  delete [] tmpEWeights;
813  }
814  }
815 
816  // Vertex weights
817  nWeightsPerVertex_ = ia->getNumWeightsPerID();
818 
819  if (nWeightsPerVertex_ > 0){
820  input_t *weightInfo = new input_t [nWeightsPerVertex_];
821  env_->localMemoryAssertion(__FILE__, __LINE__, nWeightsPerVertex_,
822  weightInfo);
823 
824  for (int idx=0; idx < nWeightsPerVertex_; idx++){
825  bool useNumNZ = ia->useDegreeAsWeight(idx);
826  if (useNumNZ){
827  scalar_t *wgts = new scalar_t [nLocalVertices_];
828  env_->localMemoryAssertion(__FILE__, __LINE__, nLocalVertices_, wgts);
829  ArrayRCP<const scalar_t> wgtArray = arcp(wgts,
830  0, nLocalVertices_, true);
831  for (size_t i=0; i < nLocalVertices_; i++)
832  wgts[i] = eOffsets_[i+1] - eOffsets_[i];
833  weightInfo[idx] = input_t(wgtArray, 1);
834  }
835  else{
836  const scalar_t *weights=NULL;
837  int stride=0;
838  ia->getWeightsView(weights, stride, idx);
839  ArrayRCP<const scalar_t> wgtArray = arcp(weights, 0,
840  stride*nLocalVertices_,
841  false);
842  weightInfo[idx] = input_t(wgtArray, stride);
843  }
844  }
845 
846  vWeights_ = arcp<input_t>(weightInfo, 0, nWeightsPerVertex_, true);
847  }
848 
849  // Compute global values
850  if (localGraph_) {
851  nGlobalVertices_ = nLocalVertices_;
852  nGlobalEdges_ = nLocalEdges_;
853  }
854  else {
855  reduceAll<int, size_t>(*comm_, Teuchos::REDUCE_SUM, 1,
856  &nLocalVertices_, &nGlobalVertices_);
857  reduceAll<int, size_t>(*comm_, Teuchos::REDUCE_SUM, 1,
858  &nLocalEdges_, &nGlobalEdges_);
859  }
860 
861  env_->memory("After construction of graph model");
862 }
863 
865 
866 template <typename Adapter>
867 template <typename AdapterWithCoords>
868 void GraphModel<Adapter>::shared_GetVertexCoords(const AdapterWithCoords *ia)
869 {
870  // get pointers to vertex coordinates from input adapter
871 
872  vCoordDim_ = ia->getDimension();
873 
874  if (vCoordDim_ > 0){
875  input_t *coordInfo = new input_t [vCoordDim_];
876  env_->localMemoryAssertion(__FILE__, __LINE__, vCoordDim_, coordInfo);
877 
878  for (int dim=0; dim < vCoordDim_; dim++){
879  const scalar_t *coords=NULL;
880  int stride=0;
881  ia->getCoordinatesView(coords, stride, dim);
882  ArrayRCP<const scalar_t> coordArray = arcp(coords, 0,
883  stride*nLocalVertices_,
884  false);
885  coordInfo[dim] = input_t(coordArray, stride);
886  }
887 
888  vCoords_ = arcp<input_t>(coordInfo, 0, vCoordDim_, true);
889  }
890 }
891 
893  template <typename Adapter>
894 void GraphModel<Adapter>::print()
895 {
896 // if (env_->getDebugLevel() < VERBOSE_DETAILED_STATUS)
897 // return;
898 
899  std::ostream *os = env_->getDebugOStream();
900 
901  int me = comm_->getRank();
902 
903  *os << me
904  << " " << (localGraph_ ? "LOCAL GRAPH " : "GLOBAL GRAPH ")
905  << " Nvtx " << nLocalVertices_
906  << " Nedge " << nLocalEdges_
907  << " NVWgt " << nWeightsPerVertex_
908  << " NEWgt " << nWeightsPerEdge_
909  << " CDim " << vCoordDim_
910  << std::endl;
911 
912  for (size_t i = 0; i < nLocalVertices_; i++) {
913  *os << me << " " << i << " GID " << vGids_[i] << ": ";
914  for (offset_t j = eOffsets_[i]; j < eOffsets_[i+1]; j++)
915  *os << eGids_[j] << " " ;
916  *os << std::endl;
917  }
918 
919  if (nWeightsPerVertex_) {
920  for (size_t i = 0; i < nLocalVertices_; i++) {
921  *os << me << " " << i << " VWGTS " << vGids_[i] << ": ";
922  for (int j = 0; j < nWeightsPerVertex_; j++)
923  *os << vWeights_[j][i] << " ";
924  *os << std::endl;
925  }
926  }
927 
928  if (nWeightsPerEdge_) {
929  for (size_t i = 0; i < nLocalVertices_; i++) {
930  *os << me << " " << i << " EWGTS " << vGids_[i] << ": ";
931  for (offset_t j = eOffsets_[i]; j < eOffsets_[i+1]; j++) {
932  *os << eGids_[j] << " (";
933  for (int w = 0; w < nWeightsPerEdge_; w++)
934  *os << eWeights_[w][j] << " ";
935  *os << ") ";
936  }
937  *os << std::endl;
938  }
939  }
940 
941  if (vCoordDim_) {
942  for (size_t i = 0; i < nLocalVertices_; i++) {
943  *os << me << " " << i << " COORDS " << vGids_[i] << ": ";
944  for (int j = 0; j < vCoordDim_; j++)
945  *os << vCoords_[j][i] << " ";
946  *os << std::endl;
947  }
948  }
949  else
950  *os << me << " NO COORDINATES AVAIL " << std::endl;
951 }
952 
953 } // namespace Zoltan2
954 
955 
956 #endif
957 
use columns as graph vertices
size_t getLocalNumEdges() const
Returns the number of edges on this process. In global or subset graphs, includes off-process edges...
void get2ndAdjsViewFromAdjs(const Teuchos::RCP< const MeshAdapter< User > > &ia, const RCP< const Comm< int > > comm, Zoltan2::MeshEntityType sourcetarget, Zoltan2::MeshEntityType through, Teuchos::ArrayRCP< typename MeshAdapter< User >::offset_t > &offsets, Teuchos::ArrayRCP< typename MeshAdapter< User >::gno_t > &adjacencyIds)
size_t getVertexList(ArrayView< const gno_t > &Ids, ArrayView< input_t > &wgts) const
Sets pointers to this process&#39; vertex Ids and their weights.
Time an algorithm (or other entity) as a whole.
IdentifierAdapter defines the interface for identifiers.
ignore invalid neighbors
int getNumWeightsPerVertex() const
Returns the number (0 or greater) of weights per vertex.
size_t getGlobalNumObjects() const
Return the global number of objects.
fast typical checks for valid arguments
#define Z2_FORWARD_EXCEPTIONS
Forward an exception back through call stack.
MatrixAdapter defines the adapter interface for matrices.
Defines the Model interface.
GraphAdapter defines the interface for graph-based user data.
algorithm requires consecutive ids
Defines the MeshAdapter interface.
MeshAdapter defines the interface for mesh input.
std::bitset< NUM_MODEL_FLAGS > modelFlag_t
model must symmetrize input
Defines the IdentifierAdapter interface.
Defines the VectorAdapter interface.
GraphModel(const RCP< const MatrixAdapter< user_t, userCoord_t > > &ia, const RCP< const Environment > &env, const RCP< const Comm< int > > &comm, modelFlag_t &modelFlags)
Constructor.
model must symmetrize input
void getVertexDist(ArrayView< size_t > &vtxdist) const
Return the vtxDist array Array of size comm-&gt;getSize() + 1 Array[n+1] - Array[n] is number of vertice...
static ArrayRCP< ArrayRCP< zscalar_t > > weights
size_t getGlobalNumVertices() const
Returns the global number vertices.
Defines helper functions for use in the models.
algorithm requires no self edges
int getCoordinateDim() const
Returns the dimension (0 to 3) of vertex coordinates.
size_t getVertexCoords(ArrayView< input_t > &xyz) const
Sets pointers to this process&#39; vertex coordinates, if available. Order of coordinate info matches tha...
use nonzeros as graph vertices
size_t getEdgeList(ArrayView< const gno_t > &edgeIds, ArrayView< const offset_t > &offsets, ArrayView< input_t > &wgts) const
Sets pointers to this process&#39; edge (neighbor) global Ids, including off-process edges.
size_t getGlobalNumEdges() const
Returns the global number edges. For local graphs, the number of global edges is the number of local ...
VectorAdapter defines the interface for vector input.
GraphModel(const RCP< const VectorAdapter< userCoord_t > > &, const RCP< const Environment > &, const RCP< const Comm< int > > &, modelFlag_t &)
The StridedData class manages lists of weights or coordinates.
Traits for application input objects.
size_t getLocalNumVertices() const
Returns the number vertices on this process.
GraphModel defines the interface required for graph models.
MeshEntityType
Enumerate entity types for meshes: Regions, Faces, Edges, or Vertices.
size_t getLocalNumObjects() const
Return the local number of objects.
Defines the MatrixAdapter interface.
The base class for all model classes.
int getNumWeightsPerEdge() const
Returns the number (0 or greater) of weights per edge.
Defines the GraphAdapter interface.
GraphModel(const RCP< const IdentifierAdapter< user_t > > &ia, const RCP< const Environment > &env, const RCP< const Comm< int > > &comm, modelFlag_t &flags)
const RCP< const Comm< int > > getComm()
Return the communicator used by the model.
model represents graph within only one rank
This file defines the StridedData class.
Zoltan2::BasicUserTypes< zscalar_t, zlno_t, zgno_t > user_t
Definition: Metric.cpp:74