42 #ifndef TPETRA_DETAILS_FILL_HPP
43 #define TPETRA_DETAILS_FILL_HPP
58 #include <type_traits>
67 memsetWrapper (
void* dest,
int ch, std::size_t count);
70 template<
class ViewType,
74 const int rank = ViewType::Rank>
78 fill (
const ExecutionSpace& execSpace,
80 const ValueType& alpha,
84 static_assert (std::is_integral<IndexType>::value,
85 "IndexType must be a built-in integer type.");
92 template<
class ViewType,
96 class Fill<ViewType, ValueType, ExecutionSpace, IndexType, 1> {
98 Fill (
const ViewType& X,
const ValueType& alpha) :
99 X_ (X), alpha_ (alpha)
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.");
107 KOKKOS_INLINE_FUNCTION
void
108 operator () (
const IndexType& i)
const
114 fill (
const ExecutionSpace& execSpace,
116 const ValueType& alpha,
117 const IndexType numRows,
120 typedef Kokkos::RangePolicy<ExecutionSpace, IndexType> range_type;
121 Kokkos::parallel_for (
"fill", range_type (0, numRows),
131 template<
class ViewType,
133 class ExecutionSpace,
135 class Fill<ViewType, ValueType, ExecutionSpace, IndexType, 2> {
137 Fill (
const ViewType& X,
138 const ValueType& alpha,
139 const IndexType numCols) :
140 X_ (X), alpha_ (alpha), numCols_ (numCols)
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.");
148 KOKKOS_INLINE_FUNCTION
void
149 operator () (
const IndexType& i)
const
151 for (IndexType j = 0; j < numCols_; ++j) {
157 fill (
const ExecutionSpace& execSpace,
159 const ValueType& alpha,
160 const IndexType numRows,
161 const IndexType numCols)
163 typedef Kokkos::RangePolicy<ExecutionSpace, IndexType> range_type;
164 Kokkos::parallel_for (
"fill", range_type (0, numRows),
165 Fill (X, alpha, numCols));
174 #if defined(KOKKOS_ENABLE_SERIAL)
175 template<
class ViewType,
180 struct Fill<ViewType,
187 fill (
const Kokkos::Serial& ,
189 const ValueType& alpha,
190 const IndexType numRows,
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;
204 if (podType && X.span_is_contiguous () && alpha == ValueType (0.0)) {
205 memsetWrapper (X.data (), 0, X.span () *
sizeof (view_value_type));
208 for (IndexType k = 0; k < numRows; ++k) {
214 #endif // defined(KOKKOS_ENABLE_SERIAL)
216 #if defined(KOKKOS_ENABLE_SERIAL)
217 template<
class ViewType,
222 struct Fill<ViewType,
229 fill (
const Kokkos::Serial& ,
231 const ValueType& alpha,
233 const IndexType numCols)
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;
245 constexpr
bool podType = BlasSupportsScalar<view_value_type>::value;
247 if (podType && alpha == ValueType (0.0)) {
248 if (X.span_is_contiguous ()) {
249 memsetWrapper (X.data (), 0, X.span () *
sizeof (view_value_type));
251 else if (std::is_same<array_layout, Kokkos::LayoutLeft>::value) {
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));
266 using execution_space =
typename ViewType::execution_space;
271 #endif // defined(KOKKOS_ENABLE_SERIAL)
286 template<
class ViewType,
289 class ExecutionSpace>
291 fill (
const ExecutionSpace& execSpace,
293 const ValueType& alpha,
294 const IndexType numRows,
295 const IndexType numCols)
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);
304 template<
class ViewType,
307 class ExecutionSpace>
309 fill (
const ExecutionSpace& execSpace,
311 const ValueType& alpha,
312 const IndexType numRows,
313 const IndexType numCols,
314 const size_t whichVectors[])
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));
335 #endif // TPETRA_DETAILS_FILL_HPP
Type traits for Tpetra'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 "scalar" (matrix ...