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