Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
Tpetra_KokkosRefactor_Details_MultiVectorDistObjectKernels_UQ_PCE.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_UQ_PCE_HPP
45 #define TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_DIST_OBJECT_KERNELS_UQ_PCE_HPP
46 
47 //----------------------------------------------------------------------------
48 // Specializations of Tpetra::MultiVector pack/unpack kernels for UQ::PCE
49 //----------------------------------------------------------------------------
50 
51 #include "Tpetra_KokkosRefactor_Details_MultiVectorDistObjectKernels.hpp"
54 
55 // For device_is_cuda<>
57 
58 namespace Tpetra {
59 namespace KokkosRefactor {
60 namespace Details {
61 
62  // Functors for implementing packAndPrepare and unpackAndCombine
63  // through parallel_for
64 
65  template <typename DS, typename ... DP,
66  typename SS, typename ... SP,
67  typename IdxView>
68  struct PackArraySingleColumn<
69  Kokkos::View<Sacado::UQ::PCE<DS>*,DP...>,
70  Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>,
71  IdxView >
72  {
73  typedef Kokkos::View<Sacado::UQ::PCE<DS>*,DP...> DstView;
74  typedef Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...> SrcView;
76  typedef typename execution_space::size_type size_type;
77 
78  static const unsigned BlockSize = 32;
79 
82  IdxView idx;
83  size_t col;
84 
86  const SrcView& src_,
87  const IdxView& idx_,
88  size_t col_) :
89  dst(dst_), src(src_), idx(idx_), col(col_) {}
90 
91  KOKKOS_INLINE_FUNCTION
92  void operator()( const size_type k ) const {
93  dst(k) = src(idx(k), col);
94  }
95 
96  KOKKOS_INLINE_FUNCTION
97  void operator()( const size_type k, const size_type tidx ) const {
98  for (size_type i=tidx; i<Kokkos::dimension_scalar(dst); i+=BlockSize)
99  dst(k).fastAccessCoeff(i) = src(idx(k), col).fastAccessCoeff(i);
100  }
101 
102  static void pack(const DstView& dst,
103  const SrcView& src,
104  const IdxView& idx,
105  size_t col) {
107  Kokkos::parallel_for(
108  Kokkos::MPVectorWorkConfig<execution_space>( idx.size(), BlockSize ),
109  PackArraySingleColumn(dst,src,idx,col) );
110  else
111  Kokkos::parallel_for( idx.size(),
112  PackArraySingleColumn(dst,src,idx,col) );
113  }
114  };
115 
116  template <typename DS, typename ... DP,
117  typename SS, typename ... SP,
118  typename IdxView>
119  struct PackArrayMultiColumn<
120  Kokkos::View<Sacado::UQ::PCE<DS>*,DP...>,
121  Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>,
122  IdxView >
123  {
124  typedef Kokkos::View<Sacado::UQ::PCE<DS>*,DP...> DstView;
125  typedef Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...> SrcView;
127  typedef typename execution_space::size_type size_type;
128 
129  static const unsigned BlockSize = 32;
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  for (size_type i=tidx; i<Kokkos::dimension_scalar(dst); i+=BlockSize)
156  dst(offset + j).fastAccessCoeff(i) =
157  src(localRow, j).fastAccessCoeff(i);
158  }
159 
160  static void pack(const DstView& dst,
161  const SrcView& src,
162  const IdxView& idx,
163  size_t numCols) {
165  Kokkos::parallel_for(
166  Kokkos::MPVectorWorkConfig<execution_space>( idx.size(), BlockSize ),
167  PackArrayMultiColumn(dst,src,idx,numCols) );
168  else
169  Kokkos::parallel_for( idx.size(),
170  PackArrayMultiColumn(dst,src,idx,numCols) );
171  }
172  };
173 
174  template <typename DS, typename ... DP,
175  typename SS, typename ... SP,
176  typename IdxView,
177  typename ColView>
178  struct PackArrayMultiColumnVariableStride<
179  Kokkos::View<Sacado::UQ::PCE<DS>*,DP...>,
180  Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>,
181  IdxView,
182  ColView>
183  {
184  typedef Kokkos::View<Sacado::UQ::PCE<DS>*,DP...> DstView;
185  typedef Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...> SrcView;
187  typedef typename execution_space::size_type size_type;
188 
189  static const unsigned BlockSize = 32;
190 
193  IdxView idx;
194  ColView col;
195  size_t numCols;
196 
198  const SrcView& src_,
199  const IdxView& idx_,
200  const ColView& col_,
201  size_t numCols_) :
202  dst(dst_), src(src_), idx(idx_), col(col_), numCols(numCols_) {}
203 
204  KOKKOS_INLINE_FUNCTION
205  void operator()( const size_type k ) const {
206  const typename IdxView::value_type localRow = idx(k);
207  const size_t offset = k*numCols;
208  for (size_t j = 0; j < numCols; ++j)
209  dst(offset + j) = src(localRow, col(j));
210  }
211 
212  KOKKOS_INLINE_FUNCTION
213  void operator()( const size_type k, const size_type tidx ) const {
214  const typename IdxView::value_type localRow = idx(k);
215  const size_t offset = k*numCols;
216  for (size_t j = 0; j < numCols; ++j)
217  for (size_type i=tidx; i<Kokkos::dimension_scalar(dst); i+=BlockSize)
218  dst(offset + j).fastAccessCoeff(i) =
219  src(localRow, col(j)).fastAccessCoeff(i);
220  }
221 
222  static void pack(const DstView& dst,
223  const SrcView& src,
224  const IdxView& idx,
225  const ColView& col,
226  size_t numCols) {
228  Kokkos::parallel_for(
229  Kokkos::MPVectorWorkConfig<execution_space>( idx.size(), BlockSize ),
230  PackArrayMultiColumnVariableStride(dst,src,idx,col,numCols) );
231  else
232  Kokkos::parallel_for( idx.size(),
233  PackArrayMultiColumnVariableStride(
234  dst,src,idx,col,numCols) );
235  }
236  };
237 
238  template <typename DS, typename ... DP,
239  typename SS, typename ... SP,
240  typename IdxView, typename Op>
241  struct UnpackArrayMultiColumn<
242  Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>,
243  Kokkos::View<const Sacado::UQ::PCE<SS>*,SP...>,
244  IdxView,
245  Op >
246  {
247  typedef Kokkos::View<Sacado::UQ::PCE<DS>**,DP...> DstView;
248  typedef Kokkos::View<const Sacado::UQ::PCE<SS>*,SP...> SrcView;
250  typedef typename execution_space::size_type size_type;
251 
252  static const unsigned BlockSize = 32;
253 
256  IdxView idx;
257  Op op;
258  size_t numCols;
259 
261  const SrcView& src_,
262  const IdxView& idx_,
263  const Op& op_,
264  size_t numCols_) :
265  dst(dst_), src(src_), idx(idx_), op(op_), numCols(numCols_) {}
266 
267  KOKKOS_INLINE_FUNCTION
268  void operator()( const size_type k ) const {
269  const typename IdxView::value_type localRow = idx(k);
270  const size_t offset = k*numCols;
271  for (size_t j = 0; j < numCols; ++j)
272  op( dst(localRow,j), src(offset+j) );
273  }
274 
275  KOKKOS_INLINE_FUNCTION
276  void operator()( const size_type k, const size_type tidx ) const {
277  const typename IdxView::value_type localRow = idx(k);
278  const size_t offset = k*numCols;
279  for (size_t j = 0; j < numCols; ++j)
280  for (size_type i=tidx; i<Kokkos::dimension_scalar(dst); i+=BlockSize)
281  op( dst(localRow,j).fastAccessCoeff(i),
282  src(offset+j).fastAccessCoeff(i) );
283  }
284 
285  static void unpack(const DstView& dst,
286  const SrcView& src,
287  const IdxView& idx,
288  const Op& op,
289  size_t numCols) {
291  Kokkos::parallel_for(
292  Kokkos::MPVectorWorkConfig<execution_space>( idx.size(), BlockSize ),
293  UnpackArrayMultiColumn(dst,src,idx,op,numCols) );
294  else
295  Kokkos::parallel_for( idx.size(),
296  UnpackArrayMultiColumn(dst,src,idx,op,numCols) );
297  }
298  };
299 
300  template <typename DS, typename ... DP,
301  typename SS, typename ... SP,
302  typename IdxView, typename ColView, typename Op>
303  struct UnpackArrayMultiColumnVariableStride<
304  Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>,
305  Kokkos::View<const Sacado::UQ::PCE<SS>*,SP...>,
306  IdxView,
307  ColView,
308  Op>
309  {
310  typedef Kokkos::View<Sacado::UQ::PCE<DS>**,DP...> DstView;
311  typedef Kokkos::View<const Sacado::UQ::PCE<SS>*,SP...> SrcView;
313  typedef typename execution_space::size_type size_type;
314 
315  static const unsigned BlockSize = 32;
316 
319  IdxView idx;
320  ColView col;
321  Op op;
322  size_t numCols;
323 
325  const SrcView& src_,
326  const IdxView& idx_,
327  const ColView& col_,
328  const Op& op_,
329  size_t numCols_) :
330  dst(dst_), src(src_), idx(idx_), col(col_), op(op_), numCols(numCols_) {}
331 
332  KOKKOS_INLINE_FUNCTION
333  void operator()( const size_type k ) const {
334  const typename IdxView::value_type localRow = idx(k);
335  const size_t offset = k*numCols;
336  for (size_t j = 0; j < numCols; ++j)
337  op( dst(localRow,col(j)), src(offset+j) );
338  }
339 
340  KOKKOS_INLINE_FUNCTION
341  void operator()( const size_type k, const size_type tidx ) const {
342  const typename IdxView::value_type localRow = idx(k);
343  const size_t offset = k*numCols;
344  for (size_t j = 0; j < numCols; ++j)
345  for (size_type i=tidx; i<Kokkos::dimension_scalar(dst); i+=BlockSize)
346  op( dst(localRow,col(j)).fastAccessCoeff(i),
347  src(offset+j).fastAccessCoeff(i) );
348  }
349 
350  static void unpack(const DstView& dst,
351  const SrcView& src,
352  const IdxView& idx,
353  const ColView& col,
354  const Op& op,
355  size_t numCols) {
357  Kokkos::parallel_for(
358  Kokkos::MPVectorWorkConfig<execution_space>( idx.size(), BlockSize ),
359  UnpackArrayMultiColumnVariableStride(dst,src,idx,col,op,numCols) );
360  else
361  Kokkos::parallel_for( idx.size(),
362  UnpackArrayMultiColumnVariableStride(
363  dst,src,idx,col,op,numCols) );
364  }
365  };
366 
367  template <typename DS, typename ... DP,
368  typename SS, typename ... SP,
369  typename DstIdxView, typename SrcIdxView>
370  struct PermuteArrayMultiColumn<
371  Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>,
372  Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>,
373  DstIdxView,
374  SrcIdxView>
375  {
376  typedef Kokkos::View<Sacado::UQ::PCE<DS>**,DP...> DstView;
377  typedef Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...> SrcView;
379  typedef typename execution_space::size_type size_type;
380 
381  static const unsigned BlockSize = 32;
382 
385  DstIdxView dst_idx;
386  SrcIdxView src_idx;
387  size_t numCols;
388 
390  const SrcView& src_,
391  const DstIdxView& dst_idx_,
392  const SrcIdxView& src_idx_,
393  size_t numCols_) :
394  dst(dst_), src(src_), dst_idx(dst_idx_), src_idx(src_idx_),
395  numCols(numCols_) {}
396 
397  KOKKOS_INLINE_FUNCTION
398  void operator()( const size_type k ) const {
399  const typename DstIdxView::value_type toRow = dst_idx(k);
400  const typename SrcIdxView::value_type fromRow = src_idx(k);
401  for (size_t j = 0; j < numCols; ++j)
402  dst(toRow, j) = src(fromRow, j);
403  }
404 
405  KOKKOS_INLINE_FUNCTION
406  void operator()( const size_type k, const size_type tidx ) const {
407  const typename DstIdxView::value_type toRow = dst_idx(k);
408  const typename SrcIdxView::value_type fromRow = src_idx(k);
409  for (size_t j = 0; j < numCols; ++j)
410  for (size_type i=tidx; i<Kokkos::dimension_scalar(dst); i+=BlockSize)
411  dst(toRow, j).fastAccessCoeff(i) =
412  src(fromRow, j).fastAccessCoeff(i);
413  }
414 
415  static void permute(const DstView& dst,
416  const SrcView& src,
417  const DstIdxView& dst_idx,
418  const SrcIdxView& src_idx,
419  size_t numCols) {
420  const size_type n = std::min( dst_idx.size(), src_idx.size() );
422  Kokkos::parallel_for(
424  PermuteArrayMultiColumn(dst,src,dst_idx,src_idx,numCols) );
425  else
426  Kokkos::parallel_for(
427  n, PermuteArrayMultiColumn(dst,src,dst_idx,src_idx,numCols) );
428  }
429  };
430 
431  template <typename DS, typename ... DP,
432  typename SS, typename ... SP,
433  typename DstIdxView, typename SrcIdxView,
434  typename DstColView, typename SrcColView>
435  struct PermuteArrayMultiColumnVariableStride<
436  Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>,
437  Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>,
438  DstIdxView, SrcIdxView,
439  DstColView, SrcColView >
440  {
441  typedef Kokkos::View<Sacado::UQ::PCE<DS>**,DP...> DstView;
442  typedef Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...> SrcView;
444  typedef typename execution_space::size_type size_type;
445 
446  static const unsigned BlockSize = 32;
447 
450  DstIdxView dst_idx;
451  SrcIdxView src_idx;
452  DstColView dst_col;
453  SrcColView src_col;
454  size_t numCols;
455 
457  const SrcView& src_,
458  const DstIdxView& dst_idx_,
459  const SrcIdxView& src_idx_,
460  const DstColView& dst_col_,
461  const SrcColView& src_col_,
462  size_t numCols_) :
463  dst(dst_), src(src_), dst_idx(dst_idx_), src_idx(src_idx_),
464  dst_col(dst_col_), src_col(src_col_),
465  numCols(numCols_) {}
466 
467  KOKKOS_INLINE_FUNCTION
468  void operator()( const size_type k ) const {
469  const typename DstIdxView::value_type toRow = dst_idx(k);
470  const typename SrcIdxView::value_type fromRow = src_idx(k);
471  for (size_t j = 0; j < numCols; ++j)
472  dst(toRow, dst_col(j)) = src(fromRow, src_col(j));
473  }
474 
475  KOKKOS_INLINE_FUNCTION
476  void operator()( const size_type k, const size_type tidx ) const {
477  const typename DstIdxView::value_type toRow = dst_idx(k);
478  const typename SrcIdxView::value_type fromRow = src_idx(k);
479  for (size_t j = 0; j < numCols; ++j)
480  for (size_type i=tidx; i<Kokkos::dimension_scalar(dst); i+=BlockSize)
481  dst(toRow, dst_col(j)).fastAccessCoeff(i) =
482  src(fromRow, src_col(j)).fastAccessCoeff(i);
483  }
484 
485  static void permute(const DstView& dst,
486  const SrcView& src,
487  const DstIdxView& dst_idx,
488  const SrcIdxView& src_idx,
489  const DstColView& dst_col,
490  const SrcColView& src_col,
491  size_t numCols) {
492  const size_type n = std::min( dst_idx.size(), src_idx.size() );
494  Kokkos::parallel_for(
496  PermuteArrayMultiColumnVariableStride(
497  dst,src,dst_idx,src_idx,dst_col,src_col,numCols) );
498  else
499  Kokkos::parallel_for(
500  n, PermuteArrayMultiColumnVariableStride(
501  dst,src,dst_idx,src_idx,dst_col,src_col,numCols) );
502  }
503  };
504 
505 } // Details namespace
506 } // KokkosRefactor namespace
507 } // Tpetra namespace
508 
509 #endif // TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_DIST_OBJECT_KERNELS_UQ_PCE_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.