Zoltan2 Version of the Day
XpetraCrsGraphInput.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 // Basic testing of Zoltan2::XpetraCrsGraphInput
00052 #include <string>
00053 
00054 #include <Zoltan2_XpetraCrsGraphInput.hpp>
00055 #include <Zoltan2_PartitioningSolution.hpp>
00056 #include <Zoltan2_TestHelpers.hpp>
00057 
00058 #include <Teuchos_GlobalMPISession.hpp>
00059 #include <Teuchos_DefaultComm.hpp>
00060 #include <Teuchos_RCP.hpp>
00061 #include <Teuchos_Comm.hpp>
00062 #include <Teuchos_CommHelpers.hpp>
00063 
00064 using namespace std;
00065 using Teuchos::RCP;
00066 using Teuchos::rcp;
00067 using Teuchos::rcp_const_cast;
00068 using Teuchos::Comm;
00069 using Teuchos::DefaultComm;
00070 using Teuchos::Array;
00071 using Teuchos::ArrayView;
00072 
00073 typedef UserInputForTests uinput_t;
00074 typedef Tpetra::CrsGraph<lno_t, gno_t, node_t> tgraph_t;
00075 typedef Xpetra::CrsGraph<lno_t, gno_t, node_t> xgraph_t;
00076 typedef Epetra_CrsGraph egraph_t;
00077 
00078 void printGraph(RCP<const Comm<int> > &comm, lno_t nvtx,
00079     const gno_t *vtxIds, const lno_t *offsets, const gno_t *edgeIds)
00080 {
00081   int rank = comm->getRank();
00082   int nprocs = comm->getSize();
00083   comm->barrier();
00084   for (int p=0; p < nprocs; p++){
00085     if (p == rank){
00086       std::cout << rank << ":" << std::endl;
00087       for (lno_t i=0; i < nvtx; i++){
00088         std::cout << " vertex " << vtxIds[i] << ": ";
00089         for (lno_t j=offsets[i]; j < offsets[i+1]; j++){
00090           std::cout << edgeIds[j] << " ";
00091         }
00092         std::cout << std::endl;
00093       }
00094       std::cout.flush();
00095     }
00096     comm->barrier();
00097   }
00098   comm->barrier();
00099 }
00100 
00101 template <typename User>
00102 int verifyInputAdapter(
00103   Zoltan2::XpetraCrsGraphInput<User> &ia, tgraph_t &graph)
00104 {
00105   RCP<const Comm<int> > comm = graph.getComm();
00106   int fail = 0, gfail=0;
00107 
00108   if (!fail && ia.getLocalNumberOfVertices() != graph.getNodeNumRows())
00109     fail = 4;
00110 
00111   if (!fail && ia.getLocalNumberOfEdges() != graph.getNodeNumEntries())
00112       fail = 6;
00113 
00114   gfail = globalFail(comm, fail);
00115 
00116   const gno_t *vtxIds=NULL, *edgeIds=NULL;
00117   const lno_t *offsets=NULL;
00118   size_t nvtx=0;
00119 
00120   if (!gfail){
00121 
00122     nvtx = ia.getVertexListView(vtxIds, offsets, edgeIds);
00123 
00124     if (nvtx != graph.getNodeNumRows())
00125       fail = 8;
00126 
00127     gfail = globalFail(comm, fail);
00128 
00129     if (gfail == 0){
00130       printGraph(comm, nvtx, vtxIds, offsets, edgeIds);
00131     }
00132     else{
00133       if (!fail) fail = 10;
00134     }
00135   }
00136   return fail;
00137 }
00138 
00139 int main(int argc, char *argv[])
00140 {
00141   Teuchos::GlobalMPISession session(&argc, &argv);
00142   RCP<const Comm<int> > comm = DefaultComm<int>::getComm();
00143   int rank = comm->getRank();
00144   int fail = 0, gfail=0;
00145 
00146   // Create an object that can give us test Tpetra, Xpetra
00147   // and Epetra graphs for testing.
00148 
00149   RCP<uinput_t> uinput;
00150 
00151   try{
00152     uinput =
00153       rcp(new uinput_t(testDataFilePath,std::string("simple"), comm, true));
00154   }
00155   catch(std::exception &e){
00156     TEST_FAIL_AND_EXIT(*comm, 0, string("input ")+e.what(), 1);
00157   }
00158 
00159   RCP<tgraph_t> tG;     // original graph (for checking)
00160   RCP<tgraph_t> newG;   // migrated graph
00161 
00162   tG = uinput->getTpetraCrsGraph();
00163   size_t nvtx = tG->getNodeNumRows();
00164   ArrayView<const gno_t> rowGids = tG->getRowMap()->getNodeElementList();
00165 
00166   // To test migration in the input adapter we need a Solution
00167   // object.  The Solution needs an IdentifierMap.
00168   // Our solution just assigns all objects to part zero.
00169 
00170   typedef Zoltan2::IdentifierMap<tgraph_t> idmap_t;
00171 
00172   RCP<const Zoltan2::Environment> env = rcp(new Zoltan2::Environment);
00173 
00174   ArrayRCP<const gno_t> gidArray = arcpFromArrayView(rowGids);
00175   RCP<const idmap_t> idMap = rcp(new idmap_t(env, comm, gidArray));
00176 
00177   int weightDim = 1;
00178 
00179   zoltan2_partId_t *p = new zoltan2_partId_t [nvtx];
00180   memset(p, 0, sizeof(zoltan2_partId_t) * nvtx);
00181   ArrayRCP<zoltan2_partId_t> solnParts(p, 0, nvtx, true);
00182 
00183   typedef Zoltan2::XpetraCrsGraphInput<tgraph_t>  adapter_t;
00184   typedef Zoltan2::PartitioningSolution<adapter_t> soln_t;
00185   soln_t solution(env, comm, idMap, weightDim);
00186   solution.setParts(gidArray, solnParts, true);
00187 
00189   // User object is Tpetra::CrsGraph
00190   if (!gfail){
00191     RCP<const tgraph_t> ctG = rcp_const_cast<const tgraph_t>(tG);
00192     RCP<Zoltan2::XpetraCrsGraphInput<tgraph_t> > tGInput;
00193 
00194     try {
00195       tGInput =
00196         rcp(new Zoltan2::XpetraCrsGraphInput<tgraph_t>(ctG));
00197     }
00198     catch (std::exception &e){
00199       TEST_FAIL_AND_EXIT(*comm, 0,
00200         string("XpetraCrsGraphInput ")+e.what(), 1);
00201     }
00202 
00203     if (rank==0)
00204       std::cout << "Input adapter for Tpetra::CrsGraph" << std::endl;
00205 
00206     fail = verifyInputAdapter<tgraph_t>(*tGInput, *tG);
00207 
00208     gfail = globalFail(comm, fail);
00209 
00210     if (!gfail){
00211       tgraph_t *mMigrate = NULL;
00212       try{
00213         tGInput->applyPartitioningSolution( *tG, mMigrate, solution);
00214         newG = rcp(mMigrate);
00215       }
00216       catch (std::exception &e){
00217         fail = 11;
00218       }
00219 
00220       gfail = globalFail(comm, fail);
00221 
00222       if (!gfail){
00223         RCP<const tgraph_t> cnewG = rcp_const_cast<const tgraph_t>(newG);
00224         RCP<Zoltan2::XpetraCrsGraphInput<tgraph_t> > newInput;
00225         try{
00226           newInput = rcp(new Zoltan2::XpetraCrsGraphInput<tgraph_t>(cnewG));
00227         }
00228         catch (std::exception &e){
00229           TEST_FAIL_AND_EXIT(*comm, 0,
00230             string("XpetraCrsGraphInput 2 ")+e.what(), 1);
00231         }
00232 
00233         if (rank==0){
00234           std::cout <<
00235            "Input adapter for Tpetra::CrsGraph migrated to proc 0" <<
00236            std::endl;
00237         }
00238         fail = verifyInputAdapter<tgraph_t>(*newInput, *newG);
00239         if (fail) fail += 100;
00240         gfail = globalFail(comm, fail);
00241       }
00242     }
00243     if (gfail){
00244       printFailureCode(comm, fail);
00245     }
00246   }
00247 
00249   // User object is Xpetra::CrsGraph
00250   if (!gfail){
00251     RCP<xgraph_t> xG = uinput->getXpetraCrsGraph();
00252     RCP<const xgraph_t> cxG = rcp_const_cast<const xgraph_t>(xG);
00253     RCP<Zoltan2::XpetraCrsGraphInput<xgraph_t> > xGInput;
00254 
00255     try {
00256       xGInput =
00257         rcp(new Zoltan2::XpetraCrsGraphInput<xgraph_t>(cxG));
00258     }
00259     catch (std::exception &e){
00260       TEST_FAIL_AND_EXIT(*comm, 0,
00261         string("XpetraCrsGraphInput 3 ")+e.what(), 1);
00262     }
00263 
00264     if (rank==0){
00265       std::cout << "Input adapter for Xpetra::CrsGraph" << std::endl;
00266     }
00267     fail = verifyInputAdapter<xgraph_t>(*xGInput, *tG);
00268 
00269     gfail = globalFail(comm, fail);
00270 
00271     if (!gfail){
00272       xgraph_t *mMigrate =NULL;
00273       try{
00274         xGInput->applyPartitioningSolution( *xG, mMigrate, solution);
00275       }
00276       catch (std::exception &e){
00277         fail = 11;
00278       }
00279 
00280       gfail = globalFail(comm, fail);
00281 
00282       if (!gfail){
00283         RCP<const xgraph_t> cnewG(mMigrate);
00284         RCP<Zoltan2::XpetraCrsGraphInput<xgraph_t> > newInput;
00285         try{
00286           newInput =
00287             rcp(new Zoltan2::XpetraCrsGraphInput<xgraph_t>(cnewG));
00288         }
00289         catch (std::exception &e){
00290           TEST_FAIL_AND_EXIT(*comm, 0,
00291             string("XpetraCrsGraphInput 4 ")+e.what(), 1);
00292         }
00293 
00294         if (rank==0){
00295           std::cout <<
00296            "Input adapter for Xpetra::CrsGraph migrated to proc 0" <<
00297            std::endl;
00298         }
00299         fail = verifyInputAdapter<xgraph_t>(*newInput, *newG);
00300         if (fail) fail += 100;
00301         gfail = globalFail(comm, fail);
00302       }
00303     }
00304     if (gfail){
00305       printFailureCode(comm, fail);
00306     }
00307   }
00308 
00309 #ifdef HAVE_EPETRA_DATA_TYPES
00310 
00311   // User object is Epetra_CrsGraph
00312   if (!gfail){
00313     RCP<egraph_t> eG = uinput->getEpetraCrsGraph();
00314     RCP<const egraph_t> ceG = rcp_const_cast<const egraph_t>(eG);
00315     RCP<Zoltan2::XpetraCrsGraphInput<egraph_t> > eGInput;
00316 
00317     try {
00318       eGInput =
00319         rcp(new Zoltan2::XpetraCrsGraphInput<egraph_t>(ceG));
00320     }
00321     catch (std::exception &e){
00322       TEST_FAIL_AND_EXIT(*comm, 0,
00323         string("XpetraCrsGraphInput 5 ")+e.what(), 1);
00324     }
00325 
00326     if (rank==0){
00327       std::cout << "Input adapter for Epetra_CrsGraph" << std::endl;
00328     }
00329     fail = verifyInputAdapter<egraph_t>(*eGInput, *tG);
00330 
00331     gfail = globalFail(comm, fail);
00332 
00333     if (!gfail){
00334       egraph_t *mMigrate =NULL;
00335       try{
00336         eGInput->applyPartitioningSolution( *eG, mMigrate, solution);
00337       }
00338       catch (std::exception &e){
00339         fail = 11;
00340       }
00341 
00342       gfail = globalFail(comm, fail);
00343 
00344       if (!gfail){
00345         RCP<const egraph_t> cnewG(mMigrate, true);
00346         RCP<Zoltan2::XpetraCrsGraphInput<egraph_t> > newInput;
00347         try{
00348           newInput =
00349             rcp(new Zoltan2::XpetraCrsGraphInput<egraph_t>(cnewG));
00350         }
00351         catch (std::exception &e){
00352           TEST_FAIL_AND_EXIT(*comm, 0,
00353             string("XpetraCrsGraphInput 6 ")+e.what(), 1);
00354         }
00355 
00356         if (rank==0){
00357           std::cout <<
00358            "Input adapter for Epetra_CrsGraph migrated to proc 0" <<
00359            std::endl;
00360         }
00361         fail = verifyInputAdapter<egraph_t>(*newInput, *newG);
00362         if (fail) fail += 100;
00363         gfail = globalFail(comm, fail);
00364       }
00365     }
00366     if (gfail){
00367       printFailureCode(comm, fail);
00368     }
00369   }
00370 #endif
00371 
00373   // DONE
00374 
00375   if (rank==0)
00376     std::cout << "PASS" << std::endl;
00377 }