Tpetra parallel linear algebra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator 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  Kokkos::deep_copy (execSpace, X, alpha);
87  }
88 };
89 
91 template<class ViewType,
92  class ValueType,
93  class ExecutionSpace,
94  class IndexType>
95 class Fill<ViewType, ValueType, ExecutionSpace, IndexType, 1> {
96 public:
97  Fill (const ViewType& X, const ValueType& alpha) :
98  X_ (X), alpha_ (alpha)
99  {
100  static_assert (ViewType::Rank == 1,
101  "ViewType must be a rank-1 Kokkos::View.");
102  static_assert (std::is_integral<IndexType>::value,
103  "IndexType must be a built-in integer type.");
104  }
105 
106  KOKKOS_INLINE_FUNCTION void
107  operator () (const IndexType& i) const
108  {
109  X_[i] = alpha_;
110  }
111 
112  static void
113  fill (const ExecutionSpace& execSpace,
114  const ViewType& X,
115  const ValueType& alpha,
116  const IndexType numRows,
117  const IndexType /* numCols */)
118  {
119  typedef Kokkos::RangePolicy<ExecutionSpace, IndexType> range_type;
120  Kokkos::parallel_for ("fill", range_type (0, numRows),
121  Fill (X, alpha));
122  }
123 
124 private:
125  ViewType X_;
126  ValueType alpha_;
127 };
128 
130 template<class ViewType,
131  class ValueType,
132  class ExecutionSpace,
133  class IndexType>
134 class Fill<ViewType, ValueType, ExecutionSpace, IndexType, 2> {
135 public:
136  Fill (const ViewType& X,
137  const ValueType& alpha,
138  const IndexType numCols) :
139  X_ (X), alpha_ (alpha), numCols_ (numCols)
140  {
141  static_assert (ViewType::Rank == 2,
142  "ViewType must be a rank-2 Kokkos::View.");
143  static_assert (std::is_integral<IndexType>::value,
144  "IndexType must be a built-in integer type.");
145  }
146 
147  KOKKOS_INLINE_FUNCTION void
148  operator () (const IndexType& i) const
149  {
150  for (IndexType j = 0; j < numCols_; ++j) {
151  X_(i,j) = alpha_;
152  }
153  }
154 
155  static void
156  fill (const ExecutionSpace& execSpace,
157  const ViewType& X,
158  const ValueType& alpha,
159  const IndexType numRows,
160  const IndexType numCols)
161  {
162  typedef Kokkos::RangePolicy<ExecutionSpace, IndexType> range_type;
163  Kokkos::parallel_for ("fill", range_type (0, numRows),
164  Fill (X, alpha, numCols));
165  }
166 
167 private:
168  ViewType X_;
169  ValueType alpha_;
170  IndexType numCols_;
171 };
172 
173 #if defined(KOKKOS_ENABLE_SERIAL)
174 template<class ViewType,
177  class ValueType,
178  class IndexType>
179 struct Fill<ViewType,
180  ValueType,
181  Kokkos::Serial,
182  IndexType,
183  1>
184 {
185  static void
186  fill (const Kokkos::Serial& /* execSpace */,
187  const ViewType& X,
188  const ValueType& alpha,
189  const IndexType numRows,
190  const IndexType /* numCols */ )
191  {
192  static_assert (ViewType::Rank == 1,
193  "ViewType must be a rank-1 Kokkos::View.");
194  static_assert (std::is_integral<IndexType>::value,
195  "IndexType must be a built-in integer type.");
196  using ::Tpetra::Details::Blas::BlasSupportsScalar;
197  typedef typename ViewType::non_const_value_type view_value_type;
198 
199  // Do sizeof(view_value_type) and taking the address of a
200  // value_type instance work correctly with memset?
201  constexpr bool podType = BlasSupportsScalar<view_value_type>::value;
202 
203  if (podType && X.span_is_contiguous () && alpha == ValueType (0.0)) {
204  memsetWrapper (X.data (), 0, X.span () * sizeof (view_value_type));
205  }
206  else {
207  for (IndexType k = 0; k < numRows; ++k) {
208  X[k] = alpha;
209  }
210  }
211  }
212 };
213 #endif // defined(KOKKOS_ENABLE_SERIAL)
214 
215 #if defined(KOKKOS_ENABLE_SERIAL)
216 template<class ViewType,
219  class ValueType,
220  class IndexType>
221 struct Fill<ViewType,
222  ValueType,
223  Kokkos::Serial,
224  IndexType,
225  2>
226 {
227  static void
228  fill (const Kokkos::Serial& /* execSpace */,
229  const ViewType& X,
230  const ValueType& alpha,
231  const IndexType /* numRows */,
232  const IndexType numCols)
233  {
234  static_assert (ViewType::Rank == 2,
235  "ViewType must be a rank-2 Kokkos::View.");
236  static_assert (std::is_integral<IndexType>::value,
237  "IndexType must be a built-in integer type.");
238  using ::Tpetra::Details::Blas::BlasSupportsScalar;
239  typedef typename ViewType::non_const_value_type view_value_type;
240  typedef typename ViewType::array_layout array_layout;
241 
242  // Do sizeof(view_value_type) and taking the address of a
243  // value_type instance work correctly with memset?
244  constexpr bool podType = BlasSupportsScalar<view_value_type>::value;
245 
246  if (podType && alpha == ValueType (0.0)) {
247  if (X.span_is_contiguous ()) {
248  memsetWrapper (X.data (), 0, X.span () * sizeof (view_value_type));
249  }
250  else if (std::is_same<array_layout, Kokkos::LayoutLeft>::value) {
251  // Tpetra::MultiVector needs to optimize for LayoutLeft.
252  for (IndexType j = 0; j < numCols; ++j) {
253  auto X_j = Kokkos::subview (X, Kokkos::ALL (), j);
254  memsetWrapper (X_j.data (), 0,
255  X_j.extent (0) * sizeof (view_value_type));
256  }
257  }
258  else {
259  Kokkos::deep_copy (X, view_value_type (0.0));
260  }
261  }
262  else {
263  Kokkos::deep_copy (X, alpha);
264  }
265  }
266 };
267 #endif // defined(KOKKOS_ENABLE_SERIAL)
268 
269 } // namespace Impl
270 
271 //
272 // SKIP TO HERE FOR THE ACTUAL INTERFACE
273 //
274 
282 template<class ViewType,
283  class ValueType,
284  class IndexType,
285  class ExecutionSpace>
286 void
287 fill (const ExecutionSpace& execSpace,
288  const ViewType& X,
289  const ValueType& alpha,
290  const IndexType numRows,
291  const IndexType numCols)
292 {
293  static_assert (std::is_integral<IndexType>::value,
294  "IndexType must be a built-in integer type.");
295  typedef Impl::Fill<ViewType, ValueType, ExecutionSpace,
296  IndexType, ViewType::Rank> impl_type;
297  impl_type::fill (execSpace, X, alpha, numRows, numCols);
298 }
299 
300 template<class ViewType,
301  class ValueType,
302  class IndexType,
303  class ExecutionSpace>
304 void
305 fill (const ExecutionSpace& execSpace,
306  const ViewType& X,
307  const ValueType& alpha,
308  const IndexType numRows,
309  const IndexType numCols,
310  const size_t whichVectors[])
311 {
312  static_assert (ViewType::Rank == 2, "ViewType must be a rank-2 "
313  "Kokkos::View in order to call the \"whichVectors\" "
314  "specialization of fill.");
315  static_assert (std::is_integral<IndexType>::value,
316  "IndexType must be a built-in integer type.");
317  for (IndexType k = 0; k < numCols; ++k) {
318  const IndexType j = whichVectors[k];
319  auto X_j = Kokkos::subview (X, Kokkos::ALL (), j);
320  typedef decltype (X_j) one_d_view_type;
321  typedef Impl::Fill<one_d_view_type, ValueType, ExecutionSpace,
322  IndexType, 1> impl_type;
323  impl_type::fill (execSpace, X_j, alpha, numRows, IndexType (1));
324  }
325 }
326 
327 } // namespace Blas
328 } // namespace Details
329 } // namespace Tpetra
330 
331 #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 ...