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 ExecutionSpace,
239  typename DS,
240  typename ... DP,
241  typename SS,
242  typename ... SP,
243  typename IdxView,
244  typename Op>
245  struct UnpackArrayMultiColumn<
246  ExecutionSpace,
247  Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>,
248  Kokkos::View<const Sacado::UQ::PCE<SS>*,SP...>,
249  IdxView,
250  Op >
251  {
252  typedef Kokkos::View<Sacado::UQ::PCE<DS>**,DP...> DstView;
253  typedef Kokkos::View<const Sacado::UQ::PCE<SS>*,SP...> SrcView;
255  typedef typename execution_space::size_type size_type;
256 
257  static const unsigned BlockSize = 32;
258 
261  IdxView idx;
262  Op op;
263  size_t numCols;
264 
265  UnpackArrayMultiColumn (const ExecutionSpace& /* execSpace */,
266  const DstView& dst_,
267  const SrcView& src_,
268  const IdxView& idx_,
269  const Op& op_,
270  const size_t numCols_) :
271  dst (dst_),
272  src (src_),
273  idx (idx_),
274  op (op_),
275  numCols (numCols_)
276  {}
277 
278  KOKKOS_INLINE_FUNCTION void
279  operator() (const size_type k) const
280  {
281  const typename IdxView::value_type localRow = idx(k);
282  const size_t offset = k*numCols;
283  for (size_t j = 0; j < numCols; ++j) {
284  op (dst(localRow, j), src(offset+j));
285  }
286  }
287 
288  KOKKOS_INLINE_FUNCTION void
289  operator() (const size_type k, const size_type tidx) const
290  {
291  const typename IdxView::value_type localRow = idx(k);
292  const size_t offset = k*numCols;
293  for (size_t j = 0; j < numCols; ++j) {
294  for (size_type i=tidx; i<Kokkos::dimension_scalar(dst); i+=BlockSize) {
295  op (dst(localRow, j).fastAccessCoeff(i),
296  src(offset+j).fastAccessCoeff(i));
297  }
298  }
299  }
300 
301  static void
302  unpack (const ExecutionSpace& execSpace,
303  const DstView& dst,
304  const SrcView& src,
305  const IdxView& idx,
306  const Op& op,
307  const size_t numCols)
308  {
310  Kokkos::parallel_for
311  ("Tpetra::MultiVector (Sacado::UQ::PCE) unpack (constant stride)",
312  Kokkos::MPVectorWorkConfig<execution_space> (idx.size (), BlockSize),
313  UnpackArrayMultiColumn (execSpace, dst, src, idx, op, numCols));
314  }
315  else {
316  Kokkos::parallel_for
317  ("Tpetra::MultiVector (Sacado::UQ::PCE) unpack (constant stride)",
318  Kokkos::RangePolicy<execution_space, size_type> (0, idx.size ()),
319  UnpackArrayMultiColumn (execSpace, dst, src, idx, op, numCols));
320  }
321  }
322  };
323 
324  template <typename ExecutionSpace,
325  typename DS,
326  typename ... DP,
327  typename SS,
328  typename ... SP,
329  typename IdxView,
330  typename ColView,
331  typename Op>
332  struct UnpackArrayMultiColumnVariableStride<
333  ExecutionSpace,
334  Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>,
335  Kokkos::View<const Sacado::UQ::PCE<SS>*,SP...>,
336  IdxView,
337  ColView,
338  Op>
339  {
340  typedef Kokkos::View<Sacado::UQ::PCE<DS>**,DP...> DstView;
341  typedef Kokkos::View<const Sacado::UQ::PCE<SS>*,SP...> SrcView;
343  typedef typename execution_space::size_type size_type;
344 
345  static const unsigned BlockSize = 32;
346 
349  IdxView idx;
350  ColView col;
351  Op op;
352  size_t numCols;
353 
354  UnpackArrayMultiColumnVariableStride (const ExecutionSpace& /* execSpace */,
355  const DstView& dst_,
356  const SrcView& src_,
357  const IdxView& idx_,
358  const ColView& col_,
359  const Op& op_,
360  size_t numCols_) :
361  dst(dst_), src(src_), idx(idx_), col(col_), op(op_), numCols(numCols_) {}
362 
363  KOKKOS_INLINE_FUNCTION
364  void operator()( const size_type k ) const {
365  const typename IdxView::value_type localRow = idx(k);
366  const size_t offset = k*numCols;
367  for (size_t j = 0; j < numCols; ++j)
368  op( dst(localRow,col(j)), src(offset+j) );
369  }
370 
371  KOKKOS_INLINE_FUNCTION
372  void operator()( const size_type k, const size_type tidx ) const {
373  const typename IdxView::value_type localRow = idx(k);
374  const size_t offset = k*numCols;
375  for (size_t j = 0; j < numCols; ++j)
376  for (size_type i=tidx; i<Kokkos::dimension_scalar(dst); i+=BlockSize)
377  op( dst(localRow,col(j)).fastAccessCoeff(i),
378  src(offset+j).fastAccessCoeff(i) );
379  }
380 
381  static void
382  unpack (const ExecutionSpace& execSpace,
383  const DstView& dst,
384  const SrcView& src,
385  const IdxView& idx,
386  const ColView& col,
387  const Op& op,
388  const size_t numCols)
389  {
391  Kokkos::parallel_for
392  ("Tpetra::MultiVector unpack (Sacado::UQ::PCE) (nonconstant stride)",
393  Kokkos::MPVectorWorkConfig<execution_space> (idx.size (), BlockSize),
394  UnpackArrayMultiColumnVariableStride (execSpace, dst, src,
395  idx, col, op, numCols));
396  }
397  else {
398  Kokkos::parallel_for
399  ("Tpetra::MultiVector unpack (Sacado::UQ::PCE) (nonconstant stride)",
400  Kokkos::RangePolicy<execution_space, size_type> (0, idx.size ()),
401  UnpackArrayMultiColumnVariableStride (execSpace, dst, src,
402  idx, col, op, numCols));
403  }
404  }
405  };
406 
407  template <typename DS, typename ... DP,
408  typename SS, typename ... SP,
409  typename DstIdxView, typename SrcIdxView>
410  struct PermuteArrayMultiColumn<
411  Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>,
412  Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>,
413  DstIdxView,
414  SrcIdxView>
415  {
416  typedef Kokkos::View<Sacado::UQ::PCE<DS>**,DP...> DstView;
417  typedef Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...> SrcView;
419  typedef typename execution_space::size_type size_type;
420 
421  static const unsigned BlockSize = 32;
422 
425  DstIdxView dst_idx;
426  SrcIdxView src_idx;
427  size_t numCols;
428 
430  const SrcView& src_,
431  const DstIdxView& dst_idx_,
432  const SrcIdxView& src_idx_,
433  size_t numCols_) :
434  dst(dst_), src(src_), dst_idx(dst_idx_), src_idx(src_idx_),
435  numCols(numCols_) {}
436 
437  KOKKOS_INLINE_FUNCTION
438  void operator()( const size_type k ) const {
439  const typename DstIdxView::value_type toRow = dst_idx(k);
440  const typename SrcIdxView::value_type fromRow = src_idx(k);
441  for (size_t j = 0; j < numCols; ++j)
442  dst(toRow, j) = src(fromRow, j);
443  }
444 
445  KOKKOS_INLINE_FUNCTION
446  void operator()( const size_type k, const size_type tidx ) const {
447  const typename DstIdxView::value_type toRow = dst_idx(k);
448  const typename SrcIdxView::value_type fromRow = src_idx(k);
449  for (size_t j = 0; j < numCols; ++j)
450  for (size_type i=tidx; i<Kokkos::dimension_scalar(dst); i+=BlockSize)
451  dst(toRow, j).fastAccessCoeff(i) =
452  src(fromRow, j).fastAccessCoeff(i);
453  }
454 
455  static void permute(const DstView& dst,
456  const SrcView& src,
457  const DstIdxView& dst_idx,
458  const SrcIdxView& src_idx,
459  size_t numCols) {
460  const size_type n = std::min( dst_idx.size(), src_idx.size() );
462  Kokkos::parallel_for(
464  PermuteArrayMultiColumn(dst,src,dst_idx,src_idx,numCols) );
465  else
466  Kokkos::parallel_for(
467  n, PermuteArrayMultiColumn(dst,src,dst_idx,src_idx,numCols) );
468  }
469  };
470 
471  template <typename DS, typename ... DP,
472  typename SS, typename ... SP,
473  typename DstIdxView, typename SrcIdxView,
474  typename DstColView, typename SrcColView>
475  struct PermuteArrayMultiColumnVariableStride<
476  Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>,
477  Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>,
478  DstIdxView, SrcIdxView,
479  DstColView, SrcColView >
480  {
481  typedef Kokkos::View<Sacado::UQ::PCE<DS>**,DP...> DstView;
482  typedef Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...> SrcView;
484  typedef typename execution_space::size_type size_type;
485 
486  static const unsigned BlockSize = 32;
487 
490  DstIdxView dst_idx;
491  SrcIdxView src_idx;
492  DstColView dst_col;
493  SrcColView src_col;
494  size_t numCols;
495 
497  const SrcView& src_,
498  const DstIdxView& dst_idx_,
499  const SrcIdxView& src_idx_,
500  const DstColView& dst_col_,
501  const SrcColView& src_col_,
502  size_t numCols_) :
503  dst(dst_), src(src_), dst_idx(dst_idx_), src_idx(src_idx_),
504  dst_col(dst_col_), src_col(src_col_),
505  numCols(numCols_) {}
506 
507  KOKKOS_INLINE_FUNCTION
508  void operator()( const size_type k ) 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)) = src(fromRow, src_col(j));
513  }
514 
515  KOKKOS_INLINE_FUNCTION
516  void operator()( const size_type k, const size_type tidx ) const {
517  const typename DstIdxView::value_type toRow = dst_idx(k);
518  const typename SrcIdxView::value_type fromRow = src_idx(k);
519  for (size_t j = 0; j < numCols; ++j)
520  for (size_type i=tidx; i<Kokkos::dimension_scalar(dst); i+=BlockSize)
521  dst(toRow, dst_col(j)).fastAccessCoeff(i) =
522  src(fromRow, src_col(j)).fastAccessCoeff(i);
523  }
524 
525  static void permute(const DstView& dst,
526  const SrcView& src,
527  const DstIdxView& dst_idx,
528  const SrcIdxView& src_idx,
529  const DstColView& dst_col,
530  const SrcColView& src_col,
531  size_t numCols) {
532  const size_type n = std::min( dst_idx.size(), src_idx.size() );
534  Kokkos::parallel_for(
536  PermuteArrayMultiColumnVariableStride(
537  dst,src,dst_idx,src_idx,dst_col,src_col,numCols) );
538  else
539  Kokkos::parallel_for(
540  n, PermuteArrayMultiColumnVariableStride(
541  dst,src,dst_idx,src_idx,dst_col,src_col,numCols) );
542  }
543  };
544 
545 } // Details namespace
546 } // KokkosRefactor namespace
547 } // Tpetra namespace
548 
549 #endif // TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_DIST_OBJECT_KERNELS_UQ_PCE_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.