Amesos2 - Direct Sparse Solver Interfaces Version of the Day
Amesos2_TpetraMultiVecAdapter_def.hpp
Go to the documentation of this file.
00001 // @HEADER
00002 //
00003 // ***********************************************************************
00004 //
00005 //           Amesos2: Templated Direct Sparse Solver Package
00006 //                  Copyright 2011 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 Michael A. Heroux (maherou@sandia.gov)
00039 //
00040 // ***********************************************************************
00041 //
00042 // @HEADER
00043 
00053 #ifndef AMESOS2_TPETRA_MULTIVEC_ADAPTER_DEF_HPP
00054 #define AMESOS2_TPETRA_MULTIVEC_ADAPTER_DEF_HPP
00055 
00056 #include "Amesos2_TpetraMultiVecAdapter_decl.hpp"
00057 
00058 
00059 namespace Amesos2 {
00060 
00061   using Tpetra::MultiVector;
00062 
00063   template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, class Node >
00064   MultiVecAdapter<
00065     MultiVector<Scalar,
00066                 LocalOrdinal,
00067                 GlobalOrdinal,
00068                 Node> >::MultiVecAdapter( const Teuchos::RCP<multivec_t>& m )
00069   : mv_(m)
00070   {}
00071 
00072 
00073   template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, class Node >
00074   void
00075   MultiVecAdapter<
00076     MultiVector<Scalar,
00077                 LocalOrdinal,
00078                 GlobalOrdinal,
00079                 Node> >::get1dCopy(const Teuchos::ArrayView<scalar_t>& av,
00080                                    size_t lda,
00081                                    Teuchos::Ptr<
00082                                      const Tpetra::Map<LocalOrdinal,
00083                                                        GlobalOrdinal,
00084                                                        Node> > distribution_map ) const
00085   {
00086     using Teuchos::as;
00087     using Teuchos::RCP;
00088     typedef Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> map_type;
00089     const size_t num_vecs = getGlobalNumVectors ();
00090 
00091     TEUCHOS_TEST_FOR_EXCEPTION(
00092       distribution_map.getRawPtr () == NULL, std::invalid_argument,
00093       "Amesos2::MultiVecAdapter::get1dCopy: distribution_map argument is null.");
00094     TEUCHOS_TEST_FOR_EXCEPTION(
00095       mv_.is_null (), std::logic_error,
00096       "Amesos2::MultiVecAdapter::get1dCopy: mv_ is null.");
00097     // Check mv_ before getMap(), because the latter calls mv_->getMap().
00098     TEUCHOS_TEST_FOR_EXCEPTION(
00099       this->getMap ().is_null (), std::logic_error,
00100       "Amesos2::MultiVecAdapter::get1dCopy: this->getMap() returns null.");
00101 
00102 #ifdef HAVE_AMESOS2_DEBUG
00103     const size_t requested_vector_length = distribution_map->getNodeNumElements ();
00104     TEUCHOS_TEST_FOR_EXCEPTION(
00105       lda < requested_vector_length, std::invalid_argument,
00106       "Amesos2::MultiVecAdapter::get1dCopy: On process " <<
00107       distribution_map->getComm ()->getRank () << " of the distribution Map's "
00108       "communicator, the given stride lda = " << lda << " is not large enough "
00109       "for the local vector length " << requested_vector_length << ".");
00110     TEUCHOS_TEST_FOR_EXCEPTION(
00111       as<size_t> (av.size ()) < as<size_t> ((num_vecs - 1) * lda + requested_vector_length),
00112       std::invalid_argument, "Amesos2::MultiVector::get1dCopy: MultiVector "
00113       "storage not large enough given leading dimension and number of vectors." );
00114 #endif // HAVE_AMESOS2_DEBUG
00115 
00116     // (Re)compute the Export object if necessary.  If not, then we
00117     // don't need to clone distribution_map; we can instead just get
00118     // the previously cloned target Map from the Export object.
00119     RCP<const map_type> distMap;
00120     if (exporter_.is_null () ||
00121         ! exporter_->getSourceMap ()->isSameAs (* (this->getMap ())) ||
00122         ! exporter_->getTargetMap ()->isSameAs (* distribution_map)) {
00123 
00124       // Since we're caching the Export object, and since the Export
00125       // needs to keep the distribution Map, we have to make a copy of
00126       // the latter in order to ensure that it will stick around past
00127       // the scope of this function call.  (Ptr is not reference
00128       // counted.)  Map's clone() method suffices, even though it only
00129       // makes a shallow copy of some of Map's data, because Map is
00130       // immutable and those data are reference-counted (e.g.,
00131       // ArrayRCP or RCP).
00132       distMap = distribution_map->template clone<Node> (distribution_map->getNode ());
00133 
00134       // (Re)create the Export object.
00135       exporter_ = rcp (new export_type (this->getMap (), distMap));
00136     }
00137     else {
00138       distMap = exporter_->getTargetMap ();
00139     }
00140 
00141     multivec_t redist_mv (distMap, num_vecs);
00142 
00143     // Redistribute the input (multi)vector.
00144     redist_mv.doExport (*mv_, *exporter_, Tpetra::REPLACE);
00145 
00146     // Copy the imported (multi)vector's data into the ArrayView.
00147     redist_mv.get1dCopy (av, lda);
00148   }
00149 
00150   template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, class Node >
00151   Teuchos::ArrayRCP<Scalar>
00152   MultiVecAdapter<
00153     MultiVector<Scalar,
00154                 LocalOrdinal,
00155                 GlobalOrdinal,
00156                 Node> >::get1dViewNonConst (bool local)
00157   {
00158     // FIXME (mfh 22 Jan 2014) When I first found this routine, all of
00159     // its code was commented out, and it didn't return anything.  The
00160     // latter is ESPECIALLY dangerous, given that the return value is
00161     // an ArrayRCP.  Thus, I added the exception throw below.
00162     TEUCHOS_TEST_FOR_EXCEPTION(
00163       true, std::logic_error, "Amesos2::MultiVecAdapter::get1dViewNonConst: "
00164       "Not implemented.");
00165 
00166     // if( local ){
00167     //   this->localize();
00168     //   /* Use the global element list returned by
00169     //    * mv_->getMap()->getNodeElementList() to get a subCopy of mv_ which we
00170     //    * assign to l_l_mv_, then return get1dViewNonConst() of l_l_mv_
00171     //    */
00172     //   if(l_l_mv_.is_null() ){
00173     //  Teuchos::ArrayView<const GlobalOrdinal> nodeElements_go
00174     //    = mv_->getMap()->getNodeElementList();
00175     //  Teuchos::Array<size_t> nodeElements_st(nodeElements_go.size());
00176 
00177     //  // Convert the node element to a list of size_t type objects
00178     //  typename Teuchos::ArrayView<const GlobalOrdinal>::iterator it_go;
00179     //  Teuchos::Array<size_t>::iterator it_st = nodeElements_st.begin();
00180     //  for( it_go = nodeElements_go.begin(); it_go != nodeElements_go.end(); ++it_go ){
00181     //    *(it_st++) = Teuchos::as<size_t>(*it_go);
00182     //  }
00183 
00184     //  // To be consistent with the globalize() function, get a view of the local mv
00185     //  l_l_mv_ = l_mv_->subViewNonConst(nodeElements_st);
00186 
00187     //  return(l_l_mv_->get1dViewNonConst());
00188     //   } else {
00189     //  // We need to re-import values to the local, since changes may have been
00190     //  // made to the global structure that are not reflected in the local
00191     //  // view.
00192     //  Teuchos::ArrayView<const GlobalOrdinal> nodeElements_go
00193     //    = mv_->getMap()->getNodeElementList();
00194     //  Teuchos::Array<size_t> nodeElements_st(nodeElements_go.size());
00195 
00196     //  // Convert the node element to a list of size_t type objects
00197     //  typename Teuchos::ArrayView<const GlobalOrdinal>::iterator it_go;
00198     //  Teuchos::Array<size_t>::iterator it_st = nodeElements_st.begin();
00199     //  for( it_go = nodeElements_go.begin(); it_go != nodeElements_go.end(); ++it_go ){
00200     //    *(it_st++) = Teuchos::as<size_t>(*it_go);
00201     //  }
00202 
00203     //  return l_l_mv_->get1dViewNonConst();
00204     //   }
00205     // } else {
00206     //   if( mv_->isDistributed() ){
00207     //  this->localize();
00208 
00209     //  return l_mv_->get1dViewNonConst();
00210     //   } else {                      // not distributed, no need to import
00211     //  return mv_->get1dViewNonConst();
00212     //   }
00213     // }
00214   }
00215 
00216 
00217   template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, class Node>
00218   void
00219   MultiVecAdapter<
00220     MultiVector<Scalar,
00221                 LocalOrdinal,
00222                 GlobalOrdinal,
00223                 Node> >::put1dData(const Teuchos::ArrayView<const scalar_t>& new_data,
00224                                    size_t lda,
00225                                    Teuchos::Ptr<
00226                                      const Tpetra::Map<LocalOrdinal,
00227                                                        GlobalOrdinal,
00228                                                        Node> > source_map)
00229   {
00230     using Teuchos::RCP;
00231     typedef Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> map_type;
00232 
00233     TEUCHOS_TEST_FOR_EXCEPTION(
00234       source_map.getRawPtr () == NULL, std::invalid_argument,
00235       "Amesos2::MultiVecAdapter::put1dData: source_map argument is null.");
00236     TEUCHOS_TEST_FOR_EXCEPTION(
00237       mv_.is_null (), std::logic_error,
00238       "Amesos2::MultiVecAdapter::put1dData: the internal MultiVector mv_ is null.");
00239     // getMap() calls mv_->getMap(), so test first whether mv_ is null.
00240     TEUCHOS_TEST_FOR_EXCEPTION(
00241       this->getMap ().is_null (), std::logic_error,
00242       "Amesos2::MultiVecAdapter::put1dData: this->getMap() returns null.");
00243 
00244     const size_t num_vecs = getGlobalNumVectors ();
00245 
00246     // (Re)compute the Import object if necessary.  If not, then we
00247     // don't need to clone source_map; we can instead just get the
00248     // previously cloned source Map from the Import object.
00249     RCP<const map_type> srcMap;
00250     if (importer_.is_null () ||
00251         ! importer_->getSourceMap ()->isSameAs (* source_map) ||
00252         ! importer_->getTargetMap ()->isSameAs (* (this->getMap ()))) {
00253 
00254       // Since we're caching the Import object, and since the Import
00255       // needs to keep the source Map, we have to make a copy of the
00256       // latter in order to ensure that it will stick around past the
00257       // scope of this function call.  (Ptr is not reference counted.)
00258       // Map's clone() method suffices, even though it only makes a
00259       // shallow copy of some of Map's data, because Map is immutable
00260       // and those data are reference-counted (e.g., ArrayRCP or RCP).
00261       srcMap = source_map->template clone<Node> (source_map->getNode ());
00262       importer_ = rcp (new import_type (srcMap, this->getMap ()));
00263     }
00264     else {
00265       srcMap = importer_->getSourceMap ();
00266     }
00267 
00268     const multivec_t source_mv (srcMap, new_data, lda, num_vecs);
00269 
00270     // Redistribute the output (multi)vector.
00271     mv_->doImport (source_mv, *importer_, Tpetra::REPLACE);
00272   }
00273 
00274 
00275   template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, class Node >
00276   std::string
00277   MultiVecAdapter<
00278     MultiVector<Scalar,
00279                 LocalOrdinal,
00280                 GlobalOrdinal,
00281                 Node> >::description() const
00282   {
00283     std::ostringstream oss;
00284     oss << "Amesos2 adapter wrapping: ";
00285     oss << mv_->description();
00286     return oss.str();
00287   }
00288 
00289 
00290   template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, class Node >
00291   void
00292   MultiVecAdapter<
00293     MultiVector<Scalar,
00294                 LocalOrdinal,
00295                 GlobalOrdinal,
00296                 Node> >::describe (Teuchos::FancyOStream& os,
00297                                    const Teuchos::EVerbosityLevel verbLevel) const
00298   {
00299     mv_->describe (os, verbLevel);
00300   }
00301 
00302 
00303   template <typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, class Node >
00304   const char* MultiVecAdapter<
00305     MultiVector<Scalar,
00306                 LocalOrdinal,
00307                 GlobalOrdinal,
00308                 Node> >::name = "Amesos2 adapter for Tpetra::MultiVector";
00309 
00310 } // end namespace Amesos2
00311 
00312 #endif // AMESOS2_TPETRA_MULTIVEC_ADAPTER_DEF_HPP