Anasazi Version of the Day
TsqrRandomizer_Tpetra_MultiVector.hpp
Go to the documentation of this file.
00001 // @HEADER
00002 // ***********************************************************************
00003 //
00004 //                 Anasazi: Block Eigensolvers Package
00005 //                 Copyright (2010) Sandia Corporation
00006 //
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 //
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
00025 //
00026 // ***********************************************************************
00027 // @HEADER
00028 
00029 #ifndef __TSQR_Trilinos_Randomizer_Tpetra_MultiVector_hpp
00030 #define __TSQR_Trilinos_Randomizer_Tpetra_MultiVector_hpp
00031 
00035 
00036 #include "TsqrTypeAdaptor_Tpetra_MultiVector_SerialNode.hpp"
00037 #ifdef HAVE_KOKKOS_TBB
00038 #  include "TsqrTypeAdaptor_Tpetra_MultiVector_TBBNode.hpp"
00039 #endif // HAVE_KOKKOS_TBB
00040 #include "TsqrCommFactory_Tpetra.hpp"
00041 
00042 #include <stdexcept>
00043 #include <sstream>
00044 
00047 
00048 namespace TSQR {
00049   namespace Trilinos {
00050 
00051     template< class S, class LO, class GO, class NodeType, class Gen >
00052     class TpetraRandomizer :
00053       public Randomizer< S, LO, GO, 
00054        Tpetra::MultiVector< S, LO, GO, NodeType >,
00055        Gen >
00056     {
00057     public:
00058       // mfh 14 Jul 2010: I'm uncomfortable naming this "MV" because
00059       // the base class already has MV defined in it.  I don't want to
00060       // go all circular (MV -> multivector_type -> MV).  TMV here is
00061       // just a way to make the base_type typedef fit on one line.
00062       typedef Tpetra::MultiVector< S, LO, GO, NodeType > TMV;
00063 
00064       // The C++ compiler needs help inheriting typedefs, when both
00065       // the base class and the derived class are templated.
00066       typedef Randomizer< double, int, int, TMV, Gen > base_type;
00067       typedef typename base_type::multivector_type multivector_type;
00068       typedef typename base_type::local_ordinal_type local_ordinal_type;
00069       typedef typename base_type::scalar_type scalar_type;
00070 
00071       typedef typename base_type::normalgen_ptr normalgen_ptr;
00072       typedef typename base_type::comm_ptr comm_ptr;
00073       typedef typename base_type::ordinal_messenger_ptr ordinal_messenger_ptr;
00074       typedef typename base_type::scalar_messenger_ptr scalar_messenger_ptr;
00075 
00085       TpetraRandomizer (const multivector_type& mv,
00086       const normalgen_ptr& pGen)
00087       {
00088   init (mv, pGen);
00089       }
00090 
00091     private:
00092       virtual void
00093       fetchDims (const multivector_type& A,
00094      local_ordinal_type& nrowsLocal, 
00095      local_ordinal_type& ncols, 
00096      local_ordinal_type& LDA) const
00097       {
00098   nrowsLocal = A.getLocalLength();
00099   ncols = A.getNumVectors();
00100   if (nrowsLocal < ncols)
00101     {
00102       std::ostringstream os;
00103       os << "The local component of the input matrix has fewer row"
00104         "s (" << nrowsLocal << ") than columns (" << ncols << ").  "
00105         "TSQR::Trilinos::TpetraRandomizer does not support this case.";
00106       throw std::runtime_error (os.str());
00107     }
00108   if (! A.isConstantStride())
00109     {
00110       // FIXME (mfh 14 June 2010) Storage of A uses nonconstant
00111       // stride internally, but that doesn't necessarily mean we
00112       // can't run TSQR.  It depends on what get1dViewNonConst()
00113       // returns.  If it's copied and packed into a matrix with
00114       // constant stride, then we are free to run TSQR.
00115       std::ostringstream os;
00116       os << "TSQR::Trilinos::TpetraRandomizer does not support "
00117         "Tpetra::MultiVector inputs that do not have constant stride.";
00118       throw std::runtime_error (os.str());
00119     }
00120   LDA = A.getStride();
00121       }
00122 
00123       virtual Teuchos::ArrayRCP< scalar_type > 
00124       fetchNonConstView (multivector_type& A) const
00125       {
00126   // This won't be efficient if the entries of A live in a
00127   // different memory space (e.g., GPU node), but it will be
00128   // correct for any Tpetra::MultiVector type.
00129   return A.get1dViewNonConst();
00130       }
00131 
00132       virtual void
00133       fetchMessengers (const multivector_type& mv,
00134            scalar_messenger_ptr& pScalarMessenger,
00135            ordinal_messenger_ptr& pOrdinalMessenger) const
00136       {
00137   typedef TpetraCommFactory< S, LO, GO, NodeType > comm_factory_type;
00138   comm_ptr pComm = mv.getMap()->getComm();
00139   comm_factory_type factory;
00140   factory.makeMessengers (pComm, pScalarMessenger, pOrdinalMessenger);
00141       }
00142     };
00143   } // namespace Trilinos
00144 } // namespace TSQR
00145 
00146 #endif // __TSQR_Trilinos_Randomizer_Tpetra_MultiVector_hpp
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends