Zoltan2 Version of the Day
GraphModel.cpp
Go to the documentation of this file.
00001 // @HEADER
00002 //
00003 // ***********************************************************************
00004 //
00005 //   Zoltan2: A package of combinatorial algorithms for scientific computing
00006 //                  Copyright 2012 Sandia Corporation
00007 //
00008 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00009 // the U.S. Government retains certain rights in this software.
00010 //
00011 // Redistribution and use in source and binary forms, with or without
00012 // modification, are permitted provided that the following conditions are
00013 // met:
00014 //
00015 // 1. Redistributions of source code must retain the above copyright
00016 // notice, this list of conditions and the following disclaimer.
00017 //
00018 // 2. Redistributions in binary form must reproduce the above copyright
00019 // notice, this list of conditions and the following disclaimer in the
00020 // documentation and/or other materials provided with the distribution.
00021 //
00022 // 3. Neither the name of the Corporation nor the names of the
00023 // contributors may be used to endorse or promote products derived from
00024 // this software without specific prior written permission.
00025 //
00026 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00027 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00028 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00029 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00030 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00031 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00032 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00033 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00034 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00035 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00036 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00037 //
00038 // Questions? Contact Karen Devine      (kddevin@sandia.gov)
00039 //                    Erik Boman        (egboman@sandia.gov)
00040 //                    Siva Rajamanickam (srajama@sandia.gov)
00041 //
00042 // ***********************************************************************
00043 //
00044 // @HEADER
00045 //
00046 // Testing of GraphModel built from Xpetra matrix input adapters.
00047 //
00048 
00058 #include <Zoltan2_GraphModel.hpp>
00059 #include <Zoltan2_XpetraCrsMatrixInput.hpp>
00060 #include <Zoltan2_TestHelpers.hpp>
00061 
00062 #include <string>
00063 #include <vector>
00064 #include <iostream>
00065 #include <bitset>
00066 
00067 #include <Teuchos_Comm.hpp>
00068 #include <Teuchos_DefaultComm.hpp>
00069 #include <Teuchos_ArrayView.hpp>
00070 
00071 using namespace std;
00072 using Teuchos::RCP;
00073 using Teuchos::rcp;
00074 using Teuchos::Comm;
00075 using Teuchos::DefaultComm;
00076 using Teuchos::ArrayView;
00077 
00078 typedef Tpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t> tcrsMatrix_t;
00079 typedef Tpetra::Map<lno_t, gno_t, node_t> tmap_t;
00080 typedef Zoltan2::StridedData<lno_t, scalar_t> input_t;
00081 
00082 using std::string;
00083 using std::vector;
00084 
00085 void printGraph(lno_t nrows, const gno_t *v, const lno_t *elid, 
00086     const gno_t *egid,
00087     const int *owner, const lno_t *idx, const RCP<const Comm<int> > &comm)
00088 {
00089   int rank = comm->getRank();
00090   int nprocs = comm->getSize();
00091   comm->barrier();
00092   if (rank == 0){
00093     if (owner)
00094       std::cout << "Global graph:" << std::endl;
00095     else
00096       std::cout << "Local graph:" << std::endl;
00097   }
00098 
00099   for (int p=0; p < nprocs; p++){
00100     if (p == rank){
00101       std::cout << "Rank " << p << std::endl;
00102       if (owner){
00103         for (lno_t i=0; i < nrows; i++){
00104           std::cout << "  Vtx " << *v++ << ": ";
00105           if (elid)
00106             for (lno_t j=idx[i]; j < idx[i+1]; j++)
00107               std::cout << *elid++ << " (" << *owner++ << ") ";
00108           else
00109             for (lno_t j=idx[i]; j < idx[i+1]; j++)
00110               std::cout << *egid++ << " (" << *owner++ << ") ";
00111           std::cout << std::endl;
00112         }
00113         std::cout.flush();
00114       }
00115       else{
00116         for (lno_t i=0; i < nrows; i++){
00117           std::cout << "  Vtx " << i << ": ";
00118           if (elid)
00119             for (lno_t j=idx[i]; j < idx[i+1]; j++)
00120               std::cout << *elid++ << " ";
00121           else
00122             for (lno_t j=idx[i]; j < idx[i+1]; j++)
00123               std::cout << *egid++ << " ";
00124           std::cout << std::endl;
00125         }
00126         std::cout.flush();
00127       }
00128     }
00129     comm->barrier();
00130   }
00131   comm->barrier();
00132 }
00133 
00134 void testMatrixInput(RCP<const Tpetra::CrsMatrix<scalar_t, lno_t, gno_t> > &M,
00135     const RCP<const Comm<int> > &comm,
00136     bool idsAreConsecutive,
00137     int rowWeightDim, int nnzDim, int coordDim,
00138     bool consecutiveIdsRequested, bool removeSelfEdges)
00139 {
00140   int fail=0;
00141   int rank = comm->getRank();
00142   int nprocs = comm->getSize();
00143   RCP<const Zoltan2::Environment> env = rcp(new Zoltan2::Environment);
00144 
00145   lno_t nLocalRows = M->getNodeNumRows();
00146   lno_t nLocalNZ = M->getNodeNumEntries();
00147   gno_t nGlobalRows =  M->getGlobalNumRows();
00148   gno_t nGlobalNZ = M->getGlobalNumEntries();
00149 
00150   // Create a matrix input adapter.
00151 
00152   typedef Zoltan2::MatrixInput<tcrsMatrix_t> base_adapter_t;
00153   typedef Zoltan2::XpetraCrsMatrixInput<tcrsMatrix_t> adapter_t;
00154 
00155   std::bitset<Zoltan2::NUM_MODEL_FLAGS> modelFlags;
00156   if (consecutiveIdsRequested)
00157     modelFlags.set(Zoltan2::IDS_MUST_BE_GLOBALLY_CONSECUTIVE);
00158   if (removeSelfEdges)
00159     modelFlags.set(Zoltan2::SELF_EDGES_MUST_BE_REMOVED);
00160 
00161   scalar_t **coords=NULL;
00162   scalar_t **rowWeights=NULL;
00163 
00164   if (coordDim > 0){
00165     coords = new scalar_t * [coordDim];
00166     for (int i=0; i < coordDim; i++){
00167       coords[i] = new scalar_t [nLocalRows];
00168       for (lno_t j=0; j < nLocalRows; j++){
00169         coords[i][j] = 100000*i + j;
00170       }
00171     }
00172   }
00173 
00174   if (rowWeightDim > 0){
00175     rowWeights = new scalar_t * [rowWeightDim];
00176     for (int i=0; i < rowWeightDim; i++){
00177       if (nnzDim == i)
00178         rowWeights[i] = NULL;
00179       else{
00180         rowWeights[i] = new scalar_t [nLocalRows];
00181         for (lno_t j=0; j < nLocalRows; j++){
00182           rowWeights[i][j] = 200000*i + j;
00183         }
00184       }
00185     }
00186   }
00187 
00188   adapter_t tmi(M, coordDim, rowWeightDim);
00189 
00190   for (int i=0; i < rowWeightDim; i++){
00191     if (nnzDim == i)
00192       tmi.setRowWeightIsNumberOfNonZeros(i);
00193     else
00194       tmi.setRowWeights(i, rowWeights[i], 1);
00195   }
00196 
00197   for (int i=0; i < coordDim; i++){
00198     tmi.setRowCoordinates(i, coords[i], 1);
00199   }
00200 
00201   int numLocalDiags = M->getNodeNumDiags();
00202   int numGlobalDiags = M->getGlobalNumDiags();
00203 
00204   const RCP<const tmap_t> rowMap = M->getRowMap();
00205   const RCP<const tmap_t> colMap = M->getColMap();
00206 
00207   // How many neighbor vertices are not on my process.
00208 
00209   // we know GIDs are consecutive
00210   int minLocalGID = rowMap->getMinGlobalIndex();
00211   int maxLocalGID = rowMap->getMaxGlobalIndex();
00212 
00213   int *numNbors = new int [nLocalRows];
00214   int *numLocalNbors = new int [nLocalRows];
00215   bool *haveDiag = new bool [nLocalRows];
00216   gno_t totalLocalNbors = 0;
00217 
00218   for (lno_t i=0; i < nLocalRows; i++){
00219     numLocalNbors[i] = 0;
00220     haveDiag[i] = false;
00221     ArrayView<const lno_t> idx;
00222     ArrayView<const scalar_t> val;
00223     M->getLocalRowView(i, idx, val);
00224     numNbors[i] = idx.size();
00225 
00226     for (lno_t j=0; j < idx.size(); j++){
00227       if (idx[j] == i){
00228         haveDiag[i] = true;
00229       }
00230       else{
00231         gno_t gidVal = colMap->getGlobalElement(idx[j]);
00232         if (gidVal >= minLocalGID && gidVal <= maxLocalGID){
00233           numLocalNbors[i]++;
00234           totalLocalNbors++;
00235         }
00236       }
00237     }
00238   }
00239 
00240   // Create a GraphModel based on this input data.
00241 
00242   Zoltan2::GraphModel<base_adapter_t> *model = NULL;
00243   const base_adapter_t *baseTmi = &tmi;
00244 
00245   try{
00246     model = new Zoltan2::GraphModel<base_adapter_t>(baseTmi, env, 
00247       comm, modelFlags);
00248   }
00249   catch (std::exception &e){
00250     std::cerr << rank << ") " << e.what() << std::endl;
00251     fail = 1;
00252   }
00253   TEST_FAIL_AND_EXIT(*comm, !fail, "Creating xpetra graph model", 1)
00254 
00255   // Test the GraphModel interface
00256 
00257   if (model->getLocalNumVertices() != size_t(nLocalRows))
00258     fail = 1;
00259   TEST_FAIL_AND_EXIT(*comm, !fail, "getLocalNumVertices", 1)
00260 
00261   if (model->getGlobalNumVertices() != size_t(nGlobalRows))
00262     fail = 1;
00263   TEST_FAIL_AND_EXIT(*comm, !fail, "getGlobalNumVertices", 1)
00264 
00265   size_t num = (removeSelfEdges ? (nLocalNZ-numLocalDiags) : nLocalNZ);
00266 
00267   if (model->getLocalNumGlobalEdges() != num)
00268     fail = 1;
00269   TEST_FAIL_AND_EXIT(*comm, !fail, "getLocalNumGlobalEdges", 1)
00270 
00271   num = totalLocalNbors;
00272   if (!removeSelfEdges){
00273     for (lno_t i=0; i < nLocalRows; i++)
00274       if (haveDiag[i])
00275         num++;
00276   }
00277 
00278   if (model->getLocalNumLocalEdges() != num)
00279     fail = 1;
00280   TEST_FAIL_AND_EXIT(*comm, !fail, "getLocalNumLocalEdges", 1)
00281 
00282   num = (removeSelfEdges ? (nGlobalNZ-numGlobalDiags) : nGlobalNZ);
00283 
00284   if (model->getGlobalNumEdges() != num)
00285     fail = 1;
00286   TEST_FAIL_AND_EXIT(*comm, !fail, "getGlobalNumEdges", 1)
00287 
00288   if (model->getVertexWeightDim() != rowWeightDim)
00289     fail = 1;
00290   TEST_FAIL_AND_EXIT(*comm, !fail, "getVertexWeightDim", 1)
00291 
00292   if (model->getEdgeWeightDim() != 0)
00293     fail = 1;
00294   TEST_FAIL_AND_EXIT(*comm, !fail, "getEdgeWeightDim", 1)
00295 
00296   if (model->getCoordinateDim() != coordDim)
00297     fail = 1;
00298   TEST_FAIL_AND_EXIT(*comm, !fail, "getCoordinateDim", 1)
00299 
00300   ArrayView<const gno_t> vertexGids;
00301   ArrayView<input_t> crds;
00302   ArrayView<input_t> wgts;
00303 
00304   try{
00305     model->getVertexList(vertexGids, crds, wgts);
00306   }
00307   catch (std::exception &e){
00308     std::cerr << rank << ") Error " << e.what() << std::endl;
00309     fail = 1;
00310   }
00311   TEST_FAIL_AND_EXIT(*comm, !fail, "getVertexList", 1)
00312 
00313   if (vertexGids.size() != nLocalRows)
00314     fail = 1;
00315   TEST_FAIL_AND_EXIT(*comm, !fail, "getVertexList", 1)
00316 
00317   // We know model stores things in same order we gave it.
00318 
00319   if (idsAreConsecutive){
00320     for (lno_t i=0; i < nLocalRows; i++){
00321       if (vertexGids[i] != minLocalGID + i) {
00322         fail = 1;
00323         break;
00324       }
00325     }
00326   }
00327   else{  // round robin ids
00328     gno_t myGid = rank;
00329     for (lno_t i=0; i < nLocalRows; i++, myGid += nprocs){
00330       if (vertexGids[i] != myGid){
00331         fail = 1;
00332         break;
00333       }
00334     }
00335   }
00336 
00337   TEST_FAIL_AND_EXIT(*comm, !fail, "vertex gids", 1)
00338 
00339   if ((crds.size() != coordDim) || (wgts.size() != rowWeightDim))
00340     fail = 1;
00341   TEST_FAIL_AND_EXIT(*comm, !fail, "coord or weight array size", 1)
00342 
00343   for (int i=0; !fail && i < coordDim; i++){
00344     for (lno_t j=0; j < nLocalRows; j++){
00345       if (crds[i][j] != 200000*i + j){
00346         fail = 1;
00347         break;
00348       }
00349     }
00350   }
00351   TEST_FAIL_AND_EXIT(*comm, !fail, "coord values", 1)
00352 
00353   for (int i=0; !fail && i < rowWeightDim; i++){
00354     if (nnzDim == i){
00355       for (lno_t j=0; j < nLocalRows; j++){
00356         scalar_t val = numNbors[j];
00357         if (removeSelfEdges && haveDiag[j])
00358           val -= 1;
00359         if (wgts[i][j] != val){
00360           fail = 1;
00361           break;
00362         }
00363       }
00364     }
00365     else{
00366       for (lno_t j=0; j < nLocalRows; j++){
00367         if (wgts[i][j] != 100000*i + j){
00368           fail = 1;
00369           break;
00370         }
00371       }
00372     }
00373   }
00374 
00375   TEST_FAIL_AND_EXIT(*comm, !fail, "row weight values", 1)
00376   
00377   ArrayView<const gno_t> edgeGids;
00378   ArrayView<const int> procIds;
00379   ArrayView<const lno_t> offsets;
00380   size_t numEdges=0;
00381 
00382   try{
00383     numEdges = model->getEdgeList(edgeGids, procIds, offsets, wgts);
00384   }
00385   catch(std::exception &e){
00386     std::cerr << rank << ") Error " << e.what() << std::endl;
00387     fail = 1;
00388   }
00389   TEST_FAIL_AND_EXIT(*comm, !fail, "getEdgeList", 1)
00390 
00391   TEST_FAIL_AND_EXIT(*comm, wgts.size() == 0, "edge weights present", 1)
00392 
00393   num = 0;
00394   for (ArrayView<const lno_t>::size_type i=0; i < offsets.size()-1; i++){
00395     size_t edgeListSize = offsets[i+1] - offsets[i];
00396     num += edgeListSize;
00397     size_t val = numNbors[i];
00398     if (removeSelfEdges && haveDiag[i])
00399       val--;
00400     if (edgeListSize != val){
00401       fail = 1;
00402       break;
00403     }
00404   }
00405 
00406   TEST_FAIL_AND_EXIT(*comm, numEdges==num, "getEdgeList size", 1)
00407 
00408   if (nGlobalRows < 200){
00409     printGraph(nLocalRows, vertexGids.getRawPtr(), NULL,
00410       edgeGids.getRawPtr(), procIds.getRawPtr(), offsets.getRawPtr(), comm);
00411   }
00412   else{
00413     if (rank==0) 
00414       std::cout << "    " << nGlobalRows << " total rows" << std::endl;
00415   }
00416 
00417   // Get graph restricted to this process
00418 
00419   ArrayView<const lno_t> localEdges;
00420   ArrayView<const lno_t> localOffsets;
00421   size_t numLocalNeighbors=0;
00422 
00423   try{
00424     numLocalNeighbors= model->getLocalEdgeList(localEdges, localOffsets, wgts);
00425   }
00426   catch(std::exception &e){
00427     std::cerr << rank << ") Error " << e.what() << std::endl;
00428     fail = 1;
00429   }
00430   TEST_FAIL_AND_EXIT(*comm, !fail, "getLocalEdgeList", 1)
00431 
00432 #ifdef THIS_CODE_DOESNT_FAIL
00433   num = 0;
00434   for (size_t i=0; i < localOffsets.size()-1; i++){
00435     size_t edgeListSize = localOffsets[i+1] - localOffsets[i];
00436     num += edgeListSize;
00437     size_t val = numLocalNbors[i];
00438     if (removeSelfEdges && haveDiag[i])
00439       val--;
00440     if (edgeListSize != val){
00441       fail = 1;
00442       break;
00443     }
00444   }
00445   TEST_FAIL_AND_EXIT(*comm, !fail, "getLocalEdgeList list sizes", 1)
00446 
00447   TEST_FAIL_AND_EXIT(*comm, numLocalNeighbors==num, "getLocalEdgeList size", 1)
00448 
00449   if (size_t(totalLocalNbors) != numLocalNeighbors)
00450     fail = 1;
00451 
00452   TEST_FAIL_AND_EXIT(*comm, !fail, "getLocalEdgeList size", 1)
00453 
00454   if (nGlobalRows < 200){
00455     if (totalLocalNbors == 0){
00456       if (rank == 0)
00457         std::cout << "  Graph of local edges is empty" << std::endl; 
00458     }
00459     else{
00460       printGraph(nLocalRows, vertexGids.getRawPtr(), 
00461         localEdges.getRawPtr(), NULL, NULL, localOffsets.getRawPtr(), comm);
00462     }
00463   }
00464 #endif
00465 
00466   delete model;
00467 
00468   if (nLocalRows){
00469     delete [] numNbors;
00470     delete [] numLocalNbors;
00471     delete [] haveDiag;
00472 
00473     if (rowWeightDim > 0){
00474       for (int i=0; i < rowWeightDim; i++){
00475         if (rowWeights[i])
00476           delete [] rowWeights[i];
00477       }
00478       delete [] rowWeights;
00479     }
00480 
00481     if (coordDim > 0){
00482       for (int i=0; i < coordDim; i++){
00483         if (coords[i])
00484           delete [] coords[i];
00485       }
00486       delete [] coords;
00487     }
00488   }
00489 
00490   if (rank==0) std::cout << "    OK" << std::endl;
00491 }
00492 
00493 void testGraphModel(string fname, gno_t xdim, gno_t ydim, gno_t zdim,
00494     const RCP<const Comm<int> > &comm,
00495     int rowWeightDim, int nnzDim, int coordDim,
00496     bool consecutiveIdsRequested, bool removeSelfEdges)
00497 {
00498   int rank = comm->getRank();
00499 
00500   if (rank==0){
00501     cout << endl << "=======================" << endl;
00502     if (fname.size() > 0)
00503       cout << endl << "Test parameters: file name " << fname << endl;
00504     else{
00505       cout << endl << "Test parameters: dimension ";
00506       cout  << xdim << "x" << ydim << "x" << zdim << endl;
00507     }
00508 
00509     cout << "Vertex weight dim: " << rowWeightDim << endl;
00510     if (nnzDim >= 0)
00511      cout << "  Dimension " << nnzDim << " is number of neighbors" << endl;
00512 
00513     cout << "Coordinate dim: " << coordDim << endl;
00514     cout << "Request consecutive vertex gids: ";
00515     cout << (consecutiveIdsRequested ? "yes" : "no") << endl;
00516     cout << "Request to remove self edges: ";
00517     cout << (removeSelfEdges ? "yes" : "no") << endl;
00518   }
00519 
00520   typedef Tpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t> tcrsMatrix_t;
00521 
00522   // Input generator
00523   UserInputForTests *input;
00524 
00525   if (fname.size() > 0)
00526     input = new UserInputForTests(testDataFilePath, fname, comm, true);
00527   else
00528     input = new UserInputForTests(xdim,ydim,zdim,string(""), comm, true);
00529 
00530   RCP<tcrsMatrix_t> M = input->getTpetraCrsMatrix();
00531 
00532   // Row Ids of test input are already consecutive
00533 
00534   RCP<const tcrsMatrix_t> Mconsec = rcp_const_cast<const tcrsMatrix_t>(M);
00535 
00536   RCP<const Tpetra::CrsGraph<lno_t, gno_t> > graph = Mconsec->getCrsGraph();
00537 
00538   printTpetraGraph<lno_t, gno_t>(comm, *graph, cout, 100, 
00539     "Graph with consecutive IDs");
00540 
00541   bool idsAreConsecutive = true;
00542 
00543   testMatrixInput(Mconsec, comm,  idsAreConsecutive,
00544     rowWeightDim, nnzDim, coordDim,
00545     consecutiveIdsRequested, removeSelfEdges);
00546 
00547 #if 0
00548   testGraphInput(Mconsec, comm, idsAreConsecutive,
00549     rowWeightDim, nnzDim, coordDim,
00550     consecutiveIdsRequested, removeSelfEdges);
00551 #endif
00552 
00553   // Do a round robin migration so that global IDs are not consecutive.
00554 
00555 #ifdef TODO_THESE_HAVE_BEEN_TESTED
00556   Array<gno_t> myNewRows;
00557   for (size_t i=rank; i < Mconsec->getGlobalNumRows(); i+=nprocs)
00558     myNewRows.push_back(i);
00559 
00560   RCP<const tcrsMatrix_t> Mnonconsec = 
00561     Zoltan2::XpetraTraits<tcrsMatrix_t>::doMigration(
00562       Mconsec, myNewRows.size(), myNewRows.getRawPtr());
00563 
00564   graph = Mnonconsec->getCrsGraph();
00565 
00566   printTpetraGraph<lno_t, gno_t>(comm, *graph, cout, 100, 
00567     "Graph with non-consecutive IDs");
00568 
00569   idsAreConsecutive = false;
00570 
00571   testMatrixInput(Mnonconsec, comm, idsAreConsecutive,
00572     rowWeightDim, nnzDim, coordDim,
00573     consecutiveIdsRequested, removeSelfEdges);
00574 
00575 #if 0
00576   testGraphInput(Mnonconsec, comm, idsAreConsecutive,
00577      rowWeightDim, nnzDim, coordDim,
00578     consecutiveIdsRequested, removeSelfEdges);
00579 #endif
00580 
00581 #endif
00582 
00583   delete input;
00584 }
00585 
00586 int main(int argc, char *argv[])
00587 {
00588   Teuchos::GlobalMPISession session(&argc, &argv);
00589   Teuchos::RCP<const Teuchos::Comm<int> > comm =
00590     Teuchos::DefaultComm<int>::getComm();
00591 
00592   int rank = comm->getRank();
00593 
00594   int rowWeightDim=0;
00595   int nnzDim = -1; 
00596   int coordDim=0;
00597   bool consecutiveIdsRequested=false, removeSelfEdges=false;
00598   string fname("simple");
00599 
00600   testGraphModel(fname, 0, 0, 0, comm,
00601     rowWeightDim, nnzDim, coordDim,
00602     consecutiveIdsRequested, removeSelfEdges);
00603 
00604 #ifdef TODO_THESE_HAVE_BEEN_TESTED
00605   rowWeightDim = 1;
00606 
00607   testGraphModel(fname, 0, 0, 0, comm,
00608     rowWeightDim, nnzDim, coordDim,
00609     consecutiveIdsRequested, removeSelfEdges);
00610 
00611   nnzDim = 1;
00612 
00613   testGraphModel(fname, 0, 0, 0, comm,
00614     rowWeightDim, nnzDim, coordDim,
00615     consecutiveIdsRequested, removeSelfEdges);
00616 
00617   rowWeightDim = 2;
00618   coordDim = 3;
00619 
00620   testGraphModel(fname, 0, 0, 0, comm,
00621     rowWeightDim, nnzDim, coordDim,
00622     consecutiveIdsRequested, removeSelfEdges);
00623 
00624   consecutiveIdsRequested = true;
00625 
00626   testGraphModel(fname, 0, 0, 0, comm,
00627     rowWeightDim, nnzDim, coordDim,
00628     consecutiveIdsRequested, removeSelfEdges);
00629 
00630   removeSelfEdges = true;
00631 
00632   testGraphModel(fname, 0, 0, 0, comm,
00633     rowWeightDim, nnzDim, coordDim,
00634     consecutiveIdsRequested, removeSelfEdges);
00635 
00636   consecutiveIdsRequested = false;
00637 
00638   testGraphModel(fname, 0, 0, 0, comm,
00639     rowWeightDim, nnzDim, coordDim,
00640     consecutiveIdsRequested, removeSelfEdges);
00641 #endif
00642 
00643   if (rank==0)
00644     cout << "PASS" << endl;
00645 
00646   return 0;
00647 }
00648