Tpetra_KokkosRefactor_MultiVector_decl.hpp

// @HEADER
// ***********************************************************************
//
//          Tpetra: Templated Linear Algebra Services Package
//                 Copyright (2008) Sandia Corporation
//
// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
// the U.S. Government retains certain rights in this software.
//
// Redistribution and use in source and binary forms, with or without
// modification, are permitted provided that the following conditions are
// met:
//
// 1. Redistributions of source code must retain the above copyright
// notice, this list of conditions and the following disclaimer.
//
// 2. Redistributions in binary form must reproduce the above copyright
// notice, this list of conditions and the following disclaimer in the
// documentation and/or other materials provided with the distribution.
//
// 3. Neither the name of the Corporation nor the names of the
// contributors may be used to endorse or promote products derived from
// this software without specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
//
// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
//
// ************************************************************************
// @HEADER

#ifndef TPETRA_KOKKOS_REFACTOR_MULTIVECTOR_DECL_HPP
#define TPETRA_KOKKOS_REFACTOR_MULTIVECTOR_DECL_HPP

#include <Kokkos_DualView.hpp>
#include <KokkosCompat_ClassicNodeAPI_Wrapper.hpp>
#include <Kokkos_InnerProductSpaceTraits.hpp>
#include <Kokkos_ArithTraits.hpp>
#include <Tpetra_KokkosRefactor_DistObject.hpp>

namespace KokkosClassic {
#ifndef DOXYGEN_SHOULD_SKIP_THIS
  // forward declaration of DefaultArithmetic
  template<class KokkosMultiVectorType>
  class DefaultArithmetic;
#endif // DOXYGEN_SHOULD_SKIP_THIS
} // namespace KokkosClassic

namespace Tpetra {

#ifndef DOXYGEN_SHOULD_SKIP_THIS
  // forward declaration of Vector, needed to prevent circular inclusions
  template<class S, class LO, class GO, class N> class Vector;

  // forward declaration of Map
  template<class LO, class GO, class N> class Map;
#endif // DOXYGEN_SHOULD_SKIP_THIS

  template<class Scalar,
           class LocalOrdinal,
           class GlobalOrdinal,
           class DeviceType>
  class MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType> > :
    public DistObject<Scalar, LocalOrdinal, GlobalOrdinal, Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType> >
  {
  private:
    typedef Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType> Node;
    typedef DistObject<Scalar, LocalOrdinal, GlobalOrdinal, Node> base_type;

  public:


    typedef Scalar scalar_type;
    typedef LocalOrdinal local_ordinal_type;
    typedef GlobalOrdinal global_ordinal_type;
    typedef Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType> node_type;

    typedef typename Kokkos::Details::InnerProductSpaceTraits<Scalar>::dot_type dot_type;

    typedef typename Kokkos::Details::ArithTraits<Scalar>::mag_type mag_type;

    typedef DeviceType device_type;

    typedef typename DeviceType::host_mirror_device_type host_mirror_device_type;

    typedef Kokkos::DualView<scalar_type**, Kokkos::LayoutLeft, device_type> dual_view_type;




    MultiVector ();

    MultiVector (const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& map,
                 size_t NumVectors,
                 bool zeroOut=true);

    MultiVector (const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> &source);

    MultiVector (const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& source,
                 const Teuchos::DataAccess copyOrView);

    MultiVector (const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& map,
                 const Teuchos::ArrayView<const Scalar>& A,
                 const size_t LDA,
                 const size_t NumVectors);

    MultiVector (const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& map,
                 const Teuchos::ArrayView<const Teuchos::ArrayView<const Scalar> >&ArrayOfPtrs,
                 const size_t NumVectors);

    MultiVector (const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& map,
                 const dual_view_type& view);

    MultiVector (const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& map,
                 const dual_view_type& view,
                 const dual_view_type& origView);

    MultiVector (const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& map,
                 const dual_view_type& view,
                 const Teuchos::ArrayView<const size_t>& whichVectors);

    MultiVector (const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& map,
                 const dual_view_type& view,
                 const dual_view_type& origView,
                 const Teuchos::ArrayView<const size_t>& whichVectors);

    template <class Node2>
    Teuchos::RCP<MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node2> >
    clone (const Teuchos::RCP<Node2> &node2) const;

    virtual ~MultiVector();




    void
    replaceGlobalValue (GlobalOrdinal globalRow,
                        size_t vectorIndex,
                        const Scalar &value);

    void
    sumIntoGlobalValue (GlobalOrdinal globalRow,
                        size_t vectorIndex,
                        const Scalar &value);

    void
    replaceLocalValue (LocalOrdinal myRow,
                       size_t vectorIndex,
                       const Scalar &value);

    void
    sumIntoLocalValue (LocalOrdinal myRow,
                       size_t vectorIndex,
                       const Scalar &value);

    void putScalar (const Scalar &value);

    void randomize();

    void replaceMap (const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& map);

    void reduce();

    MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
    operator= (const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& source);




    Teuchos::RCP<MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
    subCopy (const Teuchos::Range1D &colRng) const;

    Teuchos::RCP<MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
    subCopy (const Teuchos::ArrayView<const size_t> &cols) const;

    Teuchos::RCP<const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
    subView (const Teuchos::Range1D &colRng) const;

    Teuchos::RCP<const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
    subView (const Teuchos::ArrayView<const size_t> &cols) const;

    Teuchos::RCP<MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
    subViewNonConst (const Teuchos::Range1D &colRng);

    Teuchos::RCP<MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
    subViewNonConst (const Teuchos::ArrayView<const size_t> &cols);

    Teuchos::RCP<const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
    offsetView (const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& subMap,
                size_t offset) const;

    Teuchos::RCP<MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
    offsetViewNonConst (const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > &subMap,
                        size_t offset);

    Teuchos::RCP<const Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
    getVector (size_t j) const;

    Teuchos::RCP<Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
    getVectorNonConst (size_t j);

    Teuchos::ArrayRCP<const Scalar> getData(size_t j) const;

    Teuchos::ArrayRCP<Scalar> getDataNonConst(size_t j);

    void get1dCopy (Teuchos::ArrayView<Scalar> A, size_t LDA) const;

    void get2dCopy (Teuchos::ArrayView<const Teuchos::ArrayView<Scalar> > ArrayOfPtrs) const;

    Teuchos::ArrayRCP<const Scalar> get1dView() const;

    Teuchos::ArrayRCP<Teuchos::ArrayRCP<const Scalar> > get2dView() const;

    Teuchos::ArrayRCP<Scalar> get1dViewNonConst();

    Teuchos::ArrayRCP<Teuchos::ArrayRCP<Scalar> > get2dViewNonConst();

    KokkosClassic::MultiVector<Scalar,Node> getLocalMV() const;

    TEUCHOS_DEPRECATED KokkosClassic::MultiVector<Scalar,Node>
    getLocalMVNonConst ();

    dual_view_type getDualView () const;

    template<class TargetDeviceType>
    void sync () {
      getDualView ().template sync<TargetDeviceType> ();
    }

    template<class TargetDeviceType>
    void modify () {
      getDualView ().template modify<TargetDeviceType> ();
    }

    template<class TargetDeviceType>
    typename Kokkos::Impl::if_c<
      Kokkos::Impl::is_same<
        typename device_type::memory_space,
        typename TargetDeviceType::memory_space>::value,
      typename dual_view_type::t_dev,
      typename dual_view_type::t_host>::type
    getLocalView () const {
      return getDualView ().template view<TargetDeviceType> ();
    }




    void
    dot (const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A,
         const Teuchos::ArrayView<dot_type>& dots) const;

    template <typename T>
    typename Kokkos::Impl::enable_if< !(Kokkos::Impl::is_same<dot_type, T>::value), void >::type
    dot (const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A,
         const Teuchos::ArrayView<T> &dots) const
    {
      const size_t sz = dots.size ();
      Teuchos::Array<dot_type> dts (sz);
      this->dot (A, dts);
      for (size_t i = 0; i < sz; ++i) {
        dots[i] = dts[i];
      }
    }

    void
    dot (const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A,
         const Kokkos::View<dot_type*, device_type>& dots) const;

    template <typename T>
    typename Kokkos::Impl::enable_if< !(Kokkos::Impl::is_same<dot_type, T>::value), void >::type
    dot (const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A,
         const Kokkos::View<T*, device_type>& dots) const
    {
      const size_t numDots = dots.dimension_0 ();
      Kokkos::View<dot_type*, device_type> dts ("MV::dot tmp", numDots);
      // Call overload that takes a Kokkos::View<dot_type*, device_type>.
      this->dot (A, dts);
      // FIXME (mfh 14 Jul 2014) Does this actually work if dot_type
      // and T differ?  We would need a test for this, but only the
      // Sacado and Stokhos packages are likely to care about this use
      // case.
      Kokkos::deep_copy (dots, dts);
    }

    void abs(const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> &A);

    void reciprocal(const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> &A);

    void scale (const Scalar &alpha);

    void scale (Teuchos::ArrayView<const Scalar> alpha);

    void scale (const Kokkos::View<const Scalar*, device_type> alpha);

    void
    scale (const Scalar& alpha,
           const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A);

    void
    update (const Scalar& alpha,
            const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A,
            const Scalar& beta);

    void
    update (const Scalar& alpha,
            const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A,
            const Scalar& beta,
            const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& B,
            const Scalar& gamma);

    void norm1 (const Kokkos::View<mag_type*, device_type>& norms) const;

    template <typename T>
    typename Kokkos::Impl::enable_if< !(Kokkos::Impl::is_same<mag_type, T>::value), void >::type
    norm1 (const Kokkos::View<T*, device_type>& norms) const
    {
      const size_t numNorms = norms.dimension_0 ();
      Kokkos::View<mag_type*, device_type> tmpNorms ("MV::norm1 tmp", numNorms);
      // Call overload that takes a Kokkos::View<mag_type*, device_type>.
      this->norm1 (tmpNorms);
      // FIXME (mfh 15 Jul 2014) Does this actually work if mag_type
      // and T differ?  We would need a test for this, but only the
      // Sacado and Stokhos packages are likely to care about this use
      // case.
      Kokkos::deep_copy (norms, tmpNorms);
    }


    void norm1 (const Teuchos::ArrayView<mag_type>& norms) const;

    template <typename T>
    typename Kokkos::Impl::enable_if< !(Kokkos::Impl::is_same<mag_type,T>::value), void >::type
    norm1 (const Teuchos::ArrayView<T>& norms) const
    {
      typedef typename Teuchos::ArrayView<T>::size_type size_type;
      const size_type sz = norms.size ();
      Teuchos::Array<mag_type> theNorms (sz);
      this->norm1 (theNorms);
      for (size_type i = 0; i < sz; ++i) {
        norms[i] = theNorms[i];
      }
    }

    void norm2 (const Kokkos::View<mag_type*, device_type>& norms) const;

    template<typename T>
    typename Kokkos::Impl::enable_if< !(Kokkos::Impl::is_same<mag_type, T>::value), void >::type
    norm2 (const Kokkos::View<T*, device_type>& norms) const
    {
      const size_t numNorms = norms.dimension_0 ();
      Kokkos::View<mag_type*, device_type> theNorms ("MV::norm2 tmp", numNorms);
      // Call overload that takes a Kokkos::View<mag_type*, device_type>.
      this->norm2 (theNorms);
      // FIXME (mfh 14 Jul 2014) Does this actually work if mag_type
      // and T differ?  We would need a test for this, but only the
      // Sacado and Stokhos packages are likely to care about this use
      // case.
      Kokkos::deep_copy (norms, theNorms);
    }

    void norm2 (const Teuchos::ArrayView<mag_type>& norms) const;

    template <typename T>
    typename Kokkos::Impl::enable_if< !(Kokkos::Impl::is_same<mag_type,T>::value), void >::type
    norm2 (const Teuchos::ArrayView<T>& norms) const
    {
      typedef typename Teuchos::ArrayView<T>::size_type size_type;
      const size_type sz = norms.size ();
      Teuchos::Array<mag_type> theNorms (sz);
      this->norm2 (theNorms);
      for (size_type i = 0; i < sz; ++i) {
        norms[i] = theNorms[i];
      }
    }

    void normInf (const Kokkos::View<mag_type*, device_type>& norms) const;

    template<typename T>
    typename Kokkos::Impl::enable_if< !(Kokkos::Impl::is_same<mag_type, T>::value), void >::type
    normInf (const Kokkos::View<T*, device_type>& norms) const
    {
      const size_t numNorms = norms.dimension_0 ();
      Kokkos::View<mag_type*, device_type> theNorms ("MV::normInf tmp", numNorms);
      // Call overload that takes a Kokkos::View<mag_type*, device_type>.
      this->normInf (theNorms);
      // FIXME (mfh 15 Jul 2014) Does this actually work if mag_type
      // and T differ?  We would need a test for this, but only the
      // Sacado and Stokhos packages are likely to care about this use
      // case.
      Kokkos::deep_copy (norms, theNorms);
    }

    void normInf (const Teuchos::ArrayView<mag_type>& norms) const;

    template <typename T>
    typename Kokkos::Impl::enable_if< !(Kokkos::Impl::is_same<mag_type,T>::value), void >::type
    normInf (const Teuchos::ArrayView<T>& norms) const
    {
      typedef typename Teuchos::ArrayView<T>::size_type size_type;
      const size_type sz = norms.size ();
      Teuchos::Array<mag_type> theNorms (sz);
      this->norm2 (theNorms);
      for (size_type i = 0; i < sz; ++i) {
        norms[i] = theNorms[i];
      }
    }

    void
    normWeighted (const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& weights,
                  const Teuchos::ArrayView<mag_type>& norms) const;

    template <typename T>
    typename Kokkos::Impl::enable_if< !(Kokkos::Impl::is_same<mag_type,T>::value), void >::type
    normWeighted (const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& weights,
                  const Teuchos::ArrayView<T>& norms) const
    {
      typedef typename Teuchos::ArrayView<T>::size_type size_type;
      const size_type sz = norms.size ();
      Teuchos::Array<mag_type> theNorms (sz);
      this->normWeighted (weights, theNorms);
      for (size_type i = 0; i < sz; ++i) {
        norms[i] = theNorms[i];
      }
    }

    void meanValue(const Teuchos::ArrayView<Scalar> &means) const;

    void
    multiply (Teuchos::ETransp transA,
              Teuchos::ETransp transB,
              const Scalar& alpha,
              const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A,
              const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& B,
              const Scalar& beta);

    void
    elementWiseMultiply (Scalar scalarAB,
                         const Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A,
                         const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& B,
                         Scalar scalarThis);



    size_t getNumVectors() const;

    size_t getLocalLength() const;

    global_size_t getGlobalLength() const;

    size_t getStride() const;

    bool isConstantStride() const;




    virtual std::string description() const;

    virtual void
    describe (Teuchos::FancyOStream& out,
              const Teuchos::EVerbosityLevel verbLevel =
              Teuchos::Describable::verbLevel_default) const;

    virtual void
    removeEmptyProcessesInPlace (const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> >& newMap);

    void setCopyOrView (const Teuchos::DataAccess copyOrView) {
      TEUCHOS_TEST_FOR_EXCEPTION(
        copyOrView == Teuchos::Copy, std::invalid_argument,
        "Tpetra::MultiVector::setCopyOrView: The Kokkos refactor version of "
        "MultiVector _only_ implements view semantics.  You may not call this "
        "method with copyOrView = Teuchos::Copy.  The only valid argument is "
        "Teuchos::View.");
    }

    Teuchos::DataAccess getCopyOrView () const {
      return Teuchos::View;
    }

    void
    assign (const MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, node_type>& src);

  protected:
    template <class S, class LO, class GO, class D>
    friend MultiVector<S,LO,GO,Kokkos::Compat::KokkosDeviceWrapperNode<D> >
    createCopy (const MultiVector<S,LO,GO,Kokkos::Compat::KokkosDeviceWrapperNode<D> >& src);

    template <class DS, class DL, class DG, class DD, class SS, class SL, class SG, class SD>
    friend void
    deep_copy (MultiVector<DS,DL,DG,Kokkos::Compat::KokkosDeviceWrapperNode<DD> >& dst,
               const MultiVector<SS,SL,SG,Kokkos::Compat::KokkosDeviceWrapperNode<SD> >& src);

    typedef KokkosClassic::MultiVector<Scalar,Node> KMV;
    typedef KokkosClassic::DefaultArithmetic<KMV>   MVT;

    //KMV lclMV_;

    mutable dual_view_type view_;

    mutable dual_view_type origView_;

    Array<size_t> whichVectors_;



    template <class S,class LO,class GO,class N>
    friend Teuchos::RCP<MultiVector<S,LO,GO,N> >
    createMultiVectorFromView (const Teuchos::RCP<const Map<LO,GO,N> >&, const Teuchos::ArrayRCP<S>&, const size_t, const size_t);

    bool vectorIndexOutOfRange (size_t VectorIndex) const;

    template <class T>
    ArrayRCP<T> getSubArrayRCP(ArrayRCP<T> arr, size_t j) const;

    size_t getOrigNumLocalRows () const;

    size_t getOrigNumLocalCols () const;




    virtual bool
    checkSizes (const SrcDistObject& sourceObj);

    virtual size_t constantNumberOfPackets () const;

    virtual bool useNewInterface () { return true; }

    virtual void
    copyAndPermuteNew (
      const SrcDistObject& sourceObj,
      size_t numSameIDs,
      const Kokkos::View<const LocalOrdinal*, device_type>& permuteToLIDs,
      const Kokkos::View<const LocalOrdinal*, device_type>& permuteFromLIDs);

    virtual void
    packAndPrepareNew (
      const SrcDistObject& sourceObj,
      const Kokkos::View<const LocalOrdinal*, device_type>& exportLIDs,
      Kokkos::View<Scalar*, device_type>& exports,
      const Kokkos::View<size_t*, device_type>& numPacketsPerLID,
      size_t& constantNumPackets,
      Distributor &distor);

    virtual void
    unpackAndCombineNew (
      const Kokkos::View<const LocalOrdinal*, device_type>& importLIDs,
      const Kokkos::View<const Scalar*, device_type>& imports,
      const Kokkos::View<size_t*, device_type>& numPacketsPerLID,
      size_t constantNumPackets,
      Distributor &distor,
      CombineMode CM);

    void createViews () const;
    void createViewsNonConst (KokkosClassic::ReadWriteOption rwo);
    void releaseViews () const;

    typename dual_view_type::t_dev getKokkosView() const { return view_.d_view; }

  }; // class MultiVector

  namespace Details {

    // Partial specialization of MultiVectorCloner, for when the
    // source and destination MultiVector types are both Kokkos
    // refactor types, and when they both have the same Scalar type,
    // but all their other template parameters might be different.
    template<class ScalarType,
             class DstLocalOrdinalType, class DstGlobalOrdinalType, class DstDeviceType,
             class SrcLocalOrdinalType, class SrcGlobalOrdinalType, class SrcDeviceType>
    struct MultiVectorCloner< ::Tpetra::MultiVector<ScalarType,
                                                    DstLocalOrdinalType,
                                                    DstGlobalOrdinalType,
                                                    Kokkos::Compat::KokkosDeviceWrapperNode<DstDeviceType> >,
                              ::Tpetra::MultiVector<ScalarType,
                                                    SrcLocalOrdinalType,
                                                    SrcGlobalOrdinalType,
                                                    Kokkos::Compat::KokkosDeviceWrapperNode<SrcDeviceType> > >
    {
      typedef Kokkos::Compat::KokkosDeviceWrapperNode<DstDeviceType> dst_node_type;
      typedef Kokkos::Compat::KokkosDeviceWrapperNode<SrcDeviceType> src_node_type;
      typedef ::Tpetra::MultiVector<ScalarType, DstLocalOrdinalType,
                                    DstGlobalOrdinalType,
                                    dst_node_type> dst_mv_type;
      typedef ::Tpetra::MultiVector<ScalarType, SrcLocalOrdinalType,
                                    SrcGlobalOrdinalType,
                                    src_node_type> src_mv_type;

      static Teuchos::RCP<dst_mv_type>
      clone (const src_mv_type& X, const Teuchos::RCP<dst_node_type>& node2)
      {
        typedef typename src_mv_type::map_type src_map_type;
        typedef typename dst_mv_type::map_type dst_map_type;
        typedef typename dst_mv_type::node_type dst_node_type;
        typedef typename src_mv_type::dual_view_type src_dual_view_type;
        typedef typename dst_mv_type::dual_view_type dst_dual_view_type;

        // Clone X's Map to have the new Node type.
        RCP<const src_map_type> map1 = X.getMap ();
        RCP<const dst_map_type> map2 = map1.is_null () ?
          Teuchos::null : map1->template clone<dst_node_type> (node2);

        const size_t lclNumRows = X.getLocalLength ();
        const size_t numCols = X.getNumVectors ();
        src_dual_view_type X_view = X.getDualView ();
        dst_dual_view_type Y_view ("MV::dual_view", lclNumRows, numCols);

        RCP<dst_mv_type> Y = rcp (new dst_mv_type (map2, Y_view));
        // Let deep_copy do the work for us, to avoid code duplication.
        ::Tpetra::deep_copy (Y, X);
        return Y ;
      }
    };

    // Partial specialization for the Kokkos refactor specialization
    // of Tpetra::MultiVector.
    template<class S, class LO, class GO, class DeviceType>
    struct CreateMultiVectorFromView<MultiVector<S, LO, GO, Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType> > > {
      typedef S scalar_type;
      typedef LO local_ordinal_type;
      typedef GO global_ordinal_type;
      typedef Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType> node_type;
      typedef ::Tpetra::Map<local_ordinal_type, global_ordinal_type, node_type> map_type;
      typedef MultiVector<S, LO, GO, Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType> > MultiVectorType;

      static Teuchos::RCP<MultiVectorType>
      create (const Teuchos::RCP<const map_type>& map,
              const Teuchos::ArrayRCP<scalar_type>& view,
              const size_t LDA,
              const size_t numVectors)
      {
        (void) map;
        (void) view;
        (void) LDA;
        (void) numVectors;

        TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
      }
    };
  } // namespace Details

  template <class ST, class LO, class GO, class DT>
  MultiVector<ST, LO, GO, Kokkos::Compat::KokkosDeviceWrapperNode<DT> >
  createCopy (const MultiVector<ST, LO, GO, Kokkos::Compat::KokkosDeviceWrapperNode<DT> >& src);

  namespace { // (anonymous)
    template<class DstType, class SrcType, class IndexType, class DeviceType,
             const bool DstConstStride, const bool SrcConstStride>
    struct DeepCopySelectedVectors {
      typedef DeviceType device_type;
      DstType dst_;
      SrcType src_;
      Kokkos::View<const IndexType*, DeviceType> whichVectorDst_;
      Kokkos::View<const IndexType*, DeviceType> whichVectorSrc_;
      const IndexType numVecs_;

      DeepCopySelectedVectors (DstType dst,
                               SrcType src,
                               const Kokkos::View<const IndexType*, DeviceType>& whichVectorDst,
                               const Kokkos::View<const IndexType*, DeviceType>& whichVectorSrc) :
        dst_ (dst),
        src_ (src),
        whichVectorDst_ (whichVectorDst),
        whichVectorSrc_ (whichVectorSrc),
        numVecs_ (whichVectorSrc_.dimension_0 ())
      {}

      DeepCopySelectedVectors (DstType dst, SrcType src) :
        dst_ (dst),
        src_ (src),
        numVecs_ (dst.dimension_1 ())
      {
        TEUCHOS_TEST_FOR_EXCEPTION(
                                   ! DstConstStride || ! SrcConstStride, std::logic_error,
                                   "Tpetra::DeepCopySelectedVectors: You may not use the constant-stride "
                                   "constructor if either of the Boolean template parameters is false.");
      }

      void KOKKOS_INLINE_FUNCTION operator () (const IndexType i) const {
        if (DstConstStride) {
          if (SrcConstStride) {
            for (IndexType j = 0; j < numVecs_; ++j) {
              dst_(i,j) = src_(i,j);
            }
          } else {
            for (IndexType j = 0; j < numVecs_; ++j) {
              dst_(i,j) = src_(i,whichVectorSrc_(j));
            }
          }
        } else {
          if (SrcConstStride) {
            for (IndexType j = 0; j < numVecs_; ++j) {
              dst_(i,whichVectorDst_(j)) = src_(i,j);
            }
          } else {
            for (IndexType j = 0; j < numVecs_; ++j) {
              dst_(i,whichVectorDst_(j)) = src_(i,whichVectorSrc_(j));
            }
          }
        }
      }
    };
  } // namespace (anonymous)


  // NOTE (mfh 11 Sep 2014) Even though this partial specialization
  // looks redundant with the one in Tpetra_MultiVector_decl.hpp, it
  // needs to be here, else GCC 4.8.2 gives a compiler error saying
  // that calls to Tpetra::deep_copy are ambiguous.
  template <class ST, class LO, class GO, class DT>
  void
  deep_copy (MultiVector<ST, LO, GO, Kokkos::Compat::KokkosDeviceWrapperNode<DT> >& dst,
             const MultiVector<ST, LO, GO, Kokkos::Compat::KokkosDeviceWrapperNode<DT> >& src)
  {
    // NOTE (mfh 11 Sep 2014) We can't implement deep_copy with
    // shallow-copy operator=, because that would invalidate existing
    // views of dst!
    dst.assign (src);
  }


  template <class DS, class DL, class DG, class DD,
            class SS, class SL, class SG, class SD>
  void
  deep_copy (MultiVector<DS, DL, DG, Kokkos::Compat::KokkosDeviceWrapperNode<DD> >& dst,
             const MultiVector<SS, SL, SG, Kokkos::Compat::KokkosDeviceWrapperNode<SD> >& src)
  {
    using Kokkos::parallel_for;
    typedef MultiVector<DS,DL,DG,Kokkos::Compat::KokkosDeviceWrapperNode<DD> > MVD;
    typedef const MultiVector<SS,SL,SG,Kokkos::Compat::KokkosDeviceWrapperNode<SD> > MVS;

    TEUCHOS_TEST_FOR_EXCEPTION(
      dst.getGlobalLength () != src.getGlobalLength () ||
      dst.getNumVectors () != src.getNumVectors (), std::invalid_argument,
      "Tpetra::deep_copy: Global dimensions of the two Tpetra::MultiVector "
      "objects do not match.  src has dimensions [" << src.getGlobalLength ()
      << "," << src.getNumVectors () << "], and dst has dimensions ["
      << dst.getGlobalLength () << "," << dst.getNumVectors () << "].");

    // FIXME (mfh 28 Jul 2014) Don't throw; just set a local error flag.
    TEUCHOS_TEST_FOR_EXCEPTION(
      dst.getLocalLength () != src.getLocalLength (), std::invalid_argument,
      "Tpetra::deep_copy: The local row counts of the two Tpetra::MultiVector "
      "objects do not match.  src has " << src.getLocalLength () << " row(s) "
      << " and dst has " << dst.getLocalLength () << " row(s).");

    if (src.isConstantStride () && dst.isConstantStride ()) {
      Kokkos::deep_copy (dst.getDualView (), src.getDualView ());
    }
    else {
      if (dst.isConstantStride ()) {
        const SL numWhichVecs = static_cast<SL> (src.whichVectors_.size ());
        const std::string whichVecsLabel ("MV::deep_copy::whichVecs");

        // We can't sync src, since it is only an input argument.
        // Thus, we have to use the most recently modified version of
        // src, device or host.
        if (src.getDualView ().modified_device >= src.getDualView ().modified_host) {
          // Copy from the device version of src.
          //
          // whichVecs tells the kernel which vectors (columns) of src
          // to copy.  Fill whichVecs on the host, and sync to device.
          typedef Kokkos::DualView<SL*, DD> whichvecs_type;
          whichvecs_type whichVecs (whichVecsLabel, numWhichVecs);
          whichVecs.template modify<typename DD::host_mirror_device_type> ();
          for (SL i = 0; i < numWhichVecs; ++i) {
            whichVecs.h_view(i) = static_cast<SL> (src.whichVectors_[i]);
          }
          // Sync the host version of whichVecs to the device.
          whichVecs.template sync<DD> ();

          // Mark the device version of dst's DualView as modified.
          dst.template modify<DD> ();
          // Copy from the selected vectors of src to dst, on the
          // device.  The functor ignores its 3rd arg in this case.
          typedef DeepCopySelectedVectors<typename MVD::dual_view_type::t_dev,
            typename MVS::dual_view_type::t_dev, SL, DD, true, false> functor_type;
          functor_type f (dst.getDualView ().template view<DD> (),
                          src.getDualView ().template view<DD> (),
                          whichVecs.d_view, whichVecs.d_view);
          Kokkos::parallel_for (src.getLocalLength (), f);
          // Sync dst's DualView to the host.  This is cheaper than
          // repeating the above copy from src to dst on the host.
          dst.template sync<typename DD::host_mirror_device_type> ();
        }
        else { // host version of src was the most recently modified
          // Copy from the host version of src.
          //
          // whichVecs tells the kernel which vectors (columns) of src
          // to copy.  Fill whichVecs on the host, and use it there.
          typedef typename DD::host_mirror_device_type host_dev_type;
          typedef Kokkos::View<SL*, host_dev_type> whichvecs_type;
          whichvecs_type whichVecs (whichVecsLabel, numWhichVecs);
          for (SL i = 0; i < numWhichVecs; ++i) {
            whichVecs(i) = static_cast<SL> (src.whichVectors_[i]);
          }
          // Copy from the selected vectors of src to dst, on the host.
          // The functor ignores its 3rd arg in this case.
          typedef DeepCopySelectedVectors<typename MVD::dual_view_type::t_host,
            typename MVS::dual_view_type::t_host, SL, host_dev_type,
            true, false> functor_type;
          functor_type f (dst.getDualView ().template view<host_dev_type> (),
                          src.getDualView ().template view<host_dev_type> (),
                          whichVecs, whichVecs);
          Kokkos::parallel_for (src.getLocalLength (), f);
          // Sync dst back to the device, since we only copied on the host.
          dst.template sync<DD> ();
        }
      }
      else { // dst is NOT constant stride
        if (src.isConstantStride ()) {
          if (src.getDualView ().modified_device >= src.getDualView ().modified_host) {
            // Copy from the device version of src.
            //
            // whichVecs tells the kernel which vectors (columns) of dst
            // to copy.  Fill whichVecs on the host, and sync to device.
            typedef Kokkos::DualView<DL*, DD> whichvecs_type;
            const std::string whichVecsLabel ("MV::deep_copy::whichVecs");
            const DL numWhichVecs = static_cast<DL> (dst.whichVectors_.size ());
            whichvecs_type whichVecs (whichVecsLabel, numWhichVecs);
            whichVecs.template modify<typename DD::host_mirror_device_type> ();
            for (DL i = 0; i < numWhichVecs; ++i) {
              whichVecs.h_view(i) = dst.whichVectors_[i];
            }
            // Sync the host version of whichVecs to the device.
            whichVecs.template sync<DD> ();

            // Copy src to the selected vectors of dst, on the device.
            // The functor ignores its 4th arg in this case.
            typedef DeepCopySelectedVectors<typename MVD::dual_view_type::t_dev,
              typename MVS::dual_view_type::t_dev, DL, DD, false, true> functor_type;
            functor_type f (dst.getDualView ().template view<DD> (),
                            src.getDualView ().template view<DD> (),
                            whichVecs.d_view, whichVecs.d_view);
            Kokkos::parallel_for (src.getLocalLength (), f);
            // We can't sync src and repeat the above copy on the
            // host, so sync dst back to the host.
            //
            // FIXME (mfh 29 Jul 2014) This may overwrite columns that
            // don't actually belong to dst's view.
            dst.template sync<typename DD::host_mirror_device_type> ();
          }
          else { // host version of src was the most recently modified
            // Copy from the host version of src.
            //
            // whichVecs tells the kernel which vectors (columns) of src
            // to copy.  Fill whichVecs on the host, and use it there.
            typedef typename DD::host_mirror_device_type host_dev_type;
            typedef Kokkos::View<DL*, host_dev_type> whichvecs_type;
            const DL numWhichVecs = static_cast<DL> (dst.whichVectors_.size ());
            whichvecs_type whichVecs ("MV::deep_copy::whichVecs", numWhichVecs);
            for (DL i = 0; i < numWhichVecs; ++i) {
              whichVecs(i) = static_cast<DL> (dst.whichVectors_[i]);
            }
            // Copy from src to the selected vectors of dst, on the
            // host.  The functor ignores its 4th arg in this case.
            typedef DeepCopySelectedVectors<typename MVD::dual_view_type::t_host,
              typename MVS::dual_view_type::t_host, DL, host_dev_type,
              false, true> functor_type;
            functor_type f (dst.getDualView ().template view<host_dev_type> (),
                            src.getDualView ().template view<host_dev_type> (),
                            whichVecs, whichVecs);
            Kokkos::parallel_for (src.getLocalLength (), f);
            // Sync dst back to the device, since we only copied on the host.
            //
            // FIXME (mfh 29 Jul 2014) This may overwrite columns that
            // don't actually belong to dst's view.
            dst.template sync<DD> ();
          }
        }
        else { // neither src nor dst have constant stride
          if (src.getDualView ().modified_device >= src.getDualView ().modified_host) {
            // Copy from the device version of src.
            //
            // whichVectorsDst tells the kernel which vectors
            // (columns) of dst to copy.  Fill it on the host, and
            // sync to device.
            const DL dstNumWhichVecs = static_cast<DL> (dst.whichVectors_.size ());
            Kokkos::DualView<DL*, DD> whichVecsDst ("MV::deep_copy::whichVecsDst",
                                                    dstNumWhichVecs);
            whichVecsDst.template modify<typename DD::host_mirror_device_type> ();
            for (DL i = 0; i < dstNumWhichVecs; ++i) {
              whichVecsDst.h_view(i) = static_cast<DL> (dst.whichVectors_[i]);
            }
            // Sync the host version of whichVecsDst to the device.
            whichVecsDst.template sync<DD> ();

            // whichVectorsSrc tells the kernel which vectors
            // (columns) of src to copy.  Fill it on the host, and
            // sync to device.  Use the destination MultiVector's
            // LocalOrdinal type here.
            const DL srcNumWhichVecs = static_cast<DL> (src.whichVectors_.size ());
            Kokkos::DualView<DL*, DD> whichVecsSrc ("MV::deep_copy::whichVecsSrc",
                                                    srcNumWhichVecs);
            whichVecsSrc.template modify<typename DD::host_mirror_device_type> ();
            for (DL i = 0; i < srcNumWhichVecs; ++i) {
              whichVecsSrc.h_view(i) = static_cast<DL> (src.whichVectors_[i]);
            }
            // Sync the host version of whichVecsSrc to the device.
            whichVecsSrc.template sync<DD> ();

            // Copy from the selected vectors of src to the selected
            // vectors of dst, on the device.
            typedef DeepCopySelectedVectors<typename MVD::dual_view_type::t_dev,
              typename MVS::dual_view_type::t_dev, DL, DD, false, false>
              functor_type;
            functor_type f (dst.getDualView ().template view<DD> (),
                            src.getDualView ().template view<DD> (),
                            whichVecsDst.d_view, whichVecsSrc.d_view);
            Kokkos::parallel_for (src.getLocalLength (), f);
          }
          else {
            typedef typename DD::host_mirror_device_type host_dev_type;

            const DL dstNumWhichVecs = static_cast<DL> (dst.whichVectors_.size ());
            Kokkos::View<DL*, host_dev_type> whichVectorsDst ("dstWhichVecs", dstNumWhichVecs);
            for (DL i = 0; i < dstNumWhichVecs; ++i) {
              whichVectorsDst(i) = dst.whichVectors_[i];
            }

            // Use the destination MultiVector's LocalOrdinal type here.
            const DL srcNumWhichVecs = static_cast<DL> (src.whichVectors_.size ());
            Kokkos::View<DL*, host_dev_type> whichVectorsSrc ("srcWhichVecs", srcNumWhichVecs);
            for (DL i = 0; i < srcNumWhichVecs; ++i) {
              whichVectorsSrc(i) = src.whichVectors_[i];
            }

            typedef DeepCopySelectedVectors<typename MVD::dual_view_type::t_host,
              typename MVS::dual_view_type::t_host,
              DL, host_dev_type, false, false> functor_type;
            functor_type f (dst.getDualView ().template view<host_dev_type> (),
                            src.getDualView ().template view<host_dev_type> (),
                            whichVectorsDst, whichVectorsSrc);
            Kokkos::parallel_for (src.getLocalLength (), f);

            // We can't sync src and repeat the above copy on the
            // host, so sync dst back to the host.
            //
            // FIXME (mfh 29 Jul 2014) This may overwrite columns that
            // don't actually belong to dst's view.
            dst.template sync<typename DD::host_mirror_device_type> ();
          }
        }
      }
    }
  }
} // namespace Tpetra

#endif // TPETRA_KOKKOS_REFACTOR_MULTIVECTOR_DECL_HPP