Zoltan2 Version of the Day
XpetraVectorInput.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 
00052 #include <string>
00053 
00054 #include <Zoltan2_XpetraVectorInput.hpp>
00055 #include <Zoltan2_InputTraits.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 
00071 typedef UserInputForTests uinput_t;
00072 typedef Tpetra::Vector<scalar_t, lno_t, gno_t, node_t> tvector_t;
00073 typedef Xpetra::Vector<scalar_t, lno_t, gno_t, node_t> xvector_t;
00074 typedef Epetra_Vector evector_t;
00075 
00076 void printVector(RCP<const Comm<int> > &comm, lno_t vlen,
00077     const gno_t *vtxIds, const scalar_t *vals)
00078 {
00079   int rank = comm->getRank();
00080   int nprocs = comm->getSize();
00081   comm->barrier();
00082   for (int p=0; p < nprocs; p++){
00083     if (p == rank){
00084       std::cout << rank << ":" << std::endl;
00085       for (lno_t i=0; i < vlen; i++){
00086         std::cout << " " << vtxIds[i] << ": " << vals[i] << std::endl;
00087       }
00088       std::cout.flush();
00089     }
00090     comm->barrier();
00091   }
00092   comm->barrier();
00093 }
00094 
00095 template <typename User>
00096 int verifyInputAdapter(
00097   Zoltan2::XpetraVectorInput<User> &ia, tvector_t &vector, int wdim, 
00098     scalar_t **weights, int *strides)
00099 {
00100   RCP<const Comm<int> > comm = vector.getMap()->getComm();
00101   int fail = 0, gfail=0;
00102 
00103   if (!fail && ia.getNumberOfVectors() !=1) 
00104     fail = 42;
00105 
00106   if (!fail && ia.getNumberOfWeights() !=wdim) 
00107     fail = 41;
00108 
00109   if (!fail && ia.getLocalLength() != vector.getLocalLength())
00110     fail = 4;
00111 
00112   gfail = globalFail(comm, fail);
00113 
00114   if (!gfail){
00115     const gno_t *vtxIds=NULL;
00116     const scalar_t *vals=NULL;
00117     int stride;
00118 
00119     size_t nvals = ia.getVector(vtxIds, vals, stride);
00120 
00121     if (nvals != vector.getLocalLength())
00122       fail = 8;
00123     if (!fail && stride != 1)
00124       fail = 10;
00125 
00126     gfail = globalFail(comm, fail);
00127 
00128     if (gfail == 0){
00129       printVector(comm, nvals, vtxIds, vals);
00130     }
00131   }
00132 
00133   if (!gfail && wdim){
00134     const scalar_t *wgt =NULL;
00135     int stride;
00136 
00137     for (int w=0; !fail && w < wdim; w++){
00138       size_t nvals = ia.getVectorWeights(w, wgt, stride);
00139 
00140       if (nvals != vector.getLocalLength())
00141         fail = 100;
00142       if (!fail && stride != strides[w])
00143         fail = 101;
00144 
00145       for (size_t v=0; !fail && v < nvals; v++){
00146         if (wgt[v*stride] != weights[w][v*stride])
00147           fail=102;
00148       }
00149     }
00150 
00151     gfail = globalFail(comm, fail);
00152   }
00153 
00154   return gfail;
00155 }
00156 
00157 int main(int argc, char *argv[])
00158 {
00159   Teuchos::GlobalMPISession session(&argc, &argv);
00160   RCP<const Comm<int> > comm = DefaultComm<int>::getComm();
00161   int rank = comm->getRank();
00162   int fail = 0, gfail=0;
00163 
00164   // Create object that can give us test Tpetra, Xpetra
00165   // and Epetra vectors for testing.
00166 
00167   RCP<uinput_t> uinput;
00168 
00169   try{
00170     uinput = 
00171       rcp(new uinput_t(testDataFilePath,std::string("simple"), comm, true));
00172   }
00173   catch(std::exception &e){
00174     TEST_FAIL_AND_EXIT(*comm, 0, string("input ")+e.what(), 1);
00175   }
00176 
00177   RCP<tvector_t> tV;     // original vector (for checking)
00178   RCP<tvector_t> newV;   // migrated vector
00179 
00180   tV = uinput->getTpetraVector();
00181   size_t vlen = tV->getLocalLength();
00182   Teuchos::ArrayView<const gno_t> rowGids = tV->getMap()->getNodeElementList();
00183 
00184   // To test migration in the input adapter we need a Solution
00185   // object.  The Solution needs an IdentifierMap.
00186 
00187   typedef Zoltan2::IdentifierMap<tvector_t> idmap_t;
00188 
00189   RCP<const Zoltan2::Environment> env = rcp(new Zoltan2::Environment);
00190 
00191   ArrayRCP<const gno_t> gidArray = arcpFromArrayView(rowGids);
00192   RCP<const idmap_t> idMap = rcp(new idmap_t(env, comm, gidArray));
00193 
00194   int weightDim = 1;
00195 
00196   zoltan2_partId_t *p = new zoltan2_partId_t [vlen];
00197   memset(p, 0, sizeof(zoltan2_partId_t) * vlen);
00198   ArrayRCP<zoltan2_partId_t> solnParts(p, 0, vlen, true);
00199 
00200   std::vector<const scalar_t *> emptyWeights;
00201   std::vector<int> emptyStrides;
00202 
00203   typedef Zoltan2::XpetraVectorInput<tvector_t> adapter_t;
00204   Zoltan2::PartitioningSolution<adapter_t> solution(
00205     env, comm, idMap, weightDim);
00206   solution.setParts(gidArray, solnParts, true);
00207 
00209   // User object is Tpetra::Vector, no weights
00210   if (!gfail){ 
00211     RCP<const tvector_t> ctV = rcp_const_cast<const tvector_t>(tV);
00212     RCP<adapter_t> tVInput;
00213   
00214     try {
00215       tVInput = 
00216         rcp(new Zoltan2::XpetraVectorInput<tvector_t>(ctV, 
00217           emptyWeights, emptyStrides));
00218     }
00219     catch (std::exception &e){
00220       TEST_FAIL_AND_EXIT(*comm, 0, 
00221         string("XpetraVectorInput ")+e.what(), 1);
00222     }
00223   
00224     if (rank==0){
00225       std::cout << tVInput->inputAdapterName() << ", constructed with ";
00226       std::cout  << "Tpetra::Vector" << std::endl;
00227     }
00228     
00229     fail = verifyInputAdapter<tvector_t>(*tVInput, *tV, 0, NULL, NULL);
00230   
00231     gfail = globalFail(comm, fail);
00232   
00233     if (!gfail){
00234       tvector_t *vMigrate = NULL;
00235       try{
00236         tVInput->applyPartitioningSolution(*tV, vMigrate, solution);
00237         newV = rcp(vMigrate);
00238       }
00239       catch (std::exception &e){
00240         fail = 11;
00241       }
00242 
00243       gfail = globalFail(comm, fail);
00244   
00245       if (!gfail){
00246         RCP<const tvector_t> cnewV = rcp_const_cast<const tvector_t>(newV);
00247         RCP<Zoltan2::XpetraVectorInput<tvector_t> > newInput;
00248         try{
00249           newInput = rcp(new Zoltan2::XpetraVectorInput<tvector_t>(cnewV,
00250             emptyWeights, emptyStrides));
00251         }
00252         catch (std::exception &e){
00253           TEST_FAIL_AND_EXIT(*comm, 0, 
00254             string("XpetraVectorInput 2 ")+e.what(), 1);
00255         }
00256   
00257         if (rank==0){
00258           std::cout << tVInput->inputAdapterName() << ", constructed with ";
00259           std::cout << "Tpetra::Vector migrated to proc 0" << std::endl;
00260         }
00261         fail = verifyInputAdapter<tvector_t>(*newInput, *newV, 0, NULL, NULL);
00262         if (fail) fail += 100;
00263         gfail = globalFail(comm, fail);
00264       }
00265     }
00266     if (gfail){
00267       printFailureCode(comm, fail);
00268     }
00269   }
00270 
00272   // User object is Xpetra::Vector
00273   if (!gfail){ 
00274     RCP<xvector_t> xV = uinput->getXpetraVector();
00275     RCP<const xvector_t> cxV = rcp_const_cast<const xvector_t>(xV);
00276     RCP<Zoltan2::XpetraVectorInput<xvector_t> > xVInput;
00277   
00278     try {
00279       xVInput = 
00280         rcp(new Zoltan2::XpetraVectorInput<xvector_t>(cxV,
00281           emptyWeights, emptyStrides));
00282     }
00283     catch (std::exception &e){
00284       TEST_FAIL_AND_EXIT(*comm, 0, 
00285         string("XpetraVectorInput 3 ")+e.what(), 1);
00286     }
00287   
00288     if (rank==0){
00289       std::cout << xVInput->inputAdapterName() << ", constructed with ";
00290       std::cout << "Xpetra::Vector" << std::endl;
00291     }
00292     fail = verifyInputAdapter<xvector_t>(*xVInput, *tV, 0, NULL, NULL);
00293   
00294     gfail = globalFail(comm, fail);
00295   
00296     if (!gfail){
00297       xvector_t *vMigrate =NULL;
00298       try{
00299         xVInput->applyPartitioningSolution(*xV, vMigrate, solution);
00300       }
00301       catch (std::exception &e){
00302         fail = 11;
00303       }
00304   
00305       gfail = globalFail(comm, fail);
00306   
00307       if (!gfail){
00308         RCP<const xvector_t> cnewV(vMigrate);
00309         RCP<Zoltan2::XpetraVectorInput<xvector_t> > newInput;
00310         try{
00311           newInput = 
00312             rcp(new Zoltan2::XpetraVectorInput<xvector_t>(cnewV, 
00313               emptyWeights, emptyStrides));
00314         }
00315         catch (std::exception &e){
00316           TEST_FAIL_AND_EXIT(*comm, 0, 
00317             string("XpetraVectorInput 4 ")+e.what(), 1);
00318         }
00319   
00320         if (rank==0){
00321           std::cout << xVInput->inputAdapterName() << ", constructed with ";
00322           std::cout << "Xpetra::Vector migrated to proc 0" << std::endl;
00323         }
00324         fail = verifyInputAdapter<xvector_t>(*newInput, *newV, 0, NULL, NULL);
00325         if (fail) fail += 100;
00326         gfail = globalFail(comm, fail);
00327       }
00328     }
00329     if (gfail){
00330       printFailureCode(comm, fail);
00331     }
00332   }
00333 
00334 #ifdef HAVE_EPETRA_DATA_TYPES
00335 
00336   // User object is Epetra_Vector
00337   if (!gfail){ 
00338     RCP<evector_t> eV = uinput->getEpetraVector();
00339     RCP<const evector_t> ceV = rcp_const_cast<const evector_t>(eV);
00340     RCP<Zoltan2::XpetraVectorInput<evector_t> > eVInput;
00341   
00342     try {
00343       eVInput = 
00344         rcp(new Zoltan2::XpetraVectorInput<evector_t>(ceV,
00345           emptyWeights, emptyStrides));
00346     }
00347     catch (std::exception &e){
00348       TEST_FAIL_AND_EXIT(*comm, 0, 
00349         string("XpetraVectorInput 5 ")+e.what(), 1);
00350     }
00351   
00352     if (rank==0){
00353       std::cout << eVInput->inputAdapterName() << ", constructed with ";
00354       std::cout << "Epetra_Vector" << std::endl;
00355     }
00356     fail = verifyInputAdapter<evector_t>(*eVInput, *tV, 0, NULL, NULL);
00357   
00358     gfail = globalFail(comm, fail);
00359   
00360     if (!gfail){
00361       evector_t *vMigrate =NULL;
00362       try{
00363         eVInput->applyPartitioningSolution(*eV, vMigrate, solution);
00364       }
00365       catch (std::exception &e){
00366         fail = 11;
00367       }
00368   
00369       gfail = globalFail(comm, fail);
00370   
00371       if (!gfail){
00372         RCP<const evector_t> cnewV(vMigrate, true);
00373         RCP<Zoltan2::XpetraVectorInput<evector_t> > newInput;
00374         try{
00375           newInput = 
00376             rcp(new Zoltan2::XpetraVectorInput<evector_t>(cnewV, 
00377               emptyWeights, emptyStrides));
00378         }
00379         catch (std::exception &e){
00380           TEST_FAIL_AND_EXIT(*comm, 0, 
00381             string("XpetraVectorInput 6 ")+e.what(), 1);
00382         }
00383   
00384         if (rank==0){
00385            std::cout << eVInput->inputAdapterName() << ", constructed with ";
00386            std::cout << "Epetra_Vector migrated to proc 0" << std::endl;
00387         }
00388         fail = verifyInputAdapter<evector_t>(*newInput, *newV, 0, NULL, NULL);
00389         if (fail) fail += 100;
00390         gfail = globalFail(comm, fail);
00391       }
00392     }
00393     if (gfail){
00394       printFailureCode(comm, fail);
00395     }
00396   }
00397 #endif
00398 
00400   // DONE
00401 
00402   if (rank==0)
00403     std::cout << "PASS" << std::endl;
00404 }