50 #ifndef _ZOLTAN2_TPETRAROWGRAPHADAPTER_HPP_
51 #define _ZOLTAN2_TPETRAROWGRAPHADAPTER_HPP_
53 #include "Kokkos_DualView.hpp"
54 #include "Kokkos_UnorderedMap.hpp"
55 #include <Tpetra_RowGraph.hpp>
84 template <
typename User,
typename UserCoord = User>
88 #ifndef DOXYGEN_SHOULD_SKIP_THIS
96 using userCoord_t = UserCoord;
111 int nEdgeWeights = 0);
135 void setWeightsDevice(
typename Base::ConstWeightsDeviceView1D val,
int idx);
145 void setWeightsHost(
typename Base::ConstWeightsHostView1D val,
int idx);
263 const gno_t *&adjIds)
const override;
267 typename Base::ConstIdsDeviceView &adjIds)
const override;
270 typename Base::ConstIdsHostView &adjIds)
const override;
275 int idx)
const override;
278 int idx = 0)
const override;
281 typename Base::WeightsDeviceView &weights)
const override;
284 int idx = 0)
const override;
287 typename Base::WeightsHostView &weights)
const override;
294 int idx)
const override;
297 int idx = 0)
const override;
300 typename Base::WeightsDeviceView &weights)
const override;
303 int idx = 0)
const override;
306 typename Base::WeightsHostView &weights)
const override;
308 template <
typename Adapter>
310 const User &in, User *&out,
313 template <
typename Adapter>
315 const User &in, RCP<User> &out,
321 const RCP<const User> &graph)
342 virtual RCP<User>
doMigration(
const User &from,
size_t numLocalRows,
343 const gno_t *myNewRows)
const;
350 template <
typename User,
typename UserCoord>
352 const RCP<const User> &ingraph,
int nVtxWgts,
int nEdgeWgts)
353 : graph_(ingraph), nWeightsPerVertex_(nVtxWgts),
354 nWeightsPerEdge_(nEdgeWgts), edgeWeights_() {
356 using localInds_t =
typename User::nonconst_local_inds_host_view_type;
358 const auto nvtx =
graph_->getLocalNumRows();
359 const auto nedges =
graph_->getLocalNumEntries();
361 const auto maxNumEntries =
graph_->getLocalMaxNumRowEntries();
366 adjIdsHost_ =
typename Base::ConstIdsHostView(
"adjIdsHost_", nedges);
367 offsHost_ =
typename Base::ConstOffsetsHostView(
"offsHost_", nvtx + 1);
369 localInds_t nbors(
"nbors", maxNumEntries);
371 for (
size_t v = 0; v < nvtx; v++) {
372 size_t numColInds = 0;
373 graph_->getLocalRowCopy(v, nbors, numColInds);
393 "vertexWeightsDevice_", nvtx, nWeightsPerVertex_);
396 "vertexDegreeWeightsHost_", nWeightsPerVertex_);
413 template <
typename User,
typename UserCoord>
415 const scalar_t *weightVal,
int stride,
int idx) {
417 setVertexWeights(weightVal, stride, idx);
419 setEdgeWeights(weightVal, stride, idx);
423 template <
typename User,
typename UserCoord>
425 typename Base::ConstWeightsDeviceView1D val,
int idx) {
427 setVertexWeightsDevice(val, idx);
429 setEdgeWeightsDevice(val, idx);
433 template <
typename User,
typename UserCoord>
435 typename Base::ConstWeightsHostView1D val,
int idx) {
437 setVertexWeightsHost(val, idx);
439 setEdgeWeightsHost(val, idx);
443 template <
typename User,
typename UserCoord>
445 const scalar_t *weightVal,
int stride,
int idx) {
447 "Invalid vertex weight index: " + std::to_string(idx));
449 size_t nvtx = getLocalNumVertices();
450 ArrayRCP<const scalar_t> weightV(weightVal, 0, nvtx * stride,
false);
451 vertexWeights_[idx] = input_t(weightV, stride);
455 template <
typename User,
typename UserCoord>
457 typename Base::ConstWeightsDeviceView1D
weights,
int idx) {
460 "Invalid vertex weight index: " + std::to_string(idx));
465 Kokkos::parallel_for(
466 vertexWeightsDevice_.extent(0), KOKKOS_CLASS_LAMBDA(
const int vertexID) {
467 vertexWeightsDevice_(vertexID, idx) =
weights(vertexID);
474 template <
typename User,
typename UserCoord>
476 typename Base::ConstWeightsHostView1D weightsHost,
int idx) {
478 "Invalid vertex weight index: " + std::to_string(idx));
480 auto weightsDevice = Kokkos::create_mirror_view_and_copy(
483 setVertexWeightsDevice(weightsDevice, idx);
487 template <
typename User,
typename UserCoord>
490 "setWeightIsNumberOfNonZeros is supported only for vertices");
492 setVertexWeightIsDegree(idx);
496 template <
typename User,
typename UserCoord>
499 "Invalid vertex weight index.");
501 vertexDegreeWeightsHost_(idx) =
true;
505 template <
typename User,
typename UserCoord>
507 const scalar_t *weightVal,
int stride,
int idx) {
511 "Invalid edge weight index" + std::to_string(idx));
513 size_t nedges = getLocalNumEdges();
514 ArrayRCP<const scalar_t> weightV(weightVal, 0, nedges * stride,
false);
515 edgeWeights_[idx] = input_t(weightV, stride);
519 template <
typename User,
typename UserCoord>
521 typename Base::ConstWeightsDeviceView1D
weights,
int idx) {
523 "Invalid edge weight index.");
528 Kokkos::parallel_for(
529 edgeWeightsDevice_.extent(0), KOKKOS_CLASS_LAMBDA(
const int vertexID) {
530 edgeWeightsDevice_(vertexID, idx) =
weights(vertexID);
537 template <
typename User,
typename UserCoord>
539 typename Base::ConstWeightsHostView1D weightsHost,
int idx) {
541 "Invalid edge weight index.");
543 auto weightsDevice = Kokkos::create_mirror_view_and_copy(
546 setEdgeWeightsDevice(weightsDevice);
550 template <
typename User,
typename UserCoord>
552 return graph_->getLocalNumRows();
556 template <
typename User,
typename UserCoord>
558 const gno_t *&ids)
const {
560 if (getLocalNumVertices())
561 ids = graph_->getRowMap()->getLocalElementList().getRawPtr();
565 template <
typename User,
typename UserCoord>
567 typename Base::ConstIdsDeviceView &ids)
const {
571 auto idsDevice = graph_->getRowMap()->getMyGlobalIndices();
572 auto tmpIds =
typename Base::IdsDeviceView(
"", idsDevice.extent(0));
574 Kokkos::deep_copy(tmpIds, idsDevice);
580 template <
typename User,
typename UserCoord>
582 typename Base::ConstIdsHostView &ids)
const {
585 auto idsDevice = graph_->getRowMap()->getMyGlobalIndices();
586 auto tmpIds =
typename Base::IdsHostView(
"", idsDevice.extent(0));
588 Kokkos::deep_copy(tmpIds, idsDevice);
594 template <
typename User,
typename UserCoord>
596 return graph_->getLocalNumEntries();
600 template <
typename User,
typename UserCoord>
603 offsets = offsHost_.data();
604 adjIds = (getLocalNumEdges() ? adjIdsHost_.data() : NULL);
608 template <
typename User,
typename UserCoord>
610 typename Base::ConstOffsetsDeviceView &offsets,
611 typename Base::ConstIdsDeviceView &adjIds)
const {
613 offsets = offsDevice_;
614 adjIds = adjIdsDevice_;
618 template <
typename User,
typename UserCoord>
620 typename Base::ConstOffsetsHostView &offsets,
621 typename Base::ConstIdsHostView &adjIds)
const {
623 auto hostIDs = Kokkos::create_mirror_view(adjIdsDevice_);
624 Kokkos::deep_copy(hostIDs, adjIdsDevice_);
627 auto hostOffsets = Kokkos::create_mirror_view(offsDevice_);
628 Kokkos::deep_copy(hostOffsets, offsDevice_);
629 offsets = hostOffsets;
633 template <
typename User,
typename UserCoord>
635 return nWeightsPerVertex_;
639 template <
typename User,
typename UserCoord>
641 const scalar_t *&
weights,
int &stride,
int idx)
const {
644 "Invalid vertex weight index.");
647 vertexWeights_[idx].getStridedList(length, weights, stride);
651 template <
typename User,
typename UserCoord>
653 typename Base::WeightsDeviceView1D &
weights,
int idx)
const {
655 "Invalid vertex weight index.");
657 const auto size = vertexWeightsDevice_.extent(0);
658 weights =
typename Base::WeightsDeviceView1D(
"weights", size);
660 Kokkos::parallel_for(
661 size, KOKKOS_CLASS_LAMBDA(
const int id) {
662 weights(
id) = vertexWeightsDevice_(
id, idx);
669 template <
typename User,
typename UserCoord>
671 typename Base::WeightsDeviceView &
weights)
const {
673 weights = vertexWeightsDevice_;
677 template <
typename User,
typename UserCoord>
679 typename Base::WeightsHostView1D &
weights,
int idx)
const {
681 "Invalid vertex weight index.");
683 auto weightsDevice =
typename Base::WeightsDeviceView1D(
684 "weights", vertexWeightsDevice_.extent(0));
685 getVertexWeightsDeviceView(weightsDevice, idx);
687 weights = Kokkos::create_mirror_view(weightsDevice);
688 Kokkos::deep_copy(weights, weightsDevice);
692 template <
typename User,
typename UserCoord>
694 typename Base::WeightsHostView &
weights)
const {
696 weights = Kokkos::create_mirror_view(vertexWeightsDevice_);
697 Kokkos::deep_copy(weights, vertexWeightsDevice_);
701 template <
typename User,
typename UserCoord>
704 return vertexDegreeWeightsHost_(idx);
708 template <
typename User,
typename UserCoord>
710 return nWeightsPerEdge_;
714 template <
typename User,
typename UserCoord>
716 const scalar_t *&
weights,
int &stride,
int idx)
const {
718 "Invalid edge weight index.");
721 edgeWeights_[idx].getStridedList(length, weights, stride);
725 template <
typename User,
typename UserCoord>
727 typename Base::WeightsDeviceView1D &
weights,
int idx)
const {
729 weights = Kokkos::subview(edgeWeightsDevice_, Kokkos::ALL, idx);
733 template <
typename User,
typename UserCoord>
735 typename Base::WeightsDeviceView &
weights)
const {
737 weights = edgeWeightsDevice_;
741 template <
typename User,
typename UserCoord>
743 typename Base::WeightsHostView1D &
weights,
int idx)
const {
745 auto weightsDevice = Kokkos::subview(edgeWeightsDevice_, Kokkos::ALL, idx);
746 weights = Kokkos::create_mirror_view(weightsDevice);
747 Kokkos::deep_copy(weights, weightsDevice);
751 template <
typename User,
typename UserCoord>
753 typename Base::WeightsHostView &
weights)
const {
755 weights = Kokkos::create_mirror_view(edgeWeightsDevice_);
756 Kokkos::deep_copy(weights, edgeWeightsDevice_);
760 template <
typename User,
typename UserCoord>
761 template <
typename Adapter>
763 const User &in, User *&out,
767 ArrayRCP<gno_t> importList;
770 Zoltan2::getImportList<Adapter, TpetraRowGraphAdapter<User, UserCoord>>(
771 solution,
this, importList);
776 RCP<User> outPtr = doMigration(in, numNewVtx, importList.getRawPtr());
782 template <
typename User,
typename UserCoord>
783 template <
typename Adapter>
785 const User &in, RCP<User> &out,
789 ArrayRCP<gno_t> importList;
792 Zoltan2::getImportList<Adapter, TpetraRowGraphAdapter<User, UserCoord>>(
793 solution,
this, importList);
798 out = doMigration(in, numNewVtx, importList.getRawPtr());
802 template <
typename User,
typename UserCoord>
804 const User &from,
size_t numLocalRows,
const gno_t *myNewRows)
const {
805 typedef Tpetra::Map<lno_t, gno_t, node_t>
map_t;
806 typedef Tpetra::CrsGraph<lno_t, gno_t, node_t> tcrsgraph_t;
817 const tcrsgraph_t *pCrsGraphSrc =
dynamic_cast<const tcrsgraph_t *
>(&from);
820 throw std::logic_error(
"TpetraRowGraphAdapter cannot migrate data for "
821 "your RowGraph; it can migrate data only for "
823 "You can inherit from TpetraRowGraphAdapter and "
824 "implement migration for your RowGraph.");
828 const RCP<const map_t> &smap = from.getRowMap();
829 int oldNumElts = smap->getLocalNumElements();
830 gno_t numGlobalRows = smap->getGlobalNumElements();
831 gno_t base = smap->getMinAllGlobalIndex();
834 ArrayView<const gno_t> rowList(myNewRows, numLocalRows);
835 const RCP<const Teuchos::Comm<int>> &comm = from.getComm();
836 RCP<const map_t> tmap = rcp(
new map_t(numGlobalRows, rowList, base, comm));
839 Tpetra::Import<lno_t, gno_t, node_t> importer(smap, tmap);
842 typedef Tpetra::Vector<gno_t, lno_t, gno_t, node_t> vector_t;
843 vector_t numOld(smap);
844 vector_t numNew(tmap);
845 for (
int lid = 0; lid < oldNumElts; lid++) {
846 numOld.replaceGlobalValue(smap->getGlobalElement(lid),
847 from.getNumEntriesInLocalRow(lid));
849 numNew.doImport(numOld, importer, Tpetra::INSERT);
851 size_t numElts = tmap->getLocalNumElements();
852 ArrayRCP<const gno_t> nnz;
854 nnz = numNew.getData(0);
856 ArrayRCP<const size_t> nnz_size_t;
858 if (numElts &&
sizeof(
gno_t) !=
sizeof(
size_t)) {
859 size_t *vals =
new size_t[numElts];
860 nnz_size_t = arcp(vals, 0, numElts,
true);
861 for (
size_t i = 0; i < numElts; i++) {
862 vals[i] =
static_cast<size_t>(nnz[i]);
865 nnz_size_t = arcp_reinterpret_cast<
const size_t>(nnz);
869 RCP<tcrsgraph_t> G = rcp(
new tcrsgraph_t(tmap, nnz_size_t()));
871 G->doImport(*pCrsGraphSrc, importer, Tpetra::INSERT);
873 return Teuchos::rcp_dynamic_cast<User>(G);
TpetraRowGraphAdapter(int nVtxWgts, int nEdgeWgts, const RCP< const User > &graph)
ArrayRCP< StridedData< lno_t, scalar_t > > edgeWeights_
Base::ConstIdsHostView adjIdsHost_
void getVertexIDsHostView(typename Base::ConstIdsHostView &ids) const override
Sets pointers to this process' graph entries.
void setVertexWeightsHost(typename Base::ConstWeightsHostView1D val, int idx)
Provide a host view to vertex weights.
void getEdgeWeightsView(const scalar_t *&weights, int &stride, int idx) const override
Provide a pointer to the edge weights, if any.
Helper functions for Partitioning Problems.
void setWeightIsDegree(int idx)
Specify an index for which the weight should be the degree of the entity.
void getVertexIDsView(const gno_t *&ids) const override
Sets pointers to this process' graph entries.
#define Z2_FORWARD_EXCEPTIONS
Forward an exception back through call stack.
typename InputTraits< User >::scalar_t scalar_t
Base::WeightsDeviceView edgeWeightsDevice_
static void AssertCondition(bool condition, const std::string &message, const char *file=__FILE__, int line=__LINE__)
void getEdgesView(const offset_t *&offsets, const gno_t *&adjIds) const override
Gets adjacency lists for all vertices in a compressed sparse row (CSR) format.
void getVertexWeightsDeviceView(typename Base::WeightsDeviceView1D &weights, int idx=0) const override
Provide a device view of the vertex weights, if any.
Base::ConstOffsetsDeviceView offsDevice_
GraphAdapter defines the interface for graph-based user data.
int getNumWeightsPerVertex() const override
Returns the number (0 or greater) of weights per vertex.
void setVertexWeights(const scalar_t *val, int stride, int idx)
Provide a pointer to vertex weights.
size_t getLocalNumVertices() const override
Returns the number of vertices on this process.
ArrayRCP< StridedData< lno_t, scalar_t > > vertexWeights_
typename InputTraits< User >::part_t part_t
int getNumWeightsPerEdge() const override
Returns the number (0 or greater) of edge weights.
void getEdgesDeviceView(typename Base::ConstOffsetsDeviceView &offsets, typename Base::ConstIdsDeviceView &adjIds) const override
Gets adjacency lists for all vertices in a compressed sparse row (CSR) format.
void getEdgeWeightsDeviceView(typename Base::WeightsDeviceView1D &weights, int idx=0) const override
Provide a device view of the edge weights, if any.
void setWeightsHost(typename Base::ConstWeightsHostView1D val, int idx)
Provide a host view of weights for the primary entity type.
void getVertexWeightsView(const scalar_t *&weights, int &stride, int idx) const override
Provide a pointer to the vertex weights, if any.
bool useDegreeAsVertexWeight(int idx) const override
Indicate whether vertex weight with index idx should be the global degree of the vertex.
void applyPartitioningSolution(const User &in, User *&out, const PartitioningSolution< Adapter > &solution) const
typename InputTraits< User >::node_t node_t
A PartitioningSolution is a solution to a partitioning problem.
size_t getLocalNumEdges() const override
Returns the number of edges on this process.
void getEdgesHostView(typename Base::ConstOffsetsHostView &offsets, typename Base::ConstIdsHostView &adjIds) const override
Gets adjacency lists for all vertices in a compressed sparse row (CSR) format.
void getVertexWeightsHostView(typename Base::WeightsHostView1D &weights, int idx=0) const override
Provide a host view of the vertex weights, if any.
void getEdgeWeightsHostView(typename Base::WeightsHostView1D &weights, int idx=0) const override
Provide a host view of the edge weights, if any.
typename InputTraits< User >::gno_t gno_t
Base::ConstIdsDeviceView adjIdsDevice_
Base::VtxDegreeHostView vertexDegreeWeightsHost_
The StridedData class manages lists of weights or coordinates.
void setEdgeWeightsDevice(typename Base::ConstWeightsDeviceView1D val, int idx)
Provide a device view to edge weights.
void setWeightsDevice(typename Base::ConstWeightsDeviceView1D val, int idx)
Provide a device view of weights for the primary entity type.
map_t::local_ordinal_type lno_t
void getVertexIDsDeviceView(typename Base::ConstIdsDeviceView &ids) const override
Sets pointers to this process' graph entries.
typename InputTraits< User >::offset_t offset_t
void setVertexWeightsDevice(typename Base::ConstWeightsDeviceView1D val, int idx)
Provide a device view to vertex weights.
void setVertexWeightIsDegree(int idx)
Specify an index for which the vertex weight should be the degree of the vertex.
TpetraRowGraphAdapter(const RCP< const User > &ingraph, int nVtxWeights=0, int nEdgeWeights=0)
Constructor for graph with no weights or coordinates.
void setEdgeWeights(const scalar_t *val, int stride, int idx)
Provide a pointer to edge weights.
Base::WeightsDeviceView vertexWeightsDevice_
Defines the GraphAdapter interface.
void setWeights(const scalar_t *val, int stride, int idx)
Provide a pointer to weights for the primary entity type.
Base::ConstOffsetsHostView offsHost_
virtual RCP< User > doMigration(const User &from, size_t numLocalRows, const gno_t *myNewRows) const
void setEdgeWeightsHost(typename Base::ConstWeightsHostView1D val, int idx)
Provide a host view to edge weights.
Provides access for Zoltan2 to Tpetra::RowGraph data.
This file defines the StridedData class.