Tpetra parallel linear algebra  Version of the Day
Tpetra_KokkosRefactor_Details_MultiVectorLocalDeepCopy.hpp
1 /*
2 // @HEADER
3 // ***********************************************************************
4 //
5 // Tpetra: Templated Linear Algebra Services Package
6 // Copyright (2008) Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
39 //
40 // ************************************************************************
41 // @HEADER
42 */
43 
44 #ifndef TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_LOCAL_DEEP_COPY_HPP
45 #define TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_LOCAL_DEEP_COPY_HPP
46 
47 #include "Kokkos_Core.hpp"
48 #include "Kokkos_ArithTraits.hpp"
49 
50 namespace Tpetra {
51 namespace Details {
52 
53  // Functor for deep-copying between two 2-D Kokkos Views.
54  // This implements Tpetra::MultiVector deep copy, as in
55  // - Tpetra::deep_copy
56  // - Tpetra::MultiVector::assign
57  // - Tpetra::MultiVector::createCopy
58  // - The two-argument MultiVector copy constructor with
59  // Teuchos::Copy as the second argument
60  //
61  // DstViewType and SrcViewType must be 2-D Kokkos Views.
62  // DstWhichVecsType and SrcWhichVecsType must be 1-D Kokkos Views
63  // whose value type is an integer. They correspond to the "which
64  // vectors?" arrays for the destination and source Views,
65  // respectively.
66  //
67  // If DstConstStride is true, dstWhichVecs is ignored. If
68  // SrcConstStride is true, srcWhichVecs is ignored. These bool
69  // template parameters take advantage of compilers' ability to
70  // disable code inside an "if (false) { ... }" scope.
71  template<class DstViewType,
72  class SrcViewType,
73  class DstWhichVecsType,
74  class SrcWhichVecsType,
75  const bool DstConstStride,
76  const bool SrcConstStride>
77  struct LocalDeepCopyFunctor {
78  typedef typename DstViewType::execution_space execution_space;
79  typedef typename DstViewType::size_type index_type;
80 
81  DstViewType dst_;
82  SrcViewType src_;
83  DstWhichVecsType dstWhichVecs_;
84  SrcWhichVecsType srcWhichVecs_;
85  const index_type numVecs_;
86 
87  // You may always call the 4-argument constructor.
88  // dstWhichVecs is ignored if DstConstStride is true.
89  // srcWhichVecs is ignored if SrcConstStride is true.
90  LocalDeepCopyFunctor (const DstViewType& dst,
91  const SrcViewType& src,
92  const DstWhichVecsType& dstWhichVecs,
93  const SrcWhichVecsType& srcWhichVecs) :
94  dst_ (dst),
95  src_ (src),
96  dstWhichVecs_ (dstWhichVecs),
97  srcWhichVecs_ (srcWhichVecs),
98  numVecs_ (DstConstStride ? dst.dimension_1 () : dstWhichVecs.dimension_0 ())
99  {
100  TEUCHOS_TEST_FOR_EXCEPTION(
101  ! DstConstStride && ! SrcConstStride &&
102  dstWhichVecs.dimension_0 () != srcWhichVecs.dimension_0 (),
103  std::invalid_argument, "LocalDeepCopyFunctor (4-arg constructor): "
104  "Neither src nor dst have constant stride, but "
105  "dstWhichVecs.dimension_0() = " << dstWhichVecs.dimension_0 ()
106  << " != srcWhichVecs.dimension_0() = " << srcWhichVecs.dimension_0 ()
107  << ".");
108  TEUCHOS_TEST_FOR_EXCEPTION(
109  DstConstStride && ! SrcConstStride &&
110  srcWhichVecs.dimension_0 () != dst.dimension_1 (),
111  std::invalid_argument, "LocalDeepCopyFunctor (4-arg constructor): "
112  "src does not have constant stride, but srcWhichVecs.dimension_0() = "
113  << srcWhichVecs.dimension_0 () << " != dst.dimension_1() = "
114  << dst.dimension_1 () << ".");
115  TEUCHOS_TEST_FOR_EXCEPTION(
116  ! DstConstStride && SrcConstStride &&
117  dstWhichVecs.dimension_0 () != src.dimension_1 (),
118  std::invalid_argument, "LocalDeepCopyFunctor (4-arg constructor): "
119  "dst does not have constant stride, but dstWhichVecs.dimension_0() = "
120  << dstWhichVecs.dimension_0 () << " != src.dimension_1() = "
121  << src.dimension_1 () << ".");
122  }
123 
124  // You may only call the 2-argument constructor if DstConstStride
125  // and SrcConstStride are both true.
126  LocalDeepCopyFunctor (const DstViewType& dst, const SrcViewType& src) :
127  dst_ (dst),
128  src_ (src),
129  numVecs_ (dst.dimension_1 ())
130  {
131  TEUCHOS_TEST_FOR_EXCEPTION(
132  ! DstConstStride || ! SrcConstStride, std::logic_error,
133  "Tpetra::LocalDeepCopyFunctor: You may not use the constant-stride "
134  "constructor if either of the Boolean template parameters is false.");
135  }
136 
137  void KOKKOS_INLINE_FUNCTION operator () (const index_type i) const {
138  if (DstConstStride) {
139  if (SrcConstStride) {
140  for (index_type j = 0; j < numVecs_; ++j) {
141  dst_(i,j) = src_(i,j);
142  }
143  } else {
144  for (index_type j = 0; j < numVecs_; ++j) {
145  dst_(i,j) = src_(i,srcWhichVecs_(j));
146  }
147  }
148  } else {
149  if (SrcConstStride) {
150  for (index_type j = 0; j < numVecs_; ++j) {
151  dst_(i,dstWhichVecs_(j)) = src_(i,j);
152  }
153  } else {
154  for (index_type j = 0; j < numVecs_; ++j) {
155  dst_(i,dstWhichVecs_(j)) = src_(i,srcWhichVecs_(j));
156  }
157  }
158  }
159  }
160  };
161 
162  template<class DstViewType,
163  class SrcViewType = DstViewType,
164  class DstWhichVecsType = Kokkos::View<const typename DstViewType::size_type*, typename DstViewType::execution_space>,
165  class SrcWhichVecsType = Kokkos::View<const typename SrcViewType::size_type*, typename SrcViewType::execution_space> >
166  struct LocalDeepCopy {
167  typedef typename DstViewType::size_type index_type;
168 
169  static void
170  run (const DstViewType& dst, const SrcViewType& src,
171  const bool dstConstStride, const bool srcConstStride)
172  {
173  TEUCHOS_TEST_FOR_EXCEPTION(
174  ! dstConstStride || ! srcConstStride, std::invalid_argument,
175  "LocalDeepCopy::run: You may only call the 4-argument version of this "
176  "function if dstConstStride and srcConstStride are both true.");
177 
178  // FIXME (mfh 22 Jul 2014, 10 Dec 2014) Currently, it doesn't
179  // work to do a 2-D copy, even if both MultiVectors have
180  // constant stride. This is because Kokkos can't currently tell
181  // the difference between padding (which permits a single
182  // deep_copy for the whole 2-D View) and stride > numRows (which
183  // does NOT permit a single deep_copy for the whole 2-D View).
184  // Carter is working on this, but for now, the temporary fix is
185  // to copy one column at a time.
186  //
187  // FIXME (mfh 10 Dec 2014) Call Kokkos::deep_copy if appropriate
188  // for dst and src. See note above.
189  typedef LocalDeepCopyFunctor<DstViewType, SrcViewType,
190  DstWhichVecsType, SrcWhichVecsType, true, true> functor_type;
191  functor_type f (dst, src);
192  Kokkos::parallel_for (dst.dimension_0 (), f);
193  }
194 
195  static void
196  run (const DstViewType& dst, const SrcViewType& src,
197  const bool dstConstStride, const bool srcConstStride,
198  const DstWhichVecsType& dstWhichVecs,
199  const SrcWhichVecsType& srcWhichVecs)
200  {
201  // FIXME (mfh 22 Jul 2014, 10 Dec 2014) Currently, it doesn't
202  // work to do a 2-D copy, even if both MultiVectors have
203  // constant stride. This is because Kokkos can't currently tell
204  // the difference between padding (which permits a single
205  // deep_copy for the whole 2-D View) and stride > numRows (which
206  // does NOT permit a single deep_copy for the whole 2-D View).
207  // Carter is working on this, but for now, the temporary fix is
208  // to copy one column at a time.
209 
210  if (dstConstStride) {
211  if (srcConstStride) {
212  // FIXME (mfh 10 Dec 2014) Do a Kokkos::deep_copy if
213  // appropriate for dst and src. See note above.
214  typedef LocalDeepCopyFunctor<DstViewType, SrcViewType,
215  DstWhichVecsType, SrcWhichVecsType, true, true> functor_type;
216  functor_type f (dst, src);
217  Kokkos::parallel_for (dst.dimension_0 (), f);
218  }
219  else { // ! srcConstStride
220  typedef LocalDeepCopyFunctor<DstViewType, SrcViewType,
221  DstWhichVecsType, SrcWhichVecsType, true, false> functor_type;
222  functor_type f (dst, src, srcWhichVecs, srcWhichVecs);
223  Kokkos::parallel_for (dst.dimension_0 (), f);
224  }
225  }
226  else { // ! dstConstStride
227  if (srcConstStride) {
228  typedef LocalDeepCopyFunctor<DstViewType, SrcViewType,
229  DstWhichVecsType, SrcWhichVecsType, false, true> functor_type;
230  functor_type f (dst, src, dstWhichVecs, dstWhichVecs);
231  Kokkos::parallel_for (dst.dimension_0 (), f);
232  }
233  else { // ! srcConstStride
234  typedef LocalDeepCopyFunctor<DstViewType, SrcViewType,
235  DstWhichVecsType, SrcWhichVecsType, false, false> functor_type;
236  functor_type f (dst, src, dstWhichVecs, srcWhichVecs);
237  Kokkos::parallel_for (dst.dimension_0 (), f);
238  }
239  }
240  }
241  };
242 
243 
244  template<class DstViewType,
245  class SrcViewType,
246  class DstWhichVecsType,
247  class SrcWhichVecsType>
248  void
249  localDeepCopy (const DstViewType& dst, const SrcViewType& src,
250  const bool dstConstStride, const bool srcConstStride,
251  const DstWhichVecsType& dstWhichVecs,
252  const SrcWhichVecsType& srcWhichVecs)
253  {
254  typedef LocalDeepCopy<DstViewType, SrcViewType,
255  DstWhichVecsType, SrcWhichVecsType> impl_type;
256  impl_type::run (dst, src, dstConstStride, srcConstStride,
257  dstWhichVecs, srcWhichVecs);
258  }
259 
260  template<class DstViewType, class SrcViewType>
261  void
262  localDeepCopyConstStride (const DstViewType& dst, const SrcViewType& src)
263  {
264  typedef LocalDeepCopy<DstViewType, SrcViewType> impl_type;
265  impl_type::run (dst, src, true, true);
266  }
267 
268 } // Details namespace
269 } // Tpetra namespace
270 
271 #endif // TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_LOCAL_DEEP_COPY_HPP
Namespace Tpetra contains the class and methods constituting the Tpetra library.
Implementation details of Tpetra.