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_HAVE_CUDA )
60  template< class D >
61  struct device_is_cuda : public Kokkos::Impl::is_same<D,Kokkos::Cuda> {};
62 #else
63  template< class D >
64  struct device_is_cuda : public Kokkos::Impl::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 DS, typename ... DP,
235  typename SS, typename ... SP,
236  typename IdxView, typename Op>
237  struct UnpackArrayMultiColumn<
238  Kokkos::View<Sacado::MP::Vector<DS>**,DP...>,
239  Kokkos::View<const Sacado::MP::Vector<SS>*,SP...>,
240  IdxView,
241  Op >
242  {
243  typedef Kokkos::View<Sacado::MP::Vector<DS>**,DP...> DstView;
244  typedef Kokkos::View<const Sacado::MP::Vector<SS>*,SP...> SrcView;
246  typedef typename execution_space::size_type size_type;
247 
250  IdxView idx;
251  Op op;
252  size_t numCols;
253 
255  const SrcView& src_,
256  const IdxView& idx_,
257  const Op& op_,
258  size_t numCols_) :
259  dst(dst_), src(src_), idx(idx_), op(op_), numCols(numCols_) {}
260 
261  KOKKOS_INLINE_FUNCTION
262  void operator()( const size_type k ) const {
263  const typename IdxView::value_type localRow = idx(k);
264  const size_t offset = k*numCols;
265  for (size_t j = 0; j < numCols; ++j)
266  op( dst(localRow,j), src(offset+j) );
267  }
268 
269  KOKKOS_INLINE_FUNCTION
270  void operator()( const size_type k, const size_type tidx ) const {
271  const typename IdxView::value_type localRow = idx(k);
272  const size_t offset = k*numCols;
273  for (size_t j = 0; j < numCols; ++j)
274  op( dst(localRow,j).fastAccessCoeff(tidx),
275  src(offset+j).fastAccessCoeff(tidx) );
276  }
277 
278  static void unpack(const DstView& dst,
279  const SrcView& src,
280  const IdxView& idx,
281  const Op& op,
282  size_t numCols) {
284  Kokkos::parallel_for(
286  UnpackArrayMultiColumn(dst,src,idx,op,numCols) );
287  else
288  Kokkos::parallel_for( idx.size(),
289  UnpackArrayMultiColumn(dst,src,idx,op,numCols) );
290  }
291  };
292 
293  template <typename DS, typename ... DP,
294  typename SS, typename ... SP,
295  typename IdxView, typename ColView, typename Op>
296  struct UnpackArrayMultiColumnVariableStride<
297  Kokkos::View<Sacado::MP::Vector<DS>**,DP...>,
298  Kokkos::View<const Sacado::MP::Vector<SS>*,SP...>,
299  IdxView,
300  ColView,
301  Op>
302  {
303  typedef Kokkos::View<Sacado::MP::Vector<DS>**,DP...> DstView;
304  typedef Kokkos::View<const Sacado::MP::Vector<SS>*,SP...> SrcView;
306  typedef typename execution_space::size_type size_type;
307 
310  IdxView idx;
311  ColView col;
312  Op op;
313  size_t numCols;
314 
316  const SrcView& src_,
317  const IdxView& idx_,
318  const ColView& col_,
319  const Op& op_,
320  size_t numCols_) :
321  dst(dst_), src(src_), idx(idx_), col(col_), op(op_), numCols(numCols_) {}
322 
323  KOKKOS_INLINE_FUNCTION
324  void operator()( const size_type k ) const {
325  const typename IdxView::value_type localRow = idx(k);
326  const size_t offset = k*numCols;
327  for (size_t j = 0; j < numCols; ++j)
328  op( dst(localRow,col(j)), src(offset+j) );
329  }
330 
331  KOKKOS_INLINE_FUNCTION
332  void operator()( const size_type k, const size_type tidx ) const {
333  const typename IdxView::value_type localRow = idx(k);
334  const size_t offset = k*numCols;
335  for (size_t j = 0; j < numCols; ++j)
336  op( dst(localRow,col(j)).fastAccessCoeff(tidx),
337  src(offset+j).fastAccessCoeff(tidx) );
338  }
339 
340  static void unpack(const DstView& dst,
341  const SrcView& src,
342  const IdxView& idx,
343  const ColView& col,
344  const Op& op,
345  size_t numCols) {
347  Kokkos::parallel_for(
349  UnpackArrayMultiColumnVariableStride(dst,src,idx,col,op,numCols) );
350  else
351  Kokkos::parallel_for( idx.size(),
352  UnpackArrayMultiColumnVariableStride(
353  dst,src,idx,col,op,numCols) );
354  }
355  };
356 
357  template <typename DS, typename ... DP,
358  typename SS, typename ... SP,
359  typename DstIdxView, typename SrcIdxView>
360  struct PermuteArrayMultiColumn<
361  Kokkos::View<Sacado::MP::Vector<DS>**,DP...>,
362  Kokkos::View<const Sacado::MP::Vector<SS>**,SP...>,
363  DstIdxView,
364  SrcIdxView>
365  {
366  typedef Kokkos::View<Sacado::MP::Vector<DS>**,DP...> DstView;
367  typedef Kokkos::View<const Sacado::MP::Vector<SS>**,SP...> SrcView;
369  typedef typename execution_space::size_type size_type;
370 
373  DstIdxView dst_idx;
374  SrcIdxView src_idx;
375  size_t numCols;
376 
378  const SrcView& src_,
379  const DstIdxView& dst_idx_,
380  const SrcIdxView& src_idx_,
381  size_t numCols_) :
382  dst(dst_), src(src_), dst_idx(dst_idx_), src_idx(src_idx_),
383  numCols(numCols_) {}
384 
385  KOKKOS_INLINE_FUNCTION
386  void operator()( const size_type k ) const {
387  const typename DstIdxView::value_type toRow = dst_idx(k);
388  const typename SrcIdxView::value_type fromRow = src_idx(k);
389  for (size_t j = 0; j < numCols; ++j)
390  dst(toRow, j) = src(fromRow, j);
391  }
392 
393  KOKKOS_INLINE_FUNCTION
394  void operator()( const size_type k, const size_type tidx ) const {
395  const typename DstIdxView::value_type toRow = dst_idx(k);
396  const typename SrcIdxView::value_type fromRow = src_idx(k);
397  for (size_t j = 0; j < numCols; ++j)
398  dst(toRow, j).fastAccessCoeff(tidx) =
399  src(fromRow, j).fastAccessCoeff(tidx);
400  }
401 
402  static void permute(const DstView& dst,
403  const SrcView& src,
404  const DstIdxView& dst_idx,
405  const SrcIdxView& src_idx,
406  size_t numCols) {
407  const size_type n = std::min( dst_idx.size(), src_idx.size() );
409  Kokkos::parallel_for(
411  PermuteArrayMultiColumn(dst,src,dst_idx,src_idx,numCols) );
412  else
413  Kokkos::parallel_for(
414  n, PermuteArrayMultiColumn(dst,src,dst_idx,src_idx,numCols) );
415  }
416  };
417 
418  template <typename DS, typename ... DP,
419  typename SS, typename ... SP,
420  typename DstIdxView, typename SrcIdxView,
421  typename DstColView, typename SrcColView>
422  struct PermuteArrayMultiColumnVariableStride<
423  Kokkos::View<Sacado::MP::Vector<DS>**,DP...>,
424  Kokkos::View<const Sacado::MP::Vector<SS>**,SP...>,
425  DstIdxView, SrcIdxView,
426  DstColView, SrcColView >
427  {
428  typedef Kokkos::View<Sacado::MP::Vector<DS>**,DP...> DstView;
429  typedef Kokkos::View<const Sacado::MP::Vector<SS>**,SP...> SrcView;
431  typedef typename execution_space::size_type size_type;
432 
435  DstIdxView dst_idx;
436  SrcIdxView src_idx;
437  DstColView dst_col;
438  SrcColView src_col;
439  size_t numCols;
440 
442  const SrcView& src_,
443  const DstIdxView& dst_idx_,
444  const SrcIdxView& src_idx_,
445  const DstColView& dst_col_,
446  const SrcColView& src_col_,
447  size_t numCols_) :
448  dst(dst_), src(src_), dst_idx(dst_idx_), src_idx(src_idx_),
449  dst_col(dst_col_), src_col(src_col_),
450  numCols(numCols_) {}
451 
452  KOKKOS_INLINE_FUNCTION
453  void operator()( const size_type k ) const {
454  const typename DstIdxView::value_type toRow = dst_idx(k);
455  const typename SrcIdxView::value_type fromRow = src_idx(k);
456  for (size_t j = 0; j < numCols; ++j)
457  dst(toRow, dst_col(j)) = src(fromRow, src_col(j));
458  }
459 
460  KOKKOS_INLINE_FUNCTION
461  void operator()( const size_type k, const size_type tidx ) const {
462  const typename DstIdxView::value_type toRow = dst_idx(k);
463  const typename SrcIdxView::value_type fromRow = src_idx(k);
464  for (size_t j = 0; j < numCols; ++j)
465  dst(toRow, dst_col(j)).fastAccessCoeff(tidx) =
466  src(fromRow, src_col(j)).fastAccessCoeff(tidx);
467  }
468 
469  static void permute(const DstView& dst,
470  const SrcView& src,
471  const DstIdxView& dst_idx,
472  const SrcIdxView& src_idx,
473  const DstColView& dst_col,
474  const SrcColView& src_col,
475  size_t numCols) {
476  const size_type n = std::min( dst_idx.size(), src_idx.size() );
478  Kokkos::parallel_for(
480  PermuteArrayMultiColumnVariableStride(
481  dst,src,dst_idx,src_idx,dst_col,src_col,numCols) );
482  else
483  Kokkos::parallel_for(
484  n, PermuteArrayMultiColumnVariableStride(
485  dst,src,dst_idx,src_idx,dst_col,src_col,numCols) );
486  }
487  };
488 
489 } // Details namespace
490 } // KokkosRefactor namespace
491 } // Tpetra namespace
492 
493 #endif // TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_DIST_OBJECT_KERNELS_MP_VECTOR_HPP
Kokkos::DefaultExecutionSpace execution_space
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)
expr expr expr expr j
Team-based parallel work configuration for Sacado::MP::Vector.