44 #ifndef TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_LOCAL_DEEP_COPY_HPP 45 #define TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_LOCAL_DEEP_COPY_HPP 47 #include "Kokkos_Core.hpp" 48 #include "Kokkos_ArithTraits.hpp" 71 template<
class DstViewType,
73 class DstWhichVecsType,
74 class SrcWhichVecsType,
75 const bool DstConstStride,
76 const bool SrcConstStride>
77 struct LocalDeepCopyFunctor {
78 typedef typename DstViewType::execution_space execution_space;
79 typedef typename DstViewType::size_type index_type;
83 DstWhichVecsType dstWhichVecs_;
84 SrcWhichVecsType srcWhichVecs_;
85 const index_type numVecs_;
90 LocalDeepCopyFunctor (
const DstViewType& dst,
91 const SrcViewType& src,
92 const DstWhichVecsType& dstWhichVecs,
93 const SrcWhichVecsType& srcWhichVecs) :
96 dstWhichVecs_ (dstWhichVecs),
97 srcWhichVecs_ (srcWhichVecs),
98 numVecs_ (DstConstStride ? dst.dimension_1 () : dstWhichVecs.dimension_0 ())
100 TEUCHOS_TEST_FOR_EXCEPTION(
101 ! DstConstStride && ! SrcConstStride &&
102 dstWhichVecs.dimension_0 () != srcWhichVecs.dimension_0 (),
103 std::invalid_argument,
"LocalDeepCopyFunctor (4-arg constructor): " 104 "Neither src nor dst have constant stride, but " 105 "dstWhichVecs.dimension_0() = " << dstWhichVecs.dimension_0 ()
106 <<
" != srcWhichVecs.dimension_0() = " << srcWhichVecs.dimension_0 ()
108 TEUCHOS_TEST_FOR_EXCEPTION(
109 DstConstStride && ! SrcConstStride &&
110 srcWhichVecs.dimension_0 () != dst.dimension_1 (),
111 std::invalid_argument,
"LocalDeepCopyFunctor (4-arg constructor): " 112 "src does not have constant stride, but srcWhichVecs.dimension_0() = " 113 << srcWhichVecs.dimension_0 () <<
" != dst.dimension_1() = " 114 << dst.dimension_1 () <<
".");
115 TEUCHOS_TEST_FOR_EXCEPTION(
116 ! DstConstStride && SrcConstStride &&
117 dstWhichVecs.dimension_0 () != src.dimension_1 (),
118 std::invalid_argument,
"LocalDeepCopyFunctor (4-arg constructor): " 119 "dst does not have constant stride, but dstWhichVecs.dimension_0() = " 120 << dstWhichVecs.dimension_0 () <<
" != src.dimension_1() = " 121 << src.dimension_1 () <<
".");
126 LocalDeepCopyFunctor (
const DstViewType& dst,
const SrcViewType& src) :
129 numVecs_ (dst.dimension_1 ())
131 TEUCHOS_TEST_FOR_EXCEPTION(
132 ! DstConstStride || ! SrcConstStride, std::logic_error,
133 "Tpetra::LocalDeepCopyFunctor: You may not use the constant-stride " 134 "constructor if either of the Boolean template parameters is false.");
137 void KOKKOS_INLINE_FUNCTION operator () (
const index_type i)
const {
138 if (DstConstStride) {
139 if (SrcConstStride) {
140 for (index_type j = 0; j < numVecs_; ++j) {
141 dst_(i,j) = src_(i,j);
144 for (index_type j = 0; j < numVecs_; ++j) {
145 dst_(i,j) = src_(i,srcWhichVecs_(j));
149 if (SrcConstStride) {
150 for (index_type j = 0; j < numVecs_; ++j) {
151 dst_(i,dstWhichVecs_(j)) = src_(i,j);
154 for (index_type j = 0; j < numVecs_; ++j) {
155 dst_(i,dstWhichVecs_(j)) = src_(i,srcWhichVecs_(j));
162 template<
class DstViewType,
163 class SrcViewType = DstViewType,
164 class DstWhichVecsType = Kokkos::View<const typename DstViewType::size_type*, typename DstViewType::execution_space>,
165 class SrcWhichVecsType = Kokkos::View<const typename SrcViewType::size_type*, typename SrcViewType::execution_space> >
166 struct LocalDeepCopy {
167 typedef typename DstViewType::size_type index_type;
170 run (
const DstViewType& dst,
const SrcViewType& src,
171 const bool dstConstStride,
const bool srcConstStride)
173 TEUCHOS_TEST_FOR_EXCEPTION(
174 ! dstConstStride || ! srcConstStride, std::invalid_argument,
175 "LocalDeepCopy::run: You may only call the 4-argument version of this " 176 "function if dstConstStride and srcConstStride are both true.");
189 typedef LocalDeepCopyFunctor<DstViewType, SrcViewType,
190 DstWhichVecsType, SrcWhichVecsType,
true,
true> functor_type;
191 functor_type f (dst, src);
192 Kokkos::parallel_for (dst.dimension_0 (), f);
196 run (
const DstViewType& dst,
const SrcViewType& src,
197 const bool dstConstStride,
const bool srcConstStride,
198 const DstWhichVecsType& dstWhichVecs,
199 const SrcWhichVecsType& srcWhichVecs)
210 if (dstConstStride) {
211 if (srcConstStride) {
214 typedef LocalDeepCopyFunctor<DstViewType, SrcViewType,
215 DstWhichVecsType, SrcWhichVecsType,
true,
true> functor_type;
216 functor_type f (dst, src);
217 Kokkos::parallel_for (dst.dimension_0 (), f);
220 typedef LocalDeepCopyFunctor<DstViewType, SrcViewType,
221 DstWhichVecsType, SrcWhichVecsType,
true,
false> functor_type;
222 functor_type f (dst, src, srcWhichVecs, srcWhichVecs);
223 Kokkos::parallel_for (dst.dimension_0 (), f);
227 if (srcConstStride) {
228 typedef LocalDeepCopyFunctor<DstViewType, SrcViewType,
229 DstWhichVecsType, SrcWhichVecsType,
false,
true> functor_type;
230 functor_type f (dst, src, dstWhichVecs, dstWhichVecs);
231 Kokkos::parallel_for (dst.dimension_0 (), f);
234 typedef LocalDeepCopyFunctor<DstViewType, SrcViewType,
235 DstWhichVecsType, SrcWhichVecsType,
false,
false> functor_type;
236 functor_type f (dst, src, dstWhichVecs, srcWhichVecs);
237 Kokkos::parallel_for (dst.dimension_0 (), f);
244 template<
class DstViewType,
246 class DstWhichVecsType,
247 class SrcWhichVecsType>
249 localDeepCopy (
const DstViewType& dst,
const SrcViewType& src,
250 const bool dstConstStride,
const bool srcConstStride,
251 const DstWhichVecsType& dstWhichVecs,
252 const SrcWhichVecsType& srcWhichVecs)
254 typedef LocalDeepCopy<DstViewType, SrcViewType,
255 DstWhichVecsType, SrcWhichVecsType> impl_type;
256 impl_type::run (dst, src, dstConstStride, srcConstStride,
257 dstWhichVecs, srcWhichVecs);
260 template<
class DstViewType,
class SrcViewType>
262 localDeepCopyConstStride (
const DstViewType& dst,
const SrcViewType& src)
264 typedef LocalDeepCopy<DstViewType, SrcViewType> impl_type;
265 impl_type::run (dst, src,
true,
true);
271 #endif // TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_LOCAL_DEEP_COPY_HPP Namespace Tpetra contains the class and methods constituting the Tpetra library.
Implementation details of Tpetra.