MueLu Version of the Day
Loading...
Searching...
No Matches
MueLu_Zoltan2GraphAdapter.hpp
Go to the documentation of this file.
1// @HEADER
2//
3// ***********************************************************************
4//
5// MueLu: A package for multigrid based preconditioning
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
39// Jonathan Hu (jhu@sandia.gov)
40// Andrey Prokopenko (aprokop@sandia.gov)
41// Ray Tuminaro (rstumin@sandia.gov)
42//
43// ***********************************************************************
44//
45// @HEADER
46
47#ifndef MUELU_ZOLTAN2GRAPHADAPTER_HPP_
48#define MUELU_ZOLTAN2GRAPHADAPTER_HPP_
49
50#include "MueLu_ConfigDefs.hpp"
51
52#if defined(HAVE_MUELU_ZOLTAN2)
53
54#include <Teuchos_RCP.hpp>
55#include <Teuchos_Comm.hpp>
56#include <Teuchos_ArrayView.hpp>
57#include <Xpetra_Map.hpp>
58#include <Zoltan2_InputTraits.hpp>
59#include <Zoltan2_GraphAdapter.hpp>
60#include <Zoltan2_StridedData.hpp>
61#include <Zoltan2_PartitioningSolution.hpp>
62#include "MueLu_GraphBase.hpp"
63
64
65
66// Zoltab2 InputTraits for MueLu Graph objects
67namespace Zoltan2 {
68
69template <typename LocalOrdinal,
70 typename GlobalOrdinal,
71 typename Node>
72struct InputTraits<MueLu::GraphBase<LocalOrdinal,GlobalOrdinal,Node> >
73{
74 typedef Zoltan2::default_scalar_t scalar_t;
77 typedef size_t offset_t;
78 typedef Zoltan2::default_part_t part_t;
79 typedef Node node_t;
80 static inline std::string name() {return "MueLu::Graph";}
81
82 Z2_STATIC_ASSERT_TYPES // validate the types
83};
84}//end namespace Zoltan2
85
86
87namespace MueLu {
88
89template <typename User, typename UserCoord=User>
90class MueLuGraphBaseAdapter : public Zoltan2::GraphAdapter<User,UserCoord> {
91public:
92
93#ifndef DOXYGEN_SHOULD_SKIP_THIS
94 typedef typename Zoltan2::InputTraits<User>::scalar_t scalar_t;
95 typedef typename Zoltan2::InputTraits<User>::offset_t offset_t;
96 typedef typename Zoltan2::InputTraits<User>::lno_t lno_t;
97 typedef typename Zoltan2::InputTraits<User>::gno_t gno_t;
98 typedef typename Zoltan2::InputTraits<User>::part_t part_t;
99 typedef typename Zoltan2::InputTraits<User>::node_t node_t;
100 typedef User xgraph_t;
101 typedef User user_t;
102 typedef UserCoord userCoord_t;
103#endif
104
106 const Teuchos::RCP< const Teuchos::Comm< int > > getComm() const { return graph_->GetComm();}
107 const Teuchos::RCP< const Xpetra::Map<lno_t, gno_t, node_t> > getRowMap() const { return graph_->GetDomainMap();}
108 const RCP< const Xpetra::Map<lno_t, gno_t, node_t> > getColMap() const {
109 // For some GraphBases' this is a ColMap, in others it is a seperate map that is
110 // only non-null in parallel.
111 Teuchos::RCP<const Xpetra::Map<lno_t,gno_t,node_t> > map = graph_->GetImportMap();
112 if(map.is_null()) map = graph_->GetDomainMap();
113 return map;
114 }
115 size_t getLocalNumEntries() const { return graph_->GetNodeNumEdges();}
116 size_t getLocalNumRows() const { return getRowMap()->getLocalNumElements();}
117 size_t getLocalNumCols() const { return getColMap()->getLocalNumElements();}
118
119 void getLocalRowView(lno_t LocalRow, Teuchos::ArrayView< const lno_t > &indices) const {
120 indices = graph_->getNeighborVertices(LocalRow);
121 }
122
123
124
128
138 MueLuGraphBaseAdapter(const RCP<const User> &ingraph,
139 int nVtxWeights=0, int nEdgeWeights=0);
140
153 void setWeights(const scalar_t *val, int stride, int idx);
154
170 void setVertexWeights(const scalar_t *val, int stride, int idx);
171
177 void setWeightIsDegree(int idx);
178
184 void setVertexWeightIsDegree(int idx);
185
208 void setEdgeWeights(const scalar_t *val, int stride, int idx);
209
212 RCP<const xgraph_t> getXpetraGraph() const { return graph_; }
213
216 RCP<const User> getUserGraph() const { return ingraph_; }
217
219 // The Adapter interface.
221
223 // The GraphAdapter interface.
225
226 // TODO: Assuming rows == objects;
227 // TODO: Need to add option for columns or nonzeros?
228 size_t getLocalNumVertices() const { return getLocalNumRows(); }
229
230 void getVertexIDsView(const gno_t *&ids) const
231 {
232 ids = NULL;
234 ids = getRowMap()->getLocalElementList().getRawPtr();
235 }
236
237 size_t getLocalNumEdges() const { return getLocalNumEntries(); }
238
239 void getEdgesView(const offset_t *&offsets, const gno_t *&adjIds) const
240 {
241 offsets = offs_.getRawPtr();
242 adjIds = (getLocalNumEdges() ? adjids_.getRawPtr() : NULL);
243 }
244
246
247 void getVertexWeightsView(const scalar_t *&weights, int &stride,
248 int idx) const
249 {
250 if(idx<0 || idx >= nWeightsPerVertex_)
251 {
252 std::ostringstream emsg;
253 emsg << __FILE__ << ":" << __LINE__
254 << " Invalid vertex weight index " << idx << std::endl;
255 throw std::runtime_error(emsg.str());
256 }
257
258
259 size_t length;
260 vertexWeights_[idx].getStridedList(length, weights, stride);
261 }
262
263 bool useDegreeAsVertexWeight(int idx) const {return vertexDegreeWeight_[idx];}
264
266
267 void getEdgeWeightsView(const scalar_t *&weights, int &stride, int idx) const
268 {
269 if(idx<0 || idx >= nWeightsPerEdge_)
270 {
271 std::ostringstream emsg;
272 emsg << __FILE__ << ":" << __LINE__
273 << " Invalid edge weight index " << idx << std::endl;
274 throw std::runtime_error(emsg.str());
275 }
276
277
278 size_t length;
279 edgeWeights_[idx].getStridedList(length, weights, stride);
280 }
281
282
283 template <typename Adapter>
284 void applyPartitioningSolution(const User &in, User *&out,
285 const Zoltan2::PartitioningSolution<Adapter> &solution) const {
286 TEUCHOS_TEST_FOR_EXCEPTION(1, std::invalid_argument,"applyPartitionlingSolution not implemeneted");
287}
288
289 template <typename Adapter>
290 void applyPartitioningSolution(const User &in, RCP<User> &out,
291 const Zoltan2::PartitioningSolution<Adapter> &solution) const {
292 TEUCHOS_TEST_FOR_EXCEPTION(1, std::invalid_argument,"applyPartitionlingSolution not implemeneted");
293 }
294
295
296private:
297
298 RCP<const User > ingraph_;
299 RCP<const xgraph_t > graph_;
300 RCP<const Teuchos::Comm<int> > comm_;
301
302 ArrayRCP<const offset_t> offs_;
303 ArrayRCP<const gno_t> adjids_;
304
306 ArrayRCP<Zoltan2::StridedData<lno_t, scalar_t> > vertexWeights_;
307 ArrayRCP<bool> vertexDegreeWeight_;
308
310 ArrayRCP<Zoltan2::StridedData<lno_t, scalar_t> > edgeWeights_;
311
313 ArrayRCP<Zoltan2::StridedData<lno_t, scalar_t> > coords_;
314
315};
316
317
319// Definitions
321
322template <typename User, typename UserCoord>
324 const RCP<const User> &ingraph, int nVtxWgts, int nEdgeWgts):
325 ingraph_(ingraph), graph_(), comm_() , offs_(), adjids_(),
326 nWeightsPerVertex_(nVtxWgts), vertexWeights_(), vertexDegreeWeight_(),
327 nWeightsPerEdge_(nEdgeWgts), edgeWeights_(),
328 coordinateDim_(0), coords_()
329{
330 typedef Zoltan2::StridedData<lno_t,scalar_t> input_t;
331 graph_ = ingraph;
332
333 comm_ = getRowMap()->getComm();
334 size_t nvtx = getLocalNumRows();
335 size_t nedges = getLocalNumEntries();
336
337 // Unfortunately we have to copy the offsets and edge Ids
338 // because edge Ids are not usually stored in vertex id order.
339 size_t n = nvtx + 1;
340 offs_.resize(n);
341 offset_t* offs = const_cast<offset_t*>(offs_.getRawPtr());
342 gno_t* adjids=0;
343 if(nedges > 0) {
344 adjids_.resize(nedges);
345 adjids = const_cast<gno_t*>(adjids_.getRawPtr());
346 }
347
348 offs[0] = 0;
349 for (size_t v=0; v < nvtx; v++){
350 ArrayView<const lno_t> nbors;
351 getLocalRowView(v, nbors);
352 offs[v+1] = offs[v] + nbors.size();
353 for (offset_t e=offs[v], i=0; e < offs[v+1]; e++) {
354 adjids[e] = getColMap()->getGlobalElement(nbors[i++]);
355 }
356 }
357
358 if (nWeightsPerVertex_ > 0) {
360 arcp(new input_t[nWeightsPerVertex_], 0, nWeightsPerVertex_, true);
362 arcp(new bool[nWeightsPerVertex_], 0, nWeightsPerVertex_, true);
363 for (int i=0; i < nWeightsPerVertex_; i++)
364 vertexDegreeWeight_[i] = false;
365 }
366
367
368}
369
371template <typename User, typename UserCoord>
373 const scalar_t *weightVal, int stride, int idx)
374{
375 if (this->getPrimaryEntityType() == Zoltan2::GRAPH_VERTEX)
376 setVertexWeights(weightVal, stride, idx);
377 else
378 setEdgeWeights(weightVal, stride, idx);
379}
380
382template <typename User, typename UserCoord>
384 const scalar_t *weightVal, int stride, int idx)
385{
386 typedef Zoltan2::StridedData<lno_t,scalar_t> input_t;
387
388 if(idx<0 || idx >= nWeightsPerVertex_)
389 {
390 std::ostringstream emsg;
391 emsg << __FILE__ << ":" << __LINE__
392 << " Invalid vertex weight index " << idx << std::endl;
393 throw std::runtime_error(emsg.str());
394 }
395
396 size_t nvtx = getLocalNumVertices();
397 ArrayRCP<const scalar_t> weightV(weightVal, 0, nvtx*stride, false);
398 vertexWeights_[idx] = input_t(weightV, stride);
399}
400
402template <typename User, typename UserCoord>
404 int idx)
405{
406 if (this->getPrimaryEntityType() == Zoltan2::GRAPH_VERTEX)
407 setVertexWeightIsDegree(idx);
408 else {
409 std::ostringstream emsg;
410 emsg << __FILE__ << "," << __LINE__
411 << " error: setWeightIsNumberOfNonZeros is supported only for"
412 << " vertices" << std::endl;
413 throw std::runtime_error(emsg.str());
414 }
415}
416
418template <typename User, typename UserCoord>
420 int idx)
421{
422 if(idx<0 || idx >= nWeightsPerVertex_)
423 {
424 std::ostringstream emsg;
425 emsg << __FILE__ << ":" << __LINE__
426 << " Invalid vertex weight index " << idx << std::endl;
427 throw std::runtime_error(emsg.str());
428 }
429
430 vertexDegreeWeight_[idx] = true;
431}
432
434template <typename User, typename UserCoord>
436 const scalar_t *weightVal, int stride, int idx)
437{
438 typedef Zoltan2::StridedData<lno_t,scalar_t> input_t;
439
440 if(idx<0 || idx >= nWeightsPerEdge_)
441 {
442 std::ostringstream emsg;
443 emsg << __FILE__ << ":" << __LINE__
444 << " Invalid edge weight index " << idx << std::endl;
445 throw std::runtime_error(emsg.str());
446 }
447
448 size_t nedges = getLocalNumEdges();
449 ArrayRCP<const scalar_t> weightV(weightVal, 0, nedges*stride, false);
450 edgeWeights_[idx] = input_t(weightV, stride);
451}
452
453
454} //namespace MueLu
455
456
457#endif// MUELU_HAVE_ZOLTAN2
458
459#endif
MueLu::DefaultLocalOrdinal LocalOrdinal
MueLu::DefaultGlobalOrdinal GlobalOrdinal
MueLu::DefaultNode Node
void setEdgeWeights(const scalar_t *val, int stride, int idx)
Provide a pointer to edge weights.
ArrayRCP< Zoltan2::StridedData< lno_t, scalar_t > > vertexWeights_
RCP< const Teuchos::Comm< int > > comm_
void applyPartitioningSolution(const User &in, RCP< User > &out, const Zoltan2::PartitioningSolution< Adapter > &solution) const
RCP< const User > getUserGraph() const
Access to user's graph.
ArrayRCP< Zoltan2::StridedData< lno_t, scalar_t > > coords_
void getEdgeWeightsView(const scalar_t *&weights, int &stride, int idx) const
MueLuGraphBaseAdapter(const RCP< const User > &ingraph, int nVtxWeights=0, int nEdgeWeights=0)
Constructor for graph with no weights or coordinates.
void getEdgesView(const offset_t *&offsets, const gno_t *&adjIds) const
const Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
MueLu::GraphBase Compatibility Layer.
ArrayRCP< Zoltan2::StridedData< lno_t, scalar_t > > edgeWeights_
void setVertexWeightIsDegree(int idx)
Specify an index for which the vertex weight should be the degree of the vertex.
void applyPartitioningSolution(const User &in, User *&out, const Zoltan2::PartitioningSolution< Adapter > &solution) const
void getLocalRowView(lno_t LocalRow, Teuchos::ArrayView< const lno_t > &indices) const
void setWeights(const scalar_t *val, int stride, int idx)
Provide a pointer to weights for the primary entity type.
void getVertexWeightsView(const scalar_t *&weights, int &stride, int idx) const
void getVertexIDsView(const gno_t *&ids) const
const Teuchos::RCP< const Xpetra::Map< lno_t, gno_t, node_t > > getRowMap() const
const RCP< const Xpetra::Map< lno_t, gno_t, node_t > > getColMap() const
void setVertexWeights(const scalar_t *val, int stride, int idx)
Provide a pointer to vertex weights.
RCP< const xgraph_t > getXpetraGraph() const
Access to Xpetra-wrapped user's graph.
void setWeightIsDegree(int idx)
Specify an index for which the weight should be the degree of the entity.
Namespace for MueLu classes and methods.