70 #include <Teuchos_Comm.hpp>
71 #include <Teuchos_DefaultComm.hpp>
72 #include <Teuchos_ArrayView.hpp>
78 using Teuchos::ArrayView;
82 typedef Tpetra::CrsMatrix<zscalar_t, zlno_t, zgno_t, znode_t>
tcrsMatrix_t;
83 typedef Tpetra::CrsGraph<zlno_t, zgno_t, znode_t>
tcrsGraph_t;
84 typedef Tpetra::Map<zlno_t, zgno_t, znode_t>
tmap_t;
98 template<
typename offset_t>
102 const RCP<
const Comm<int> > &comm)
104 int rank = comm->getRank();
105 int nprocs = comm->getSize();
108 for (
int p=0; p < nprocs; p++){
110 std::cout <<
"Rank " << p << std::endl;
111 for (
zlno_t i=0; i < nrows; i++){
112 std::cout <<
" Vtx " << i <<
": ";
114 for (offset_t j=idx[i]; j < idx[i+1]; j++)
115 std::cout << *elid++ <<
" ";
117 for (offset_t j=idx[i]; j < idx[i+1]; j++)
118 std::cout << *egid++ <<
" ";
119 std::cout << std::endl;
130 template <
typename MatrixOrGraph>
132 RCP<const MatrixOrGraph> &M,
133 size_t &numLocalDiags,
134 size_t &numGlobalDiags
142 RCP<const tcrsGraph_t> &M,
143 size_t &numLocalDiags,
144 size_t &numGlobalDiags
147 typedef typename tcrsGraph_t::global_ordinal_type
gno_t;
149 size_t maxnnz = M->getLocalMaxNumRowEntries();
150 typename tcrsGraph_t::nonconst_global_inds_host_view_type colGids(
"colGIds", maxnnz);
154 int nLocalRows = M->getLocalNumRows();
155 for (
int i = 0; i < nLocalRows; i++) {
157 gno_t rowGid = M->getRowMap()->getGlobalElement(i);
159 M->getGlobalRowCopy(rowGid, colGids, nnz);
161 for (
size_t j = 0; j < nnz; j++) {
162 if (rowGid == colGids[j]) {
168 Teuchos::reduceAll<int, size_t>(*(M->getComm()), Teuchos::REDUCE_SUM, 1,
169 &numLocalDiags, &numGlobalDiags);
174 RCP<const tcrsMatrix_t> &M,
175 size_t &numLocalDiags,
176 size_t &numGlobalDiags
179 RCP<const tcrsGraph_t> graph = M->getCrsGraph();
185 template <
typename BaseAdapter,
typename Adapter,
typename MatrixOrGraph>
187 RCP<const MatrixOrGraph> &M,
188 RCP<
const Tpetra::CrsGraph<zlno_t, zgno_t> > &Mgraph,
189 const RCP<
const Comm<int> > &comm,
190 bool idsAreConsecutive,
191 int nVtxWeights,
int nEdgeWeights,
int nnzWgtIdx,
int coordDim,
192 bool consecutiveIdsRequested,
bool removeSelfEdges,
bool buildLocalGraph)
199 int rank = comm->getRank();
200 int nprocs = comm->getSize();
203 zlno_t nLocalRows = M->getLocalNumRows();
204 zlno_t nLocalNZ = M->getLocalNumEntries();
205 zgno_t nGlobalRows = M->getGlobalNumRows();
206 zgno_t nGlobalNZ = M->getGlobalNumEntries();
208 std::bitset<Zoltan2::NUM_MODEL_FLAGS> modelFlags;
209 if (consecutiveIdsRequested)
222 for (
int i=0; i < coordDim; i++){
224 for (
zlno_t j=0; j < nLocalRows; j++){
225 coords[i][j] = 100000*i + j;
230 if (nVtxWeights > 0){
231 rowWeights =
new zscalar_t * [nVtxWeights];
232 for (
int i=0; i < nVtxWeights; i++){
234 rowWeights[i] = NULL;
236 rowWeights[i] =
new zscalar_t [nLocalRows];
237 for (
zlno_t j=0; j < nLocalRows; j++){
238 rowWeights[i][j] = 200000*i + j;
244 if (nEdgeWeights > 0 && rank == 0){
245 std::cout <<
"TODO: STILL NEED TO TEST EDGE WEIGHTS!" << std::endl;
250 Adapter tmi(M, nVtxWeights);
252 for (
int i=0; i < nVtxWeights; i++){
254 tmi.setWeightIsDegree(i);
256 tmi.setWeights(rowWeights[i], 1, i);
264 gids =
new zgno_t[nLocalRows];
265 for (
zlno_t i = 0; i < nLocalRows; i++)
266 gids[i] = M->getRowMap()->getGlobalElement(i);
268 (coordDim > 1 ? coords[1] : NULL),
269 (coordDim > 2 ? coords[2] : NULL),
271 tmi.setCoordinateInput(via);
274 size_t numLocalDiags = 0;
275 size_t numGlobalDiags = 0;
276 if (removeSelfEdges) {
277 computeNumDiags<MatrixOrGraph>(M, numLocalDiags, numGlobalDiags);
280 const RCP<const tmap_t> rowMap = M->getRowMap();
281 const RCP<const tmap_t> colMap = M->getColMap();
285 int *numNbors =
new int [nLocalRows];
286 int *numLocalNbors =
new int [nLocalRows];
287 bool *haveDiag =
new bool [nLocalRows];
288 size_t totalLocalNbors = 0;
290 for (
zlno_t i=0; i < nLocalRows; i++){
291 numLocalNbors[i] = 0;
294 typename Tpetra::CrsGraph<zlno_t, zgno_t>::local_inds_host_view_type idx;
295 Mgraph->getLocalRowView(i, idx);
296 numNbors[i] = idx.extent(0);
298 for (std::size_t j=0; j < idx.size(); j++){
305 zgno_t gidVal = colMap->getGlobalElement(idx[j]);
306 if (rowMap->getLocalElement(gidVal) !=
307 Teuchos::OrdinalTraits<zlno_t>::invalid()) {
318 if (rank == 0) std::cout <<
" Creating GraphModel" << std::endl;
320 RCP<const BaseAdapter> baseTmi = rcp(dynamic_cast<BaseAdapter *>(&tmi),
false);
326 catch (std::exception &e){
327 std::cerr << rank <<
") " << e.what() << std::endl;
334 if (rank == 0) std::cout <<
" Checking counts" << std::endl;
338 if (buildLocalGraph) {
342 size_t num = (removeSelfEdges ? totalLocalNbors - numLocalDiags
354 size_t num = (removeSelfEdges ? nLocalNZ-numLocalDiags : nLocalNZ);
358 num = (removeSelfEdges ? (nGlobalNZ-numGlobalDiags) : nGlobalNZ);
373 if (rank == 0) std::cout <<
" Checking vertices" << std::endl;
374 ArrayView<const zgno_t> vertexGids;
375 ArrayView<input_t> crds;
376 ArrayView<input_t> wgts;
382 catch (std::exception &e){
383 std::cerr << rank <<
") Error " << e.what() << std::endl;
388 if (vertexGids.size() != nLocalRows) fail = 1;
391 if (buildLocalGraph) {
393 for (
zlno_t i = 0; i < nLocalRows; i++) {
394 if (vertexGids[i] != i) {
402 if (idsAreConsecutive){
403 zgno_t minLocalGID = rowMap->getMinGlobalIndex();
404 for (
zlno_t i=0; i < nLocalRows; i++){
405 if (vertexGids[i] != minLocalGID + i) {
412 if (consecutiveIdsRequested) {
414 zgno_t tnLocalRows = nLocalRows;
415 scan(*comm, Teuchos::REDUCE_SUM, 1, &tnLocalRows, &myFirstRow);
416 myFirstRow -= nLocalRows;
417 for (
zlno_t i=0; i < nLocalRows; i++){
418 if (vertexGids[i] != myFirstRow+i){
419 std::cout << rank <<
" Row " << i <<
" of " << nLocalRows
420 <<
" myFirstRow+i " << myFirstRow+i
421 <<
" vertexGids " << vertexGids[i]
430 for (
zlno_t i=0; i < nLocalRows; i++, myGid += nprocs){
431 if (vertexGids[i] != myGid){
432 std::cout << rank <<
" Row " << i <<
" of " << nLocalRows
433 <<
" myGid " << myGid <<
" vertexGids " << vertexGids[i]
445 if ((crds.size() != coordDim) || (wgts.size() != nVtxWeights)) fail = 1;
448 for (
int i=0; !fail && i < coordDim; i++){
449 for (
zlno_t j=0; j < nLocalRows; j++){
450 if (crds[i][j] != 100000*i + j){
458 for (
int i=0; !fail && i < nVtxWeights; i++){
460 for (
zlno_t j=0; j < nLocalRows; j++){
461 zscalar_t val = (buildLocalGraph ? numLocalNbors[j] : numNbors[j]);
462 if (removeSelfEdges && haveDiag[j])
464 if (wgts[i][j] != val){
471 for (
zlno_t j=0; j < nLocalRows; j++){
472 if (wgts[i][j] != 200000*i + j){
483 if (rank == 0) std::cout <<
" Checking edges" << std::endl;
485 if (!buildLocalGraph) {
486 ArrayView<const zgno_t> edgeGids;
487 ArrayView<const offset_t> offsets;
491 numEdges = model->
getEdgeList(edgeGids, offsets, wgts);
493 catch(std::exception &e){
494 std::cerr << rank <<
") Error " << e.what() << std::endl;
502 for (
typename ArrayView<const offset_t>::size_type i=0;
503 i < offsets.size()-1; i++){
504 offset_t edgeListSize = offsets[i+1] - offsets[i];
506 size_t val = numNbors[i];
507 if (removeSelfEdges && haveDiag[i])
509 if (edgeListSize != val){
518 std::cout <<
"Printing graph now " << nGlobalRows << std::endl;
519 printGraph<offset_t>(nLocalRows, vertexGids.getRawPtr(), NULL,
520 edgeGids.getRawPtr(), offsets.getRawPtr(), comm);
524 std::cout <<
" " << nGlobalRows <<
" total rows" << std::endl;
531 if (rank == 0) std::cout <<
" Checking local edges" << std::endl;
532 ArrayView<const zgno_t> localEdges;
533 ArrayView<const offset_t> localOffsets;
534 size_t numLocalNeighbors=0;
537 numLocalNeighbors= model->
getEdgeList(localEdges, localOffsets, wgts);
539 catch(std::exception &e){
540 std::cout << rank <<
") Error " << e.what() << std::endl;
541 std::cerr << rank <<
") Error " << e.what() << std::endl;
546 "getLocalEdgeList localOffsets.size", 1)
549 for (
zlno_t i=0; i < nLocalRows; i++){
550 offset_t edgeListSize = localOffsets[i+1] - localOffsets[i];
552 size_t val = numLocalNbors[i];
553 if (removeSelfEdges && haveDiag[i])
555 if (edgeListSize != val){
556 std::cout << rank <<
"vtx " << i <<
" of " << localOffsets.size()
557 <<
" Number of local edges in model " << edgeListSize
558 <<
" not equal to expected number of edges " << val
567 "getLocalEdgeList sum size", 1)
569 fail = ((removeSelfEdges ? totalLocalNbors-numLocalDiags
571 != numLocalNeighbors);
575 if (totalLocalNbors == 0){
577 std::cout <<
" Graph of local edges is empty" << std::endl;
580 printGraph<offset_t>(nLocalRows, vertexGids.getRawPtr(),
581 localEdges.getRawPtr(), NULL, localOffsets.getRawPtr(), comm);
586 if (rank == 0) std::cout <<
" Cleaning up" << std::endl;
591 delete [] numLocalNbors;
594 if (nVtxWeights > 0){
595 for (
int i=0; i < nVtxWeights; i++){
597 delete [] rowWeights[i];
599 delete [] rowWeights;
605 for (
int i=0; i < coordDim; i++){
613 if (rank==0) std::cout <<
" OK" << std::endl;
618 const RCP<
const Comm<int> > &comm,
619 int nVtxWeights,
int nnzWgtIdx,
int coordDim,
620 bool consecutiveIdsRequested,
bool removeSelfEdges,
bool buildLocalGraph)
622 int rank = comm->getRank();
625 std::cout << std::endl <<
"=======================" << std::endl;
626 if (fname.size() > 0)
627 std::cout << std::endl <<
"Test parameters: file name " << fname << std::endl;
629 std::cout << std::endl <<
"Test parameters: dimension ";
630 std::cout << xdim <<
"x" << ydim <<
"x" << zdim << std::endl;
633 std::cout <<
"Num Vertex Weights: " << nVtxWeights << std::endl;
635 std::cout <<
" Dimension " << nnzWgtIdx <<
" is number of neighbors" << std::endl;
637 std::cout <<
"Coordinate dim: " << coordDim << std::endl;
638 std::cout <<
"Request consecutive vertex gids: ";
639 std::cout << (consecutiveIdsRequested ?
"yes" :
"no") << std::endl;
640 std::cout <<
"Request to remove self edges: ";
641 std::cout << (removeSelfEdges ?
"yes" :
"no") << std::endl;
647 if (fname.size() > 0)
656 RCP<const tcrsMatrix_t> Mconsec = rcp_const_cast<
const tcrsMatrix_t>(M);
658 RCP<const Tpetra::CrsGraph<zlno_t, zgno_t> > graph = Mconsec->getCrsGraph();
664 std::cout <<
" TEST MatrixAdapter for graph having Consecutive IDs"
666 bool idsAreConsecutive =
true;
668 testAdapter<baseMAdapter_t,xMAdapter_t,tcrsMatrix_t>(Mconsec, graph, comm,
672 consecutiveIdsRequested,
677 std::cout <<
" TEST GraphAdapter for graph having Consecutive IDs"
679 testAdapter<baseGAdapter_t,xGAdapter_t,tcrsGraph_t>(graph, graph, comm,
683 consecutiveIdsRequested,
689 Array<zgno_t> myNewRows;
690 int nprocs = comm->getSize();
691 for (
size_t i=rank; i < Mconsec->getGlobalNumRows(); i+=nprocs)
692 myNewRows.push_back(i);
694 RCP<const tcrsMatrix_t> Mnonconsec = rcp_const_cast<
const tcrsMatrix_t>(
696 *Mconsec, myNewRows.size(), myNewRows.getRawPtr()));
698 graph = Mnonconsec->getCrsGraph();
704 std::cout <<
" TEST MatrixAdapter graph having Round-Robin IDs"
706 idsAreConsecutive =
false;
708 testAdapter<baseMAdapter_t,xMAdapter_t,tcrsMatrix_t>(Mnonconsec, graph, comm,
712 consecutiveIdsRequested,
717 std::cout <<
" TEST GraphAdapter graph having Round-Robin IDs"
719 testAdapter<baseGAdapter_t,xGAdapter_t,tcrsGraph_t>(graph, graph, comm,
723 consecutiveIdsRequested,
731 int main(
int narg,
char *arg[])
733 Tpetra::ScopeGuard tscope(&narg, &arg);
734 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
736 int rank = comm->getRank();
741 bool consecutiveIdsRequested =
false;
742 bool removeSelfEdges =
false;
743 bool buildLocalGraph =
false;
744 string fname(
"simple");
747 std::cout <<
"TESTING base case (global)" << std::endl;
749 nVtxWeights, nnzWgtIdx, coordDim,
750 consecutiveIdsRequested, removeSelfEdges, buildLocalGraph);
753 std::cout <<
"TESTING with row weights (global)" << std::endl;
756 nVtxWeights, nnzWgtIdx, coordDim,
757 consecutiveIdsRequested, removeSelfEdges, buildLocalGraph);
760 std::cout <<
"TESTING with weights = nnz (global)" << std::endl;
763 nVtxWeights, nnzWgtIdx, coordDim,
764 consecutiveIdsRequested, removeSelfEdges, buildLocalGraph);
767 std::cout <<
"TESTING with multiple row weights and coords (global)"
772 nVtxWeights, nnzWgtIdx, coordDim,
773 consecutiveIdsRequested, removeSelfEdges, buildLocalGraph);
776 std::cout <<
"TESTING with consecutiveIdsRequested (global)" << std::endl;
777 consecutiveIdsRequested =
true;
779 nVtxWeights, nnzWgtIdx, coordDim,
780 consecutiveIdsRequested, removeSelfEdges, buildLocalGraph);
783 std::cout <<
"TESTING with consecutiveIdsRequested and removeSelfEdges "
784 <<
"(global)" << std::endl;
785 removeSelfEdges =
true;
787 nVtxWeights, nnzWgtIdx, coordDim,
788 consecutiveIdsRequested, removeSelfEdges, buildLocalGraph);
791 std::cout <<
"TESTING with removeSelfEdges (global)" << std::endl;
792 consecutiveIdsRequested =
false;
794 nVtxWeights, nnzWgtIdx, coordDim,
795 consecutiveIdsRequested, removeSelfEdges, buildLocalGraph);
798 std::cout <<
"TESTING base case (local)" << std::endl;
799 buildLocalGraph =
true;
800 consecutiveIdsRequested =
false;
801 removeSelfEdges =
false;
803 nVtxWeights, nnzWgtIdx, coordDim,
804 consecutiveIdsRequested, removeSelfEdges, buildLocalGraph);
807 std::cout <<
"TESTING with row weights (local)" << std::endl;
810 nVtxWeights, nnzWgtIdx, coordDim,
811 consecutiveIdsRequested, removeSelfEdges, buildLocalGraph);
814 std::cout <<
"TESTING with weights = nnz (local)" << std::endl;
817 nVtxWeights, nnzWgtIdx, coordDim,
818 consecutiveIdsRequested, removeSelfEdges, buildLocalGraph);
821 std::cout <<
"TESTING with multiple row weights and coords (local)"
826 nVtxWeights, nnzWgtIdx, coordDim,
827 consecutiveIdsRequested, removeSelfEdges, buildLocalGraph);
830 std::cout <<
"TESTING with consecutiveIdsRequested (local)" << std::endl;
831 consecutiveIdsRequested =
true;
833 nVtxWeights, nnzWgtIdx, coordDim,
834 consecutiveIdsRequested, removeSelfEdges, buildLocalGraph);
837 std::cout <<
"TESTING with consecutiveIdsRequested and removeSelfEdges"
840 removeSelfEdges =
true;
842 nVtxWeights, nnzWgtIdx, coordDim,
843 consecutiveIdsRequested, removeSelfEdges, buildLocalGraph);
846 std::cout <<
"TESTING with removeSelfEdges (local)" << std::endl;
847 consecutiveIdsRequested =
false;
849 nVtxWeights, nnzWgtIdx, coordDim,
850 consecutiveIdsRequested, removeSelfEdges, buildLocalGraph);
852 std::cout <<
"PASS" << std::endl;
size_t getLocalNumEdges() const
Returns the number of edges on this process. In global or subset graphs, includes off-process edges...
Zoltan2::GraphAdapter< tcrsGraph_t, simpleUser_t > baseGAdapter_t
size_t getVertexList(ArrayView< const gno_t > &Ids, ArrayView< input_t > &wgts) const
Sets pointers to this process' vertex Ids and their weights.
int getNumWeightsPerVertex() const
Returns the number (0 or greater) of weights per vertex.
void testAdapter(RCP< const MatrixOrGraph > &M, RCP< const Tpetra::CrsGraph< zlno_t, zgno_t > > &Mgraph, const RCP< const Comm< int > > &comm, bool idsAreConsecutive, int nVtxWeights, int nEdgeWeights, int nnzWgtIdx, int coordDim, bool consecutiveIdsRequested, bool removeSelfEdges, bool buildLocalGraph)
MatrixAdapter defines the adapter interface for matrices.
void computeNumDiags< tcrsMatrix_t >(RCP< const tcrsMatrix_t > &M, size_t &numLocalDiags, size_t &numGlobalDiags)
Provides access for Zoltan2 to Xpetra::CrsMatrix data.
GraphAdapter defines the interface for graph-based user data.
algorithm requires consecutive ids
void computeNumDiags< tcrsGraph_t >(RCP< const tcrsGraph_t > &M, size_t &numLocalDiags, size_t &numGlobalDiags)
map_t::global_ordinal_type gno_t
#define TEST_FAIL_AND_EXIT(comm, ok, s, code)
static RCP< User > doMigration(const User &from, size_t numLocalRows, const gno_t *myNewRows)
Migrate the object Given a user object and a new row distribution, create and return a new user objec...
Provides access for Zoltan2 to Xpetra::CrsGraph data.
Zoltan2::BasicUserTypes< zscalar_t, zlno_t, zgno_t > simpleUser_t
int main(int narg, char **arg)
common code used by tests
size_t getGlobalNumVertices() const
Returns the global number vertices.
void computeNumDiags(RCP< const MatrixOrGraph > &M, size_t &numLocalDiags, size_t &numGlobalDiags)
Zoltan2::XpetraCrsGraphAdapter< tcrsGraph_t, simpleUser_t > xGAdapter_t
Defines XpetraCrsGraphAdapter class.
void testGraphModel(string fname, zgno_t xdim, zgno_t ydim, zgno_t zdim, const RCP< const Comm< int > > &comm, int nVtxWeights, int nnzWgtIdx, int coordDim, bool consecutiveIdsRequested, bool removeSelfEdges, bool buildLocalGraph)
algorithm requires no self edges
Defines the XpetraCrsMatrixAdapter class.
int getCoordinateDim() const
Returns the dimension (0 to 3) of vertex coordinates.
Zoltan2::BasicVectorAdapter< simpleUser_t > simpleVAdapter_t
size_t getVertexCoords(ArrayView< input_t > &xyz) const
Sets pointers to this process' vertex coordinates, if available. Order of coordinate info matches tha...
size_t getEdgeList(ArrayView< const gno_t > &edgeIds, ArrayView< const offset_t > &offsets, ArrayView< input_t > &wgts) const
Sets pointers to this process' 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 ...
The StridedData class manages lists of weights or coordinates.
Zoltan2::MatrixAdapter< tcrsMatrix_t, simpleUser_t > baseMAdapter_t
BasicVectorAdapter represents a vector (plus optional weights) supplied by the user as pointers to st...
The user parameters, debug, timing and memory profiling output objects, and error checking methods...
const int SMALL_NUMBER_OF_ROWS
Test of GraphModel interface.
size_t getLocalNumVertices() const
Returns the number vertices on this process.
Tpetra::Map::local_ordinal_type zlno_t
static const std::string fail
Tpetra::CrsMatrix< zscalar_t, zlno_t, zgno_t, znode_t > tcrsMatrix_t
Tpetra::Map< zlno_t, zgno_t, znode_t > tmap_t
Tpetra::CrsGraph< zlno_t, zgno_t, znode_t > tcrsGraph_t
GraphModel defines the interface required for graph models.
int getNumWeightsPerEdge() const
Returns the number (0 or greater) of weights per edge.
Defines the GraphModel interface.
Zoltan2::XpetraCrsMatrixAdapter< tcrsMatrix_t, simpleUser_t > xMAdapter_t
Defines the BasicVectorAdapter class.
model represents graph within only one rank
Tpetra::Map::global_ordinal_type zgno_t
std::string testDataFilePath(".")