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