Tpetra parallel linear algebra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Tpetra_Details_fill.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Tpetra: Templated Linear Algebra Services Package
5 // Copyright (2008) Sandia Corporation
6 //
7 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8 // the U.S. Government retains certain rights in this software.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ************************************************************************
40 // @HEADER
41 
42 #ifndef TPETRA_DETAILS_FILL_HPP
43 #define TPETRA_DETAILS_FILL_HPP
44 
56 
57 #include "Tpetra_Details_Blas.hpp"
58 #include <type_traits>
59 
60 namespace Tpetra {
61 namespace Details {
62 namespace Blas {
63 namespace Impl {
64 
66 void*
67 memsetWrapper (void* dest, int ch, std::size_t count);
68 
70 template<class ViewType,
71  class ValueType,
72  class ExecutionSpace,
73  class IndexType,
74  const int rank = ViewType::Rank>
75 class Fill {
76 public:
77  static void
78  fill (const ExecutionSpace& execSpace,
79  const ViewType& X,
80  const ValueType& alpha,
81  const IndexType /* numRows */,
82  const IndexType /* numCols */ )
83  {
84  static_assert (std::is_integral<IndexType>::value,
85  "IndexType must be a built-in integer type.");
86  // DEEP_COPY REVIEW - NOT TESTED
87  Kokkos::deep_copy (execSpace, X, alpha);
88  }
89 };
90 
92 template<class ViewType,
93  class ValueType,
94  class ExecutionSpace,
95  class IndexType>
96 class Fill<ViewType, ValueType, ExecutionSpace, IndexType, 1> {
97 public:
98  Fill (const ViewType& X, const ValueType& alpha) :
99  X_ (X), alpha_ (alpha)
100  {
101  static_assert (ViewType::Rank == 1,
102  "ViewType must be a rank-1 Kokkos::View.");
103  static_assert (std::is_integral<IndexType>::value,
104  "IndexType must be a built-in integer type.");
105  }
106 
107  KOKKOS_INLINE_FUNCTION void
108  operator () (const IndexType& i) const
109  {
110  X_[i] = alpha_;
111  }
112 
113  static void
114  fill (const ExecutionSpace& execSpace,
115  const ViewType& X,
116  const ValueType& alpha,
117  const IndexType numRows,
118  const IndexType /* numCols */)
119  {
120  typedef Kokkos::RangePolicy<ExecutionSpace, IndexType> range_type;
121  Kokkos::parallel_for ("fill", range_type (0, numRows),
122  Fill (X, alpha));
123  }
124 
125 private:
126  ViewType X_;
127  ValueType alpha_;
128 };
129 
131 template<class ViewType,
132  class ValueType,
133  class ExecutionSpace,
134  class IndexType>
135 class Fill<ViewType, ValueType, ExecutionSpace, IndexType, 2> {
136 public:
137  Fill (const ViewType& X,
138  const ValueType& alpha,
139  const IndexType numCols) :
140  X_ (X), alpha_ (alpha), numCols_ (numCols)
141  {
142  static_assert (ViewType::Rank == 2,
143  "ViewType must be a rank-2 Kokkos::View.");
144  static_assert (std::is_integral<IndexType>::value,
145  "IndexType must be a built-in integer type.");
146  }
147 
148  KOKKOS_INLINE_FUNCTION void
149  operator () (const IndexType& i) const
150  {
151  for (IndexType j = 0; j < numCols_; ++j) {
152  X_(i,j) = alpha_;
153  }
154  }
155 
156  static void
157  fill (const ExecutionSpace& execSpace,
158  const ViewType& X,
159  const ValueType& alpha,
160  const IndexType numRows,
161  const IndexType numCols)
162  {
163  typedef Kokkos::RangePolicy<ExecutionSpace, IndexType> range_type;
164  Kokkos::parallel_for ("fill", range_type (0, numRows),
165  Fill (X, alpha, numCols));
166  }
167 
168 private:
169  ViewType X_;
170  ValueType alpha_;
171  IndexType numCols_;
172 };
173 
174 #if defined(KOKKOS_ENABLE_SERIAL)
175 template<class ViewType,
178  class ValueType,
179  class IndexType>
180 struct Fill<ViewType,
181  ValueType,
182  Kokkos::Serial,
183  IndexType,
184  1>
185 {
186  static void
187  fill (const Kokkos::Serial& /* execSpace */,
188  const ViewType& X,
189  const ValueType& alpha,
190  const IndexType numRows,
191  const IndexType /* numCols */ )
192  {
193  static_assert (ViewType::Rank == 1,
194  "ViewType must be a rank-1 Kokkos::View.");
195  static_assert (std::is_integral<IndexType>::value,
196  "IndexType must be a built-in integer type.");
197  using ::Tpetra::Details::Blas::BlasSupportsScalar;
198  typedef typename ViewType::non_const_value_type view_value_type;
199 
200  // Do sizeof(view_value_type) and taking the address of a
201  // value_type instance work correctly with memset?
202  constexpr bool podType = BlasSupportsScalar<view_value_type>::value;
203 
204  if (podType && X.span_is_contiguous () && alpha == ValueType (0.0)) {
205  memsetWrapper (X.data (), 0, X.span () * sizeof (view_value_type));
206  }
207  else {
208  for (IndexType k = 0; k < numRows; ++k) {
209  X[k] = alpha;
210  }
211  }
212  }
213 };
214 #endif // defined(KOKKOS_ENABLE_SERIAL)
215 
216 #if defined(KOKKOS_ENABLE_SERIAL)
217 template<class ViewType,
220  class ValueType,
221  class IndexType>
222 struct Fill<ViewType,
223  ValueType,
224  Kokkos::Serial,
225  IndexType,
226  2>
227 {
228  static void
229  fill (const Kokkos::Serial& /* execSpace */,
230  const ViewType& X,
231  const ValueType& alpha,
232  const IndexType /* numRows */,
233  const IndexType numCols)
234  {
235  static_assert (ViewType::Rank == 2,
236  "ViewType must be a rank-2 Kokkos::View.");
237  static_assert (std::is_integral<IndexType>::value,
238  "IndexType must be a built-in integer type.");
239  using ::Tpetra::Details::Blas::BlasSupportsScalar;
240  typedef typename ViewType::non_const_value_type view_value_type;
241  typedef typename ViewType::array_layout array_layout;
242 
243  // Do sizeof(view_value_type) and taking the address of a
244  // value_type instance work correctly with memset?
245  constexpr bool podType = BlasSupportsScalar<view_value_type>::value;
246 
247  if (podType && alpha == ValueType (0.0)) {
248  if (X.span_is_contiguous ()) {
249  memsetWrapper (X.data (), 0, X.span () * sizeof (view_value_type));
250  }
251  else if (std::is_same<array_layout, Kokkos::LayoutLeft>::value) {
252  // Tpetra::MultiVector needs to optimize for LayoutLeft.
253  for (IndexType j = 0; j < numCols; ++j) {
254  auto X_j = Kokkos::subview (X, Kokkos::ALL (), j);
255  memsetWrapper (X_j.data (), 0,
256  X_j.extent (0) * sizeof (view_value_type));
257  }
258  }
259  else {
260  // DEEP_COPY REVIEW - NOT TESTED
261  Kokkos::deep_copy (X, view_value_type (0.0));
262  }
263  }
264  else {
265  // DEEP_COPY REVIEW - VALUE-TO-DEVICE
266  using execution_space = typename ViewType::execution_space;
267  Kokkos::deep_copy (execution_space(), X, alpha);
268  }
269  }
270 };
271 #endif // defined(KOKKOS_ENABLE_SERIAL)
272 
273 } // namespace Impl
274 
275 //
276 // SKIP TO HERE FOR THE ACTUAL INTERFACE
277 //
278 
286 template<class ViewType,
287  class ValueType,
288  class IndexType,
289  class ExecutionSpace>
290 void
291 fill (const ExecutionSpace& execSpace,
292  const ViewType& X,
293  const ValueType& alpha,
294  const IndexType numRows,
295  const IndexType numCols)
296 {
297  static_assert (std::is_integral<IndexType>::value,
298  "IndexType must be a built-in integer type.");
299  typedef Impl::Fill<ViewType, ValueType, ExecutionSpace,
300  IndexType, ViewType::Rank> impl_type;
301  impl_type::fill (execSpace, X, alpha, numRows, numCols);
302 }
303 
304 template<class ViewType,
305  class ValueType,
306  class IndexType,
307  class ExecutionSpace>
308 void
309 fill (const ExecutionSpace& execSpace,
310  const ViewType& X,
311  const ValueType& alpha,
312  const IndexType numRows,
313  const IndexType numCols,
314  const size_t whichVectors[])
315 {
316  static_assert (ViewType::Rank == 2, "ViewType must be a rank-2 "
317  "Kokkos::View in order to call the \"whichVectors\" "
318  "specialization of fill.");
319  static_assert (std::is_integral<IndexType>::value,
320  "IndexType must be a built-in integer type.");
321  for (IndexType k = 0; k < numCols; ++k) {
322  const IndexType j = whichVectors[k];
323  auto X_j = Kokkos::subview (X, Kokkos::ALL (), j);
324  typedef decltype (X_j) one_d_view_type;
325  typedef Impl::Fill<one_d_view_type, ValueType, ExecutionSpace,
326  IndexType, 1> impl_type;
327  impl_type::fill (execSpace, X_j, alpha, numRows, IndexType (1));
328  }
329 }
330 
331 } // namespace Blas
332 } // namespace Details
333 } // namespace Tpetra
334 
335 #endif // TPETRA_DETAILS_FILL_HPP
Type traits for Tpetra&#39;s BLAS wrappers; an implementation detail of Tpetra::MultiVector.
void deep_copy(MultiVector< DS, DL, DG, DN > &dst, const MultiVector< SS, SL, SG, SN > &src)
Copy the contents of the MultiVector src into dst.
Implementation of ::Tpetra::Details::Blas::fill.
Do BLAS libraries (all that are compliant with the BLAS Standard) support the given &quot;scalar&quot; (matrix ...