Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
Tpetra_KokkosRefactor_Details_MultiVectorDistObjectKernels_MP_Vector.hpp
Go to the documentation of this file.
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_DIST_OBJECT_KERNELS_MP_VECTOR_HPP
45 #define TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_DIST_OBJECT_KERNELS_MP_VECTOR_HPP
46 
47 //----------------------------------------------------------------------------
48 // Specializations of Tpetra::MultiVector pack/unpack kernels for MP::Vector
49 //----------------------------------------------------------------------------
50 
51 #include "Tpetra_KokkosRefactor_Details_MultiVectorDistObjectKernels.hpp"
54 
55 namespace Tpetra {
56 namespace KokkosRefactor {
57 namespace Details {
58 
59 #if defined( KOKKOS_ENABLE_CUDA )
60  template< class D >
61  struct device_is_cuda : public std::is_same<D,Kokkos::Cuda> {};
62 #else
63  template< class D >
64  struct device_is_cuda : public std::false_type {};
65 #endif
66 
67  // Functors for implementing packAndPrepare and unpackAndCombine
68  // through parallel_for
69 
70  template <typename DS, typename ... DP,
71  typename SS, typename ... SP,
72  typename IdxView>
73  struct PackArraySingleColumn<
74  Kokkos::View<Sacado::MP::Vector<DS>*,DP...>,
75  Kokkos::View<const Sacado::MP::Vector<SS>**,SP...>,
76  IdxView >
77  {
78  typedef Kokkos::View<Sacado::MP::Vector<DS>*,DP...> DstView;
79  typedef Kokkos::View<const Sacado::MP::Vector<SS>**,SP...> SrcView;
81  typedef typename execution_space::size_type size_type;
82 
85  IdxView idx;
86  size_t col;
87 
89  const SrcView& src_,
90  const IdxView& idx_,
91  size_t col_) :
92  dst(dst_), src(src_), idx(idx_), col(col_) {}
93 
94  KOKKOS_INLINE_FUNCTION
95  void operator()( const size_type k ) const {
96  dst(k) = src(idx(k), col);
97  }
98 
99  KOKKOS_INLINE_FUNCTION
100  void operator()( const size_type k, const size_type tidx ) const {
101  dst(k).fastAccessCoeff(tidx) = src(idx(k), col).fastAccessCoeff(tidx);
102  }
103 
104  static void pack(const DstView& dst,
105  const SrcView& src,
106  const IdxView& idx,
107  size_t col) {
109  Kokkos::parallel_for(
111  PackArraySingleColumn(dst,src,idx,col) );
112  else
113  Kokkos::parallel_for( idx.size(),
114  PackArraySingleColumn(dst,src,idx,col) );
115  }
116  };
117 
118  template <typename DS, typename ... DP,
119  typename SS, typename ... SP,
120  typename IdxView>
121  struct PackArrayMultiColumn<
122  Kokkos::View<Sacado::MP::Vector<DS>*,DP...>,
123  Kokkos::View<const Sacado::MP::Vector<SS>**,SP...>,
124  IdxView >
125  {
126  typedef Kokkos::View<Sacado::MP::Vector<DS>*,DP...> DstView;
127  typedef Kokkos::View<const Sacado::MP::Vector<SS>**,SP...> SrcView;
129  typedef typename execution_space::size_type size_type;
130 
133  IdxView idx;
134  size_t numCols;
135 
137  const SrcView& src_,
138  const IdxView& idx_,
139  size_t numCols_) :
140  dst(dst_), src(src_), idx(idx_), numCols(numCols_) {}
141 
142  KOKKOS_INLINE_FUNCTION
143  void operator()( const size_type k ) const {
144  const typename IdxView::value_type localRow = idx(k);
145  const size_t offset = k*numCols;
146  for (size_t j = 0; j < numCols; ++j)
147  dst(offset + j) = src(localRow, j);
148  }
149 
150  KOKKOS_INLINE_FUNCTION
151  void operator()( const size_type k, const size_type tidx ) const {
152  const typename IdxView::value_type localRow = idx(k);
153  const size_t offset = k*numCols;
154  for (size_t j = 0; j < numCols; ++j)
155  dst(offset + j).fastAccessCoeff(tidx) =
156  src(localRow, j).fastAccessCoeff(tidx);
157  }
158 
159  static void pack(const DstView& dst,
160  const SrcView& src,
161  const IdxView& idx,
162  size_t numCols) {
164  Kokkos::parallel_for(
166  PackArrayMultiColumn(dst,src,idx,numCols) );
167  else
168  Kokkos::parallel_for( idx.size(),
169  PackArrayMultiColumn(dst,src,idx,numCols) );
170  }
171  };
172 
173  template <typename DS, typename ... DP,
174  typename SS, typename ... SP,
175  typename IdxView,
176  typename ColView>
177  struct PackArrayMultiColumnVariableStride<
178  Kokkos::View<Sacado::MP::Vector<DS>*,DP...>,
179  Kokkos::View<const Sacado::MP::Vector<SS>**,SP...>,
180  IdxView,
181  ColView>
182  {
183  typedef Kokkos::View<Sacado::MP::Vector<DS>*,DP...> DstView;
184  typedef Kokkos::View<const Sacado::MP::Vector<SS>**,SP...> SrcView;
186  typedef typename execution_space::size_type size_type;
187 
190  IdxView idx;
191  ColView col;
192  size_t numCols;
193 
195  const SrcView& src_,
196  const IdxView& idx_,
197  const ColView& col_,
198  size_t numCols_) :
199  dst(dst_), src(src_), idx(idx_), col(col_), numCols(numCols_) {}
200 
201  KOKKOS_INLINE_FUNCTION
202  void operator()( const size_type k ) const {
203  const typename IdxView::value_type localRow = idx(k);
204  const size_t offset = k*numCols;
205  for (size_t j = 0; j < numCols; ++j)
206  dst(offset + j) = src(localRow, col(j));
207  }
208 
209  KOKKOS_INLINE_FUNCTION
210  void operator()( const size_type k, const size_type tidx ) const {
211  const typename IdxView::value_type localRow = idx(k);
212  const size_t offset = k*numCols;
213  for (size_t j = 0; j < numCols; ++j)
214  dst(offset + j).fastAccessCoeff(tidx) =
215  src(localRow, col(j)).fastAccessCoeff(tidx);
216  }
217 
218  static void pack(const DstView& dst,
219  const SrcView& src,
220  const IdxView& idx,
221  const ColView& col,
222  size_t numCols) {
224  Kokkos::parallel_for(
226  PackArrayMultiColumnVariableStride(dst,src,idx,col,numCols) );
227  else
228  Kokkos::parallel_for( idx.size(),
229  PackArrayMultiColumnVariableStride(
230  dst,src,idx,col,numCols) );
231  }
232  };
233 
234  template <typename ExecutionSpace,
235  typename DS,
236  typename ... DP,
237  typename SS,
238  typename ... SP,
239  typename IdxView,
240  typename Op>
241  struct UnpackArrayMultiColumn<
242  ExecutionSpace,
243  Kokkos::View<Sacado::MP::Vector<DS>**,DP...>,
244  Kokkos::View<const Sacado::MP::Vector<SS>*,SP...>,
245  IdxView,
246  Op >
247  {
248  typedef Kokkos::View<Sacado::MP::Vector<DS>**,DP...> DstView;
249  typedef Kokkos::View<const Sacado::MP::Vector<SS>*,SP...> SrcView;
251  typedef typename execution_space::size_type size_type;
252 
255  IdxView idx;
256  Op op;
257  size_t numCols;
258 
259  UnpackArrayMultiColumn (const ExecutionSpace& /* execSpace */,
260  const DstView& dst_,
261  const SrcView& src_,
262  const IdxView& idx_,
263  const Op& op_,
264  const size_t numCols_) :
265  dst (dst_),
266  src (src_),
267  idx (idx_),
268  op (op_),
269  numCols (numCols_)
270  {}
271 
272  KOKKOS_INLINE_FUNCTION void
273  operator() (const size_type k) const
274  {
275  const typename IdxView::value_type localRow = idx(k);
276  const size_t offset = k*numCols;
277  for (size_t j = 0; j < numCols; ++j) {
278  op (dst(localRow, j), src(offset+j));
279  }
280  }
281 
282  KOKKOS_INLINE_FUNCTION void
283  operator() (const size_type k, const size_type tidx) const
284  {
285  const typename IdxView::value_type localRow = idx(k);
286  const size_t offset = k*numCols;
287  for (size_t j = 0; j < numCols; ++j) {
288  op (dst(localRow, j).fastAccessCoeff(tidx),
289  src(offset+j).fastAccessCoeff(tidx));
290  }
291  }
292 
293  static void
294  unpack (const ExecutionSpace& execSpace,
295  const DstView& dst,
296  const SrcView& src,
297  const IdxView& idx,
298  const Op& op,
299  const size_t numCols)
300  {
302  Kokkos::parallel_for
303  ("Tpetra::MultiVector (Sacado::MP::Vector) unpack (constant stride)",
305  UnpackArrayMultiColumn (execSpace, dst, src, idx, op, numCols));
306  }
307  else {
308  Kokkos::parallel_for
309  ("Tpetra::MultiVector (Sacado::MP::Vector) unpack (constant stride)",
310  Kokkos::RangePolicy<execution_space, size_type> (0, idx.size ()),
311  UnpackArrayMultiColumn (execSpace, dst, src, idx, op, numCols));
312  }
313  }
314  };
315 
316  template <typename ExecutionSpace,
317  typename DS,
318  typename ... DP,
319  typename SS,
320  typename ... SP,
321  typename IdxView,
322  typename ColView,
323  typename Op>
324  struct UnpackArrayMultiColumnVariableStride<
325  ExecutionSpace,
326  Kokkos::View<Sacado::MP::Vector<DS>**,DP...>,
327  Kokkos::View<const Sacado::MP::Vector<SS>*,SP...>,
328  IdxView,
329  ColView,
330  Op>
331  {
332  typedef Kokkos::View<Sacado::MP::Vector<DS>**,DP...> DstView;
333  typedef Kokkos::View<const Sacado::MP::Vector<SS>*,SP...> SrcView;
335  typedef typename execution_space::size_type size_type;
336 
339  IdxView idx;
340  ColView col;
341  Op op;
342  size_t numCols;
343 
344  UnpackArrayMultiColumnVariableStride (const ExecutionSpace& /* execSpace */,
345  const DstView& dst_,
346  const SrcView& src_,
347  const IdxView& idx_,
348  const ColView& col_,
349  const Op& op_,
350  const size_t numCols_) :
351  dst (dst_),
352  src (src_),
353  idx (idx_),
354  col (col_),
355  op (op_),
356  numCols (numCols_)
357  {}
358 
359  KOKKOS_INLINE_FUNCTION void
360  operator() (const size_type k) const
361  {
362  const typename IdxView::value_type localRow = idx(k);
363  const size_t offset = k*numCols;
364  for (size_t j = 0; j < numCols; ++j) {
365  op (dst(localRow, col(j)), src(offset+j));
366  }
367  }
368 
369  KOKKOS_INLINE_FUNCTION void
370  operator() (const size_type k, const size_type tidx) const
371  {
372  const typename IdxView::value_type localRow = idx(k);
373  const size_t offset = k*numCols;
374  for (size_t j = 0; j < numCols; ++j) {
375  op (dst(localRow, col(j)).fastAccessCoeff(tidx),
376  src(offset+j).fastAccessCoeff(tidx));
377  }
378  }
379 
380  static void
381  unpack (const ExecutionSpace& execSpace,
382  const DstView& dst,
383  const SrcView& src,
384  const IdxView& idx,
385  const ColView& col,
386  const Op& op,
387  const size_t numCols)
388  {
390  Kokkos::parallel_for
391  ("Tpetra::MultiVector unpack (Sacado::MP::Vector) (nonconstant stride)",
393  UnpackArrayMultiColumnVariableStride (execSpace, dst, src, idx, col, op, numCols));
394  }
395  else {
396  Kokkos::parallel_for
397  ("Tpetra::MultiVector unpack (Sacado::MP::Vector) (nonconstant stride)",
398  Kokkos::RangePolicy<execution_space, size_type> (0, idx.size ()),
399  UnpackArrayMultiColumnVariableStride (execSpace, dst, src, idx, col, op, numCols));
400  }
401  }
402  };
403 
404  template <typename DS, typename ... DP,
405  typename SS, typename ... SP,
406  typename DstIdxView, typename SrcIdxView>
407  struct PermuteArrayMultiColumn<
408  Kokkos::View<Sacado::MP::Vector<DS>**,DP...>,
409  Kokkos::View<const Sacado::MP::Vector<SS>**,SP...>,
410  DstIdxView,
411  SrcIdxView>
412  {
413  typedef Kokkos::View<Sacado::MP::Vector<DS>**,DP...> DstView;
414  typedef Kokkos::View<const Sacado::MP::Vector<SS>**,SP...> SrcView;
416  typedef typename execution_space::size_type size_type;
417 
420  DstIdxView dst_idx;
421  SrcIdxView src_idx;
422  size_t numCols;
423 
425  const SrcView& src_,
426  const DstIdxView& dst_idx_,
427  const SrcIdxView& src_idx_,
428  size_t numCols_) :
429  dst(dst_), src(src_), dst_idx(dst_idx_), src_idx(src_idx_),
430  numCols(numCols_) {}
431 
432  KOKKOS_INLINE_FUNCTION
433  void operator()( const size_type k ) const {
434  const typename DstIdxView::value_type toRow = dst_idx(k);
435  const typename SrcIdxView::value_type fromRow = src_idx(k);
436  for (size_t j = 0; j < numCols; ++j)
437  dst(toRow, j) = src(fromRow, j);
438  }
439 
440  KOKKOS_INLINE_FUNCTION
441  void operator()( const size_type k, const size_type tidx ) const {
442  const typename DstIdxView::value_type toRow = dst_idx(k);
443  const typename SrcIdxView::value_type fromRow = src_idx(k);
444  for (size_t j = 0; j < numCols; ++j)
445  dst(toRow, j).fastAccessCoeff(tidx) =
446  src(fromRow, j).fastAccessCoeff(tidx);
447  }
448 
449  static void permute(const DstView& dst,
450  const SrcView& src,
451  const DstIdxView& dst_idx,
452  const SrcIdxView& src_idx,
453  size_t numCols) {
454  const size_type n = std::min( dst_idx.size(), src_idx.size() );
456  Kokkos::parallel_for(
458  PermuteArrayMultiColumn(dst,src,dst_idx,src_idx,numCols) );
459  else
460  Kokkos::parallel_for(
461  n, PermuteArrayMultiColumn(dst,src,dst_idx,src_idx,numCols) );
462  }
463  };
464 
465  template <typename DS, typename ... DP,
466  typename SS, typename ... SP,
467  typename DstIdxView, typename SrcIdxView,
468  typename DstColView, typename SrcColView>
469  struct PermuteArrayMultiColumnVariableStride<
470  Kokkos::View<Sacado::MP::Vector<DS>**,DP...>,
471  Kokkos::View<const Sacado::MP::Vector<SS>**,SP...>,
472  DstIdxView, SrcIdxView,
473  DstColView, SrcColView >
474  {
475  typedef Kokkos::View<Sacado::MP::Vector<DS>**,DP...> DstView;
476  typedef Kokkos::View<const Sacado::MP::Vector<SS>**,SP...> SrcView;
478  typedef typename execution_space::size_type size_type;
479 
482  DstIdxView dst_idx;
483  SrcIdxView src_idx;
484  DstColView dst_col;
485  SrcColView src_col;
486  size_t numCols;
487 
489  const SrcView& src_,
490  const DstIdxView& dst_idx_,
491  const SrcIdxView& src_idx_,
492  const DstColView& dst_col_,
493  const SrcColView& src_col_,
494  size_t numCols_) :
495  dst(dst_), src(src_), dst_idx(dst_idx_), src_idx(src_idx_),
496  dst_col(dst_col_), src_col(src_col_),
497  numCols(numCols_) {}
498 
499  KOKKOS_INLINE_FUNCTION
500  void operator()( const size_type k ) const {
501  const typename DstIdxView::value_type toRow = dst_idx(k);
502  const typename SrcIdxView::value_type fromRow = src_idx(k);
503  for (size_t j = 0; j < numCols; ++j)
504  dst(toRow, dst_col(j)) = src(fromRow, src_col(j));
505  }
506 
507  KOKKOS_INLINE_FUNCTION
508  void operator()( const size_type k, const size_type tidx ) const {
509  const typename DstIdxView::value_type toRow = dst_idx(k);
510  const typename SrcIdxView::value_type fromRow = src_idx(k);
511  for (size_t j = 0; j < numCols; ++j)
512  dst(toRow, dst_col(j)).fastAccessCoeff(tidx) =
513  src(fromRow, src_col(j)).fastAccessCoeff(tidx);
514  }
515 
516  static void permute(const DstView& dst,
517  const SrcView& src,
518  const DstIdxView& dst_idx,
519  const SrcIdxView& src_idx,
520  const DstColView& dst_col,
521  const SrcColView& src_col,
522  size_t numCols) {
523  const size_type n = std::min( dst_idx.size(), src_idx.size() );
525  Kokkos::parallel_for(
527  PermuteArrayMultiColumnVariableStride(
528  dst,src,dst_idx,src_idx,dst_col,src_col,numCols) );
529  else
530  Kokkos::parallel_for(
531  n, PermuteArrayMultiColumnVariableStride(
532  dst,src,dst_idx,src_idx,dst_col,src_col,numCols) );
533  }
534  };
535 
536 } // Details namespace
537 } // KokkosRefactor namespace
538 } // Tpetra namespace
539 
540 #endif // TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_DIST_OBJECT_KERNELS_MP_VECTOR_HPP
Kokkos::DefaultExecutionSpace execution_space
expr1 expr1 expr1 expr2 expr1 expr1 c expr2 expr1 c fastAccessCoeff(j) - expr2.val(j)
KOKKOS_INLINE_FUNCTION PCE< Storage > min(const typename PCE< Storage >::value_type &a, const PCE< Storage > &b)
KOKKOS_INLINE_FUNCTION constexpr std::enable_if< is_view_uq_pce< View< T, P... > >::value, unsigned >::type dimension_scalar(const View< T, P... > &view)
Team-based parallel work configuration for Sacado::MP::Vector.