Sacado Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
KokkosExp_View_Fad_Contiguous.hpp
Go to the documentation of this file.
1// @HEADER
2// ***********************************************************************
3//
4// Sacado Package
5// Copyright (2006) 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// This library is free software; you can redistribute it and/or modify
11// it under the terms of the GNU Lesser General Public License as
12// published by the Free Software Foundation; either version 2.1 of the
13// License, or (at your option) any later version.
14//
15// This library is distributed in the hope that it will be useful, but
16// WITHOUT ANY WARRANTY; without even the implied warranty of
17// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18// Lesser General Public License for more details.
19//
20// You should have received a copy of the GNU Lesser General Public
21// License along with this library; if not, write to the Free Software
22// Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301
23// USA
24// Questions? Contact David M. Gay (dmgay@sandia.gov) or Eric T. Phipps
25// (etphipp@sandia.gov).
26//
27// ***********************************************************************
28// @HEADER
29
30#ifndef KOKKOS_EXPERIMENTAL_VIEW_SACADO_FAD_CONTIGUOUS_HPP
31#define KOKKOS_EXPERIMENTAL_VIEW_SACADO_FAD_CONTIGUOUS_HPP
32
33#include "Sacado_ConfigDefs.h"
34#if defined(HAVE_SACADO_KOKKOSCORE)
35
37
38// Some default traits definitions that need to be defined even if the view
39// specialization is disabled
40namespace Kokkos {
41
42template <typename ViewType, typename Enabled = void>
43struct ThreadLocalScalarType {
44 typedef typename ViewType::non_const_value_type type;
45};
46
47template <typename ViewType>
48struct ViewScalarStride {
49 static const unsigned stride =
50 Impl::LayoutScalarStride< typename ViewType::array_layout>::stride;
51 static const bool is_unit_stride =
52 Impl::LayoutScalarStride< typename ViewType::array_layout>::is_unit_stride;
53};
54
55} // namespace Kokkos
56
57namespace Sacado {
58
59 namespace Fad {
60
61 /* Define a partition of a View of Sacado Fad type */
62 template <unsigned Size = 0>
63 struct Partition {
64 static const unsigned PartitionSize = Size;
65 unsigned offset ;
66 unsigned stride ;
67
68 template< typename iType0 , typename iType1 >
69 KOKKOS_INLINE_FUNCTION
70 Partition( const iType0 & i0 , const iType1 & i1 ) :
71 offset(i0), stride(i1) {
72 }
73 };
74
75 template <typename T>
76 struct is_fad_partition {
77 static const bool value = false;
78 };
79
80 template <unsigned Stride>
81 struct is_fad_partition< Partition<Stride> > {
82 static const bool value = true;
83 };
84
85 }
86
87 // Type of local scalar type when partitioning a view
88 template <typename T, unsigned Stride = 0>
89 struct LocalScalarType {
90 typedef T type;
91 };
92 template <typename T, unsigned Stride>
93 struct LocalScalarType<const T, Stride> {
94 typedef typename LocalScalarType<T,Stride>::type lst;
95 typedef const lst type;
96 };
97
98 // For DFad, the size is not part of the type, so the default implementation
99 // is sufficient
100
101 // Type of local scalar type when partitioning a view
102 //
103 // For SLFad, divde the array size by the given stride
104 namespace Fad {
105 namespace Exp {
106 template <typename T, int N> class StaticStorage;
107 template <typename S> class GeneralFad;
108 }
109 }
110 template <typename T, int N, unsigned Stride>
111 struct LocalScalarType< Fad::Exp::GeneralFad< Fad::Exp::StaticStorage<T,N> >,
112 Stride > {
113 static const int Ns = (N+Stride-1) / Stride;
114 typedef Fad::Exp::GeneralFad< Fad::Exp::StaticStorage<T,Ns> > type;
115 };
116#ifndef SACADO_NEW_FAD_DESIGN_IS_DEFAULT
117 namespace Fad {
118 template <typename T, int N> class SLFad;
119 }
120 template <typename T, int N, unsigned Stride>
121 struct LocalScalarType< Fad::SLFad<T,N>, Stride > {
122 static const int Ns = (N+Stride-1) / Stride;
123 typedef Fad::SLFad<T,Ns> type;
124 };
125#endif
126
127 // Type of local scalar type when partitioning a view
128 //
129 // For SFad, divde the array size by the given stride. If it divides evenly,
130 // use SFad, otherwise use SLFad
131 namespace Fad {
132 namespace Exp {
133 template <typename T, int N> class StaticFixedStorage;
134 template <typename T, int N> class StaticStorage;
135 template <typename S> class GeneralFad;
136 }
137 }
138 template <typename T, int N, unsigned Stride>
139 struct LocalScalarType< Fad::Exp::GeneralFad< Fad::Exp::StaticFixedStorage<T,N> >,
140 Stride > {
141 static const int Ns = (N+Stride-1) / Stride;
142 typedef typename std::conditional<
143 Ns == N/Stride ,
144 Fad::Exp::GeneralFad< Fad::Exp::StaticFixedStorage<T,Ns> > ,
145 Fad::Exp::GeneralFad< Fad::Exp::StaticStorage<T,Ns> >
146 >::type type;
147 };
148
149#ifndef SACADO_NEW_FAD_DESIGN_IS_DEFAULT
150 namespace Fad {
151 template <typename T, int N> class SFad;
152 }
153 template <typename T, int N, unsigned Stride>
154 struct LocalScalarType< Fad::SFad<T,N>, Stride > {
155 static const int Ns = (N+Stride-1) / Stride;
156 typedef typename std::conditional< Ns == N/Stride , Fad::SFad<T,Ns> , Fad::SLFad<T,Ns> >::type type;
157 };
158#endif
159
160 template <unsigned Stride, typename T>
161 KOKKOS_INLINE_FUNCTION
162 const T& partition_scalar(const T& x) { return x; }
163
164} // namespace Sacado
165
166// Make sure the user really wants these View specializations
167#if defined(HAVE_SACADO_VIEW_SPEC) && !defined(SACADO_DISABLE_FAD_VIEW_SPEC)
168
169#include "Sacado_Traits.hpp"
170#include "Kokkos_Core.hpp"
171#include "impl/Kokkos_ViewMapping.hpp"
172
173//----------------------------------------------------------------------------
174
175namespace Sacado {
176
177#if defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__)
178 namespace Fad {
179 namespace Exp {
180 template <typename T, typename U> class DynamicStorage;
181 template <typename T, int N> class StaticFixedStorage;
182 template <typename T, int N> class StaticStorage;
183 template <typename S> class GeneralFad;
184 }
185 }
186#ifndef SACADO_VIEW_CUDA_HIERARCHICAL_DFAD
187 template <unsigned Stride, typename T, typename U>
188 KOKKOS_INLINE_FUNCTION
189 typename LocalScalarType< Fad::Exp::GeneralFad< Fad::Exp::DynamicStorage<T,U> >, Stride >::type
190 partition_scalar(const Fad::Exp::GeneralFad< Fad::Exp::DynamicStorage<T,U> >& x) {
191 typedef typename LocalScalarType< Fad::Exp::GeneralFad< Fad::Exp::DynamicStorage<T,U> >, Stride >::type ret_type;
192 const int size = (x.size()+blockDim.x-threadIdx.x-1) / blockDim.x;
193 const int offset = threadIdx.x;
194 ret_type xp(size, x.val());
195
196 // Note: we can't use x.dx(offset+i*Stride) if
197 // SACADO_VIEW_CUDA_HIERARCHICAL_DFAD_STRIDED is defined because it already
198 // uses blockDim.x in its index calculation. This approach should work
199 // regardless
200 const T* dx = x.dx();
201 for (int i=0; i<size; ++i)
202 xp.fastAccessDx(i) = dx[offset+i*Stride];
203
204 return xp;
205 }
206#endif
207 template <unsigned Stride, typename T, int N>
208 KOKKOS_INLINE_FUNCTION
209 typename LocalScalarType< Fad::Exp::GeneralFad< Fad::Exp::StaticStorage<T,N> >, Stride >::type
210 partition_scalar(const Fad::Exp::GeneralFad< Fad::Exp::StaticStorage<T,N> >& x) {
211 typedef typename LocalScalarType< Fad::Exp::GeneralFad< Fad::Exp::StaticStorage<T,N> >, Stride >::type ret_type;
212 const int size = (x.size()+blockDim.x-threadIdx.x-1) / blockDim.x;
213 const int offset = threadIdx.x;
214 ret_type xp(size, x.val());
215 for (int i=0; i<size; ++i)
216 xp.fastAccessDx(i) = x.fastAccessDx(offset+i*Stride);
217 return xp;
218 }
219 template <unsigned Stride, typename T, int N>
220 KOKKOS_INLINE_FUNCTION
221 typename LocalScalarType< Fad::Exp::GeneralFad< Fad::Exp::StaticFixedStorage<T,N> >, Stride >::type
222 partition_scalar(const Fad::Exp::GeneralFad< Fad::Exp::StaticFixedStorage<T,N> >& x) {
223 typedef typename LocalScalarType< Fad::Exp::GeneralFad< Fad::Exp::StaticFixedStorage<T,N> >, Stride >::type ret_type;
224 const int size = (x.size()+blockDim.x-threadIdx.x-1) / blockDim.x;
225 const int offset = threadIdx.x;
226 ret_type xp(size, x.val());
227 for (int i=0; i<size; ++i)
228 xp.fastAccessDx(i) = x.fastAccessDx(offset+i*Stride);
229 return xp;
230 }
231
232#ifndef SACADO_NEW_FAD_DESIGN_IS_DEFAULT
233 namespace Fad {
234 template <typename T> class DFad;
235 template <typename T, int N> class SLFad;
236 template <typename T, int N> class SFad;
237 }
238#ifndef SACADO_VIEW_CUDA_HIERARCHICAL_DFAD
239 template <unsigned Stride, typename T>
240 KOKKOS_INLINE_FUNCTION
241 typename LocalScalarType< Fad::DFad<T>, Stride >::type
242 partition_scalar(const Fad::DFad<T>& x) {
243 typedef typename LocalScalarType< Fad::DFad<T>, Stride >::type ret_type;
244 const int size = (x.size()+blockDim.x-threadIdx.x-1) / blockDim.x;
245 const int offset = threadIdx.x;
246 ret_type xp(size, x.val());
247
248 // Note: we can't use x.dx(offset+i*Stride) if
249 // SACADO_VIEW_CUDA_HIERARCHICAL_DFAD_STRIDED is defined because it already
250 // uses blockDim.x in its index calculation. This approach should work
251 // regardless
252 const T* dx = x.dx();
253 for (int i=0; i<size; ++i)
254 xp.fastAccessDx(i) = dx[offset+i*Stride];
255
256 return xp;
257 }
258#endif
259 template <unsigned Stride, typename T, int N>
260 KOKKOS_INLINE_FUNCTION
261 typename LocalScalarType< Fad::SLFad<T,N>, Stride >::type
262 partition_scalar(const Fad::SLFad<T,N>& x) {
263 typedef typename LocalScalarType< Fad::SLFad<T,N>, Stride >::type ret_type;
264 const int size = (x.size()+blockDim.x-threadIdx.x-1) / blockDim.x;
265 const int offset = threadIdx.x;
266 ret_type xp(size, x.val());
267 for (int i=0; i<size; ++i)
268 xp.fastAccessDx(i) = x.fastAccessDx(offset+i*Stride);
269 return xp;
270 }
271 template <unsigned Stride, typename T, int N>
272 KOKKOS_INLINE_FUNCTION
273 typename LocalScalarType< Fad::SFad<T,N>, Stride >::type
274 partition_scalar(const Fad::SFad<T,N>& x) {
275 typedef typename LocalScalarType< Fad::SFad<T,N>, Stride >::type ret_type;
276 const int size = (x.size()+blockDim.x-threadIdx.x-1) / blockDim.x;
277 const int offset = threadIdx.x;
278 ret_type xp(size, x.val());
279 for (int i=0; i<size; ++i)
280 xp.fastAccessDx(i) = x.fastAccessDx(offset+i*Stride);
281 return xp;
282 }
283#endif // SACADO_NEW_FAD_DESIGN_IS_DEFAULT
284
285#endif // defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__)
286
287} // namespace Sacado
288
289//----------------------------------------------------------------------------
290
291namespace Kokkos {
292
293template< unsigned Stride, typename D, typename ... P >
294KOKKOS_INLINE_FUNCTION
295typename Kokkos::Impl::ViewMapping< void, typename Kokkos::ViewTraits<D,P...>, Sacado::Fad::Partition<Stride> >::type
296partition( const Kokkos::View<D,P...> & src ,
297 const unsigned offset ,
298 const unsigned stride )
299{
300 typedef Kokkos::ViewTraits<D,P...> traits;
301 typedef typename Kokkos::Impl::ViewMapping< void, traits, Sacado::Fad::Partition<Stride> >::type DstViewType;
302 const Sacado::Fad::Partition<Stride> part( offset , stride );
303 return DstViewType(src, part);
304}
305
306template <typename ViewType>
307struct ThreadLocalScalarType<
308 ViewType,
309 typename std::enable_if< is_view_fad_contiguous<ViewType>::value >::type > {
310 typedef typename ViewType::traits TraitsType;
311 typedef Impl::ViewMapping<TraitsType, typename TraitsType::specialize> MappingType;
312 typedef typename MappingType::thread_local_scalar_type type;
313};
314
315namespace Impl {
316
317#if defined (KOKKOS_ENABLE_CUDA) && defined(SACADO_VIEW_CUDA_HIERARCHICAL)
318template< class OutputView >
319struct SacadoViewFill<
320 OutputView,
321 typename std::enable_if<
322 ( Kokkos::is_view_fad_contiguous<OutputView>::value &&
323 std::is_same<typename OutputView::execution_space, Kokkos::Cuda>::value &&
324 !Kokkos::ViewScalarStride<OutputView>::is_unit_stride )
325 >::type
326 >
327{
328 typedef typename OutputView::const_value_type const_value_type ;
329 typedef typename OutputView::execution_space execution_space ;
330 typedef Kokkos::TeamPolicy< execution_space> team_policy;
331 typedef typename team_policy::member_type team_impl_handle;
332 typedef typename Kokkos::ThreadLocalScalarType<OutputView>::type local_scalar_type;
333 static const unsigned stride = Kokkos::ViewScalarStride<OutputView>::stride;
334
335 const OutputView output ;
336 const_value_type input ;
337
338 KOKKOS_INLINE_FUNCTION
339 void operator()( const size_t i0 ) const
340 {
341 local_scalar_type input_stride = Sacado::partition_scalar<stride>(input);
342
343 const size_t n1 = output.extent(1);
344 const size_t n2 = output.extent(2);
345 const size_t n3 = output.extent(3);
346 const size_t n4 = output.extent(4);
347 const size_t n5 = output.extent(5);
348 const size_t n6 = output.extent(6);
349
350 for ( size_t i1 = 0 ; i1 < n1 ; ++i1 ) {
351 for ( size_t i2 = 0 ; i2 < n2 ; ++i2 ) {
352 for ( size_t i3 = 0 ; i3 < n3 ; ++i3 ) {
353 for ( size_t i4 = 0 ; i4 < n4 ; ++i4 ) {
354 for ( size_t i5 = 0 ; i5 < n5 ; ++i5 ) {
355 for ( size_t i6 = 0 ; i6 < n6 ; ++i6 ) {
356 output.access(i0,i1,i2,i3,i4,i5,i6) = input_stride ;
357 }}}}}}
358 }
359
360 KOKKOS_INLINE_FUNCTION
361 void operator()( const team_impl_handle& team ) const
362 {
363 const size_t i0 = team.league_rank()*team.team_size() + team.team_rank();
364 if (i0 < output.extent(0))
365 (*this)(i0);
366 }
367
368 SacadoViewFill( const OutputView & arg_out , const_value_type & arg_in )
369 : output( arg_out ), input( arg_in )
370 {
371 const size_t team_size = 256 / stride;
372 team_policy policy( (output.extent(0)+team_size-1)/team_size ,
373 team_size , stride );
374 Kokkos::parallel_for( policy, *this );
375 }
376};
377#endif
378
379#if defined (KOKKOS_ENABLE_HIP) && defined(SACADO_VIEW_CUDA_HIERARCHICAL)
380template< class OutputView >
381struct SacadoViewFill<
382 OutputView,
383 typename std::enable_if<
384 ( Kokkos::is_view_fad_contiguous<OutputView>::value &&
385 std::is_same<typename OutputView::execution_space, Kokkos::Experimental::HIP>::value &&
386 !Kokkos::ViewScalarStride<OutputView>::is_unit_stride )
387 >::type
388 >
389{
390 typedef typename OutputView::const_value_type const_value_type ;
391 typedef typename OutputView::execution_space execution_space ;
392 typedef Kokkos::TeamPolicy< execution_space> team_policy;
393 typedef typename team_policy::member_type team_impl_handle;
394 typedef typename Kokkos::ThreadLocalScalarType<OutputView>::type local_scalar_type;
395 static const unsigned stride = Kokkos::ViewScalarStride<OutputView>::stride;
396
397 const OutputView output ;
398 const_value_type input ;
399
400 KOKKOS_INLINE_FUNCTION
401 void operator()( const size_t i0 ) const
402 {
403 local_scalar_type input_stride = Sacado::partition_scalar<stride>(input);
404
405 const size_t n1 = output.extent(1);
406 const size_t n2 = output.extent(2);
407 const size_t n3 = output.extent(3);
408 const size_t n4 = output.extent(4);
409 const size_t n5 = output.extent(5);
410 const size_t n6 = output.extent(6);
411 const size_t n7 = output.extent(7);
412
413 for ( size_t i1 = 0 ; i1 < n1 ; ++i1 ) {
414 for ( size_t i2 = 0 ; i2 < n2 ; ++i2 ) {
415 for ( size_t i3 = 0 ; i3 < n3 ; ++i3 ) {
416 for ( size_t i4 = 0 ; i4 < n4 ; ++i4 ) {
417 for ( size_t i5 = 0 ; i5 < n5 ; ++i5 ) {
418 for ( size_t i6 = 0 ; i6 < n6 ; ++i6 ) {
419 output.access(i0,i1,i2,i3,i4,i5,i6) = input_stride ;
420 }}}}}}
421 }
422
423 KOKKOS_INLINE_FUNCTION
424 void operator()( const team_impl_handle& team ) const
425 {
426 const size_t i0 = team.league_rank()*team.team_size() + team.team_rank();
427 if (i0 < output.extent(0))
428 (*this)(i0);
429 }
430
431 SacadoViewFill( const OutputView & arg_out , const_value_type & arg_in )
432 : output( arg_out ), input( arg_in )
433 {
434 const size_t team_size = 256 / stride;
435 team_policy policy( (output.extent(0)+team_size-1)/team_size ,
436 team_size , stride );
437 Kokkos::parallel_for( policy, *this );
438 }
439};
440#endif
441
442} // namespace Impl
443} // namespace Kokkos
444
445//----------------------------------------------------------------------------
446
447namespace Kokkos {
448namespace Impl {
449
450template< class ... Args >
451struct is_ViewSpecializeSacadoFadContiguous { enum { value = false }; };
452
453template< class D , class ... P , class ... Args >
454struct is_ViewSpecializeSacadoFadContiguous< Kokkos::View<D,P...> , Args... > {
455 enum { value =
456 std::is_same< typename Kokkos::ViewTraits<D,P...>::specialize
457 , ViewSpecializeSacadoFadContiguous >::value
458 &&
459 ( ( sizeof...(Args) == 0 ) ||
460 is_ViewSpecializeSacadoFadContiguous< Args... >::value ) };
461};
462
463} // namespace Impl
464} // namespace Kokkos
465
466//----------------------------------------------------------------------------
467
468namespace Kokkos {
469namespace Impl {
470
471// Compute a partitioned fad size given a stride. Return 0 if the stride
472// does not evenly divide the size
473KOKKOS_INLINE_FUNCTION
474constexpr unsigned computeFadPartitionSize(unsigned size, unsigned stride)
475{
476 return
477 ((size+stride-1)/stride) == (size/stride) ? ((size+stride-1)/stride) : 0;
478}
479
480// Create new Layout with Fad dimension set to the last
481// Note: we use enable_if<> here to handle both LayoutLeft and
482// LayoutContiguous<LayoutLeft>
483template <unsigned rank, unsigned static_dim, typename Layout>
484KOKKOS_INLINE_FUNCTION
485typename std::enable_if< !std::is_same<Layout, LayoutLeft>::value &&
486 !std::is_same<Layout, LayoutStride>::value,
487 Layout>::type
488create_fad_array_layout(const Layout& layout)
489{
490 size_t dims[8];
491 for (int i=0; i<8; ++i)
492 dims[i] = layout.dimension[i];
493 if (static_dim > 0)
494 dims[rank] = static_dim+1;
495 return Layout( dims[0], dims[1], dims[2], dims[3],
496 dims[4], dims[5], dims[6], dims[7] );
497}
498
499// Create new Layout with Fad dimension set to the last
500// Note: we use enable_if<> here to handle both LayoutStride and
501// LayoutContiguous<LayoutStride>
502template <unsigned rank, unsigned static_dim, typename Layout>
503KOKKOS_INLINE_FUNCTION
504typename std::enable_if< std::is_same<Layout, LayoutStride>::value, Layout>::type
505create_fad_array_layout(const Layout& layout)
506{
507 size_t dims[8], strides[8];
508 for (int i=0; i<8; ++i) {
509 dims[i] = layout.dimension[i];
510 strides[i] = layout.stride[i];
511 }
512 if (static_dim > 0) {
513 dims[rank] = static_dim+1;
514 strides[rank] = 1;
515 }
516 return Layout( dims[0], strides[0],
517 dims[1], strides[1],
518 dims[2], strides[2],
519 dims[3], strides[3],
520 dims[4], strides[4],
521 dims[5], strides[5],
522 dims[6], strides[6],
523 dims[7], strides[7] );
524}
525
526// Create new LayoutLeft with Fad dimension shuffled to the first
527// Note: we use enable_if<> here to handle both LayoutLeft and
528// LayoutContiguous<LayoutLeft>
529 template <unsigned rank, unsigned static_dim, typename Layout>
530KOKKOS_INLINE_FUNCTION
531typename std::enable_if< std::is_same<Layout, LayoutLeft>::value, Layout >::type
532create_fad_array_layout(const Layout& layout)
533{
534 size_t dims[8];
535 for (int i=0; i<8; ++i)
536 dims[i] = layout.dimension[i];
537 size_t fad_dim = static_dim == 0 ? dims[rank] : static_dim+1;
538 for (int i=rank; i>=1; --i)
539 dims[i] = dims[i-1];
540 dims[0] = fad_dim;
541 return Layout( dims[0], dims[1], dims[2], dims[3],
542 dims[4], dims[5], dims[6], dims[7] );
543}
544
545template <unsigned Rank, typename Dimension, typename Layout>
546KOKKOS_INLINE_FUNCTION
547typename std::enable_if< !std::is_same<Layout, LayoutLeft>::value, size_t>::type
548getFadDimension(const ViewOffset<Dimension,Layout,void>& offset)
549{
550 return
551 ( Rank == 0 ? offset.dimension_0() :
552 ( Rank == 1 ? offset.dimension_1() :
553 ( Rank == 2 ? offset.dimension_2() :
554 ( Rank == 3 ? offset.dimension_3() :
555 ( Rank == 4 ? offset.dimension_4() :
556 ( Rank == 5 ? offset.dimension_5() :
557 ( Rank == 6 ? offset.dimension_6() :
558 offset.dimension_7() )))))));
559}
560
561template <unsigned Rank, typename Dimension, typename Layout>
562KOKKOS_INLINE_FUNCTION
563typename std::enable_if< std::is_same<Layout, LayoutLeft>::value, size_t >::type
564getFadDimension(const ViewOffset<Dimension,Layout,void>& offset)
565{
566 return offset.dimension_0();
567}
568
569template< class Traits >
570class ViewMapping< Traits , /* View internal mapping */
571 typename std::enable_if<
572 ( std::is_same< typename Traits::specialize
573 , ViewSpecializeSacadoFadContiguous >::value
574 &&
575 ( std::is_same< typename Traits::array_layout
576 , Kokkos::LayoutLeft >::value
577 ||
578 std::is_same< typename Traits::array_layout
579 , Kokkos::LayoutRight >::value
580 ||
581 std::is_same< typename Traits::array_layout
582 , Kokkos::LayoutStride >::value
583 )
584 )
585 , typename Traits::specialize
586 >::type >
587{
588private:
589
590 template< class , class ... > friend class ViewMapping ;
591 template< class , class ... > friend class Kokkos::View ;
592
593 typedef typename Traits::value_type fad_type ;
594 typedef typename Sacado::ValueType< fad_type >::type fad_value_type ;
595 typedef typename
596 std::add_const< fad_value_type >::type const_fad_value_type ;
597
598public:
599 enum { is_assignable_data_type = true };
600
601 enum { FadStaticDimension = Sacado::StaticSize< fad_type >::value };
602 enum { PartitionedFadStride = Traits::array_layout::scalar_stride };
603
604 // The partitioned static size -- this will be 0 if ParitionedFadStride
605 // does not evenly divide FadStaticDimension
606 enum { PartitionedFadStaticDimension =
607 computeFadPartitionSize(FadStaticDimension,PartitionedFadStride) };
608
609#ifdef KOKKOS_ENABLE_CUDA
610 typedef typename Sacado::LocalScalarType< fad_type, unsigned(PartitionedFadStride) >::type strided_scalar_type;
611 typedef typename std::conditional< std::is_same<typename Traits::execution_space, Kokkos::Cuda>::value, strided_scalar_type, fad_type >::type thread_local_scalar_type;
612#elif defined(KOKKOS_ENABLE_HIP)
613 typedef typename Sacado::LocalScalarType< fad_type, unsigned(PartitionedFadStride) >::type strided_scalar_type;
614 typedef typename std::conditional< std::is_same<typename Traits::execution_space, Kokkos::Experimental::HIP>::value, strided_scalar_type, fad_type >::type thread_local_scalar_type;
615#else
616 typedef fad_type thread_local_scalar_type;
617#endif
618
619private:
621
622 typedef fad_value_type * handle_type ;
623
624 typedef ViewArrayAnalysis< typename Traits::data_type > array_analysis ;
625
626 typedef ViewOffset< typename Traits::dimension
627 , typename Traits::array_layout
628 , void
629 > offset_type ;
630
631 // Prepend/append the fad dimension for the internal offset mapping.
632 static constexpr bool is_layout_left =
633 std::is_same< typename Traits::array_layout, Kokkos::LayoutLeft>::value;
634 typedef ViewOffset<
635 typename std::conditional<
636 is_layout_left,
637 typename array_analysis::dimension::
638 template prepend<0>::type,
639 typename array_analysis::dimension::
640 template append<0>::type >::type,
641 typename Traits::array_layout,
642 void >
643 array_offset_type ;
644
645 handle_type m_impl_handle ;
646 offset_type m_impl_offset ;
647 array_offset_type m_array_offset ;
648 sacado_size_type m_fad_size ;
649
650 // These are for manual partitioning, and will likely be removed
651 unsigned m_original_fad_size ;
652 unsigned m_fad_stride ;
653 unsigned m_fad_index ;
654
655public:
656
657 //----------------------------------------
658 // Domain dimensions
659
660 enum { Rank = Traits::dimension::rank };
661
662 // Rank corresponding to the sacado dimension
663 enum { Sacado_Rank = std::is_same< typename Traits::array_layout, Kokkos::LayoutLeft >::value ? 0 : Rank+1 };
664
665 // Using the internal offset mapping so limit to public rank:
666 template< typename iType >
667 KOKKOS_INLINE_FUNCTION constexpr size_t extent( const iType & r ) const
668 { return m_impl_offset.m_dim.extent(r); }
669
670 KOKKOS_INLINE_FUNCTION constexpr
671 typename Traits::array_layout layout() const
672 { return m_impl_offset.layout(); }
673
674 KOKKOS_INLINE_FUNCTION constexpr size_t dimension_0() const
675 { return m_impl_offset.dimension_0(); }
676 KOKKOS_INLINE_FUNCTION constexpr size_t dimension_1() const
677 { return m_impl_offset.dimension_1(); }
678 KOKKOS_INLINE_FUNCTION constexpr size_t dimension_2() const
679 { return m_impl_offset.dimension_2(); }
680 KOKKOS_INLINE_FUNCTION constexpr size_t dimension_3() const
681 { return m_impl_offset.dimension_3(); }
682 KOKKOS_INLINE_FUNCTION constexpr size_t dimension_4() const
683 { return m_impl_offset.dimension_4(); }
684 KOKKOS_INLINE_FUNCTION constexpr size_t dimension_5() const
685 { return m_impl_offset.dimension_5(); }
686 KOKKOS_INLINE_FUNCTION constexpr size_t dimension_6() const
687 { return m_impl_offset.dimension_6(); }
688 KOKKOS_INLINE_FUNCTION constexpr size_t dimension_7() const
689 { return m_impl_offset.dimension_7(); }
690
691 // Is a regular layout with uniform striding for each index.
692 // Since we allow for striding within the data type, we can't guarantee
693 // regular striding
694 using is_regular = std::false_type ;
695
696 // FIXME: Adjust these for m_stride
697 KOKKOS_INLINE_FUNCTION constexpr size_t stride_0() const
698 { return m_impl_offset.stride_0(); }
699 KOKKOS_INLINE_FUNCTION constexpr size_t stride_1() const
700 { return m_impl_offset.stride_1(); }
701 KOKKOS_INLINE_FUNCTION constexpr size_t stride_2() const
702 { return m_impl_offset.stride_2(); }
703 KOKKOS_INLINE_FUNCTION constexpr size_t stride_3() const
704 { return m_impl_offset.stride_3(); }
705 KOKKOS_INLINE_FUNCTION constexpr size_t stride_4() const
706 { return m_impl_offset.stride_4(); }
707 KOKKOS_INLINE_FUNCTION constexpr size_t stride_5() const
708 { return m_impl_offset.stride_5(); }
709 KOKKOS_INLINE_FUNCTION constexpr size_t stride_6() const
710 { return m_impl_offset.stride_6(); }
711 KOKKOS_INLINE_FUNCTION constexpr size_t stride_7() const
712 { return m_impl_offset.stride_7(); }
713
714 template< typename iType >
715 KOKKOS_INLINE_FUNCTION void stride( iType * const s ) const
716 { m_impl_offset.stride(s); }
717
718 // Size of sacado scalar dimension
719 KOKKOS_FORCEINLINE_FUNCTION constexpr unsigned dimension_scalar() const
720#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && ( defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
721 { return PartitionedFadStaticDimension ? PartitionedFadStaticDimension+1 : (m_fad_size.value+blockDim.x-threadIdx.x-1) / blockDim.x + 1; }
722#else
723 { return m_fad_size.value+1; }
724#endif
725
726 // trode of sacado scalar dimension
727 KOKKOS_FORCEINLINE_FUNCTION constexpr unsigned stride_scalar() const
728 { return m_fad_stride; }
729
730 //----------------------------------------
731 // Range of mapping
732
733#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && ( defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
734 // Return type of reference operators
735 // this only works if you are using a team-parallel operation on Cuda or HIP!
736 // typedef typename
737 // Sacado::ViewFadType< thread_local_scalar_type , PartitionedFadStaticDimension , (unsigned(ParitionedFadStride) > 1 ? PartitionedFadStride : 0) >::type reference_type ;
738 typedef typename
740#else
741 // Return type of reference operators
742 typedef typename
744#endif
745
747 typedef fad_value_type * pointer_type ;
748
750 KOKKOS_INLINE_FUNCTION constexpr size_t span() const
751 { return m_array_offset.span(); }
752
754 KOKKOS_INLINE_FUNCTION constexpr bool span_is_contiguous() const
755 { return m_array_offset.span_is_contiguous() && (m_fad_stride == 1); }
756
758 KOKKOS_INLINE_FUNCTION constexpr pointer_type data() const
759#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && (defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
760 { return m_impl_handle + threadIdx.x; }
761#else
762 { return m_impl_handle + m_fad_index; }
763#endif
764
765 //----------------------------------------
766
767 KOKKOS_FORCEINLINE_FUNCTION
768 reference_type
769 reference() const
770 {
771#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && ( defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
772 const unsigned index = threadIdx.x;
773 const unsigned strd = blockDim.x;
774 const unsigned size = (m_fad_size.value+blockDim.x-threadIdx.x-1) / blockDim.x;
775#else
776 const unsigned index = m_fad_index;
777 const unsigned strd = m_fad_stride;
778 const unsigned size = m_fad_size.value;
779#endif
780 return reference_type( m_impl_handle + index
781 , m_impl_handle + m_original_fad_size
782 , size
783 , strd ); }
784
785 template< typename I0 >
786 KOKKOS_FORCEINLINE_FUNCTION
787 typename std::enable_if< Kokkos::Impl::are_integral<I0>::value &&
788 is_layout_left, reference_type>::type
789 reference( const I0 & i0 ) const
790 { pointer_type beg = m_impl_handle + m_array_offset(0,i0);
791#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && ( defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
792 const unsigned index = threadIdx.x;
793 const unsigned strd = blockDim.x;
794 const unsigned size = (m_fad_size.value+blockDim.x-threadIdx.x-1) / blockDim.x;
795#else
796 const unsigned index = m_fad_index;
797 const unsigned strd = m_fad_stride;
798 const unsigned size = m_fad_size.value;
799#endif
800 return reference_type( beg + index
801 , beg + m_original_fad_size
802 , size
803 , strd ); }
804
805 template< typename I0 >
806 KOKKOS_FORCEINLINE_FUNCTION
807 typename std::enable_if< Kokkos::Impl::are_integral<I0>::value &&
808 !is_layout_left, reference_type>::type
809 reference( const I0 & i0 ) const
810 { pointer_type beg = m_impl_handle + m_array_offset(i0,0);
811#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && ( defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
812 const unsigned index = threadIdx.x;
813 const unsigned strd = blockDim.x;
814 const unsigned size = (m_fad_size.value+blockDim.x-threadIdx.x-1) / blockDim.x;
815#else
816 const unsigned index = m_fad_index;
817 const unsigned strd = m_fad_stride;
818 const unsigned size = m_fad_size.value;
819#endif
820 return reference_type( beg + index
821 , beg + m_original_fad_size
822 , size
823 , strd ); }
824
825 template< typename I0 , typename I1 >
826 KOKKOS_FORCEINLINE_FUNCTION
827 typename std::enable_if< Kokkos::Impl::are_integral<I0,I1>::value &&
828 is_layout_left, reference_type>::type
829 reference( const I0 & i0 , const I1 & i1 ) const
830 { pointer_type beg = m_impl_handle + m_array_offset(0,i0,i1);
831#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && (defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
832 const unsigned index = threadIdx.x;
833 const unsigned strd = blockDim.x;
834 const unsigned size = (m_fad_size.value+blockDim.x-threadIdx.x-1) / blockDim.x;
835#else
836 const unsigned index = m_fad_index;
837 const unsigned strd = m_fad_stride;
838 const unsigned size = m_fad_size.value;
839#endif
840 return reference_type( beg + index
841 , beg + m_original_fad_size
842 , size
843 , strd ); }
844
845 template< typename I0 , typename I1 >
846 KOKKOS_FORCEINLINE_FUNCTION
847 typename std::enable_if< Kokkos::Impl::are_integral<I0,I1>::value &&
848 !is_layout_left, reference_type>::type
849 reference( const I0 & i0 , const I1 & i1 ) const
850 { pointer_type beg = m_impl_handle + m_array_offset(i0,i1,0);
851#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && ( defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
852 const unsigned index = threadIdx.x;
853 const unsigned strd = blockDim.x;
854 const unsigned size = (m_fad_size.value+blockDim.x-threadIdx.x-1) / blockDim.x;
855#else
856 const unsigned index = m_fad_index;
857 const unsigned strd = m_fad_stride;
858 const unsigned size = m_fad_size.value;
859#endif
860 return reference_type( beg + index
861 , beg + m_original_fad_size
862 , size
863 , strd ); }
864
865
866 template< typename I0 , typename I1 , typename I2 >
867 KOKKOS_FORCEINLINE_FUNCTION
868 typename std::enable_if< Kokkos::Impl::are_integral<I0,I1,I2>::value &&
869 is_layout_left, reference_type>::type
870 reference( const I0 & i0 , const I1 & i1 , const I2 & i2 ) const
871 { pointer_type beg = m_impl_handle + m_array_offset(0,i0,i1,i2);
872#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && ( defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
873 const unsigned index = threadIdx.x;
874 const unsigned strd = blockDim.x;
875 const unsigned size = (m_fad_size.value+blockDim.x-threadIdx.x-1) / blockDim.x;
876#else
877 const unsigned index = m_fad_index;
878 const unsigned strd = m_fad_stride;
879 const unsigned size = m_fad_size.value;
880#endif
881 return reference_type( beg + index
882 , beg + m_original_fad_size
883 , size
884 , strd ); }
885
886 template< typename I0 , typename I1 , typename I2 >
887 KOKKOS_FORCEINLINE_FUNCTION
888 typename std::enable_if< Kokkos::Impl::are_integral<I0,I1,I2>::value &&
889 !is_layout_left, reference_type>::type
890 reference( const I0 & i0 , const I1 & i1 , const I2 & i2 ) const
891 { pointer_type beg = m_impl_handle + m_array_offset(i0,i1,i2,0);
892#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && ( defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
893 const unsigned index = threadIdx.x;
894 const unsigned strd = blockDim.x;
895 const unsigned size = (m_fad_size.value+blockDim.x-threadIdx.x-1) / blockDim.x;
896#else
897 const unsigned index = m_fad_index;
898 const unsigned strd = m_fad_stride;
899 const unsigned size = m_fad_size.value;
900#endif
901 return reference_type( beg + index
902 , beg + m_original_fad_size
903 , size
904 , strd ); }
905
906 template< typename I0 , typename I1 , typename I2 , typename I3 >
907 KOKKOS_FORCEINLINE_FUNCTION
908 typename std::enable_if< Kokkos::Impl::are_integral<I0,I1,I2,I3>::value &&
909 is_layout_left, reference_type>::type
910 reference( const I0 & i0 , const I1 & i1 , const I2 & i2 , const I3 & i3 ) const
911 { pointer_type beg = m_impl_handle + m_array_offset(0,i0,i1,i2,i3);
912#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && ( defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
913 const unsigned index = threadIdx.x;
914 const unsigned strd = blockDim.x;
915 const unsigned size = (m_fad_size.value+blockDim.x-threadIdx.x-1) / blockDim.x;
916#else
917 const unsigned index = m_fad_index;
918 const unsigned strd = m_fad_stride;
919 const unsigned size = m_fad_size.value;
920#endif
921 return reference_type( beg + index
922 , beg + m_original_fad_size
923 , size
924 , strd ); }
925
926 template< typename I0 , typename I1 , typename I2 , typename I3 >
927 KOKKOS_FORCEINLINE_FUNCTION
928 typename std::enable_if< Kokkos::Impl::are_integral<I0,I1,I2,I3>::value &&
929 !is_layout_left, reference_type>::type
930 reference( const I0 & i0 , const I1 & i1 , const I2 & i2 , const I3 & i3 ) const
931 { pointer_type beg = m_impl_handle + m_array_offset(i0,i1,i2,i3,0);
932#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && ( defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
933 const unsigned index = threadIdx.x;
934 const unsigned strd = blockDim.x;
935 const unsigned size = (m_fad_size.value+blockDim.x-threadIdx.x-1) / blockDim.x;
936#else
937 const unsigned index = m_fad_index;
938 const unsigned strd = m_fad_stride;
939 const unsigned size = m_fad_size.value;
940#endif
941 return reference_type( beg + index
942 , beg + m_original_fad_size
943 , size
944 , strd ); }
945
946 template< typename I0 , typename I1 , typename I2 , typename I3
947 , typename I4 >
948 KOKKOS_FORCEINLINE_FUNCTION
949 typename std::enable_if< Kokkos::Impl::are_integral<I0,I1,I2,I3,I4>::value &&
950 is_layout_left, reference_type>::type
951 reference( const I0 & i0 , const I1 & i1 , const I2 & i2 , const I3 & i3
952 , const I4 & i4 ) const
953 { pointer_type beg = m_impl_handle + m_array_offset(0,i0,i1,i2,i3,i4);
954#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && ( defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
955 const unsigned index = threadIdx.x;
956 const unsigned strd = blockDim.x;
957 const unsigned size = (m_fad_size.value+blockDim.x-threadIdx.x-1) / blockDim.x;
958#else
959 const unsigned index = m_fad_index;
960 const unsigned strd = m_fad_stride;
961 const unsigned size = m_fad_size.value;
962#endif
963 return reference_type( beg + index
964 , beg + m_original_fad_size
965 , size
966 , strd ); }
967
968 template< typename I0 , typename I1 , typename I2 , typename I3
969 , typename I4 >
970 KOKKOS_FORCEINLINE_FUNCTION
971 typename std::enable_if< Kokkos::Impl::are_integral<I0,I1,I2,I3,I4>::value &&
972 !is_layout_left, reference_type>::type
973 reference( const I0 & i0 , const I1 & i1 , const I2 & i2 , const I3 & i3
974 , const I4 & i4 ) const
975 { pointer_type beg = m_impl_handle + m_array_offset(i0,i1,i2,i3,i4,0);
976#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && ( defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
977 const unsigned index = threadIdx.x;
978 const unsigned strd = blockDim.x;
979 const unsigned size = (m_fad_size.value+blockDim.x-threadIdx.x-1) / blockDim.x;
980#else
981 const unsigned index = m_fad_index;
982 const unsigned strd = m_fad_stride;
983 const unsigned size = m_fad_size.value;
984#endif
985 return reference_type( beg + index
986 , beg + m_original_fad_size
987 , size
988 , strd ); }
989
990 template< typename I0 , typename I1 , typename I2 , typename I3
991 , typename I4 , typename I5 >
992 KOKKOS_FORCEINLINE_FUNCTION
993 typename std::enable_if< Kokkos::Impl::are_integral<I0,I1,I2,I3,I4,I5>::value &&
994 is_layout_left, reference_type>::type
995 reference( const I0 & i0 , const I1 & i1 , const I2 & i2 , const I3 & i3
996 , const I4 & i4 , const I5 & i5 ) const
997 { pointer_type beg = m_impl_handle + m_array_offset(0,i0,i1,i2,i3,i4,i5);
998#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && ( defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
999 const unsigned index = threadIdx.x;
1000 const unsigned strd = blockDim.x;
1001 const unsigned size = (m_fad_size.value+blockDim.x-threadIdx.x-1) / blockDim.x;
1002#else
1003 const unsigned index = m_fad_index;
1004 const unsigned strd = m_fad_stride;
1005 const unsigned size = m_fad_size.value;
1006#endif
1007 return reference_type( beg + index
1008 , beg + m_original_fad_size
1009 , size
1010 , strd ); }
1011
1012 template< typename I0 , typename I1 , typename I2 , typename I3
1013 , typename I4 , typename I5 >
1014 KOKKOS_FORCEINLINE_FUNCTION
1015 typename std::enable_if< Kokkos::Impl::are_integral<I0,I1,I2,I3,I4,I5>::value &&
1016 !is_layout_left, reference_type>::type
1017 reference( const I0 & i0 , const I1 & i1 , const I2 & i2 , const I3 & i3
1018 , const I4 & i4 , const I5 & i5 ) const
1019 { pointer_type beg = m_impl_handle + m_array_offset(i0,i1,i2,i3,i4,i5,0);
1020#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && ( defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
1021 const unsigned index = threadIdx.x;
1022 const unsigned strd = blockDim.x;
1023 const unsigned size = (m_fad_size.value+blockDim.x-threadIdx.x-1) / blockDim.x;
1024#else
1025 const unsigned index = m_fad_index;
1026 const unsigned strd = m_fad_stride;
1027 const unsigned size = m_fad_size.value;
1028#endif
1029 return reference_type( beg + index
1030 , beg + m_original_fad_size
1031 , size
1032 , strd ); }
1033
1034 template< typename I0 , typename I1 , typename I2 , typename I3
1035 , typename I4 , typename I5 , typename I6 >
1036 KOKKOS_FORCEINLINE_FUNCTION
1037 typename std::enable_if< Kokkos::Impl::are_integral<I0,I1,I2,I3,I4,I5,I6>::value &&
1038 is_layout_left, reference_type>::type
1039 reference( const I0 & i0 , const I1 & i1 , const I2 & i2 , const I3 & i3
1040 , const I4 & i4 , const I5 & i5 , const I6 & i6 ) const
1041 { pointer_type beg = m_impl_handle + m_array_offset(0,i0,i1,i2,i3,i4,i5,i6);
1042#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && ( defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
1043 const unsigned index = threadIdx.x;
1044 const unsigned strd = blockDim.x;
1045 const unsigned size = (m_fad_size.value+blockDim.x-threadIdx.x-1) / blockDim.x;
1046#else
1047 const unsigned index = m_fad_index;
1048 const unsigned strd = m_fad_stride;
1049 const unsigned size = m_fad_size.value;
1050#endif
1051 return reference_type( beg + index
1052 , beg + m_original_fad_size
1053 , size
1054 , strd ); }
1055
1056 template< typename I0 , typename I1 , typename I2 , typename I3
1057 , typename I4 , typename I5 , typename I6 >
1058 KOKKOS_FORCEINLINE_FUNCTION
1059 typename std::enable_if< Kokkos::Impl::are_integral<I0,I1,I2,I3,I4,I5,I6>::value &&
1060 !is_layout_left, reference_type>::type
1061 reference( const I0 & i0 , const I1 & i1 , const I2 & i2 , const I3 & i3
1062 , const I4 & i4 , const I5 & i5 , const I6 & i6 ) const
1063 { pointer_type beg = m_impl_handle + m_array_offset(i0,i1,i2,i3,i4,i5,i6,0);
1064#if defined(SACADO_VIEW_CUDA_HIERARCHICAL) && ( defined(__CUDA_ARCH__) || defined(__HIP_DEVICE_COMPILE__) )
1065 const unsigned index = threadIdx.x;
1066 const unsigned strd = blockDim.x;
1067 const unsigned size = (m_fad_size.value+blockDim.x-threadIdx.x-1) / blockDim.x;
1068#else
1069 const unsigned index = m_fad_index;
1070 const unsigned strd = m_fad_stride;
1071 const unsigned size = m_fad_size.value;
1072#endif
1073 return reference_type( beg + index
1074 , beg + m_original_fad_size
1075 , size
1076 , strd ); }
1077
1078 //----------------------------------------
1079
1081 KOKKOS_INLINE_FUNCTION
1082 static constexpr size_t memory_span( typename Traits::array_layout const & layout )
1083 {
1084 // Do not introduce padding...
1085 typedef std::integral_constant< unsigned , 0 > padding ;
1086 return array_offset_type(
1087 padding() ,
1088 create_fad_array_layout<unsigned(Rank), unsigned(FadStaticDimension)>( layout ) ).span() * sizeof(fad_value_type);
1089 }
1090
1091 //----------------------------------------
1092
1093 KOKKOS_DEFAULTED_FUNCTION ~ViewMapping() = default ;
1094 KOKKOS_INLINE_FUNCTION ViewMapping() : m_impl_handle(0) , m_impl_offset() , m_array_offset() , m_fad_size(0) , m_original_fad_size(0) , m_fad_stride(1) , m_fad_index(0) {}
1095
1096 KOKKOS_DEFAULTED_FUNCTION ViewMapping( const ViewMapping & ) = default ;
1097 KOKKOS_DEFAULTED_FUNCTION ViewMapping & operator = ( const ViewMapping & ) = default ;
1098
1099 KOKKOS_DEFAULTED_FUNCTION ViewMapping( ViewMapping && ) = default ;
1100 KOKKOS_DEFAULTED_FUNCTION ViewMapping & operator = ( ViewMapping && ) = default ;
1101
1102 template< class ... P >
1103 KOKKOS_INLINE_FUNCTION
1104 ViewMapping
1105 ( ViewCtorProp< P ... > const & prop
1106 , typename Traits::array_layout const & local_layout
1107 )
1108 : m_impl_handle( ( (ViewCtorProp<void,pointer_type> const &) prop ).value )
1109 , m_impl_offset( std::integral_constant< unsigned , 0 >()
1110 , local_layout )
1111 , m_array_offset(
1112 std::integral_constant< unsigned , 0 >()
1113 , create_fad_array_layout<unsigned(Rank), unsigned(FadStaticDimension)>( local_layout ) )
1114 , m_fad_size( getFadDimension<unsigned(Rank)>( m_array_offset ) - 1 )
1115 , m_original_fad_size( m_fad_size.value )
1116 , m_fad_stride( 1 )
1117 , m_fad_index( 0 )
1118 {
1119 const unsigned fad_dim =
1120 getFadDimension<unsigned(Rank)>( m_array_offset );
1121 if (unsigned(FadStaticDimension) == 0 && fad_dim == 0)
1122 Kokkos::abort("invalid fad dimension (0) supplied!");
1123 }
1124
1125 //----------------------------------------
1126 /* Allocate and construct mapped array.
1127 * Allocate via shared allocation record and
1128 * return that record for allocation tracking.
1129 */
1130 template< class ... P >
1131 SharedAllocationRecord<> *
1132 allocate_shared( ViewCtorProp< P... > const & prop
1133 , typename Traits::array_layout const & local_layout
1134 , bool execution_space_specified)
1135 {
1136 typedef ViewCtorProp< P... > ctor_prop ;
1137
1138 typedef typename ctor_prop::execution_space execution_space ;
1139 typedef typename Traits::memory_space memory_space ;
1140 typedef ViewValueFunctor< execution_space , fad_value_type > functor_type ;
1141 typedef SharedAllocationRecord< memory_space , functor_type > record_type ;
1142
1143 // Disallow padding
1144 typedef std::integral_constant< unsigned , 0 > padding ;
1145
1146 // Check if ViewCtorProp has CommonViewAllocProp - if so, retrieve the fad_size and append to layout
1147 enum { test_traits_check = Kokkos::Impl::check_has_common_view_alloc_prop< P... >::value };
1148
1149 typename Traits::array_layout internal_layout =
1150 (test_traits_check == true)
1151 ? Kokkos::Impl::appendFadToLayoutViewAllocHelper< Traits, P... >::returnNewLayoutPlusFad(prop, local_layout)
1152 : local_layout;
1153
1154 m_impl_offset = offset_type( padding(), internal_layout );
1155
1156 m_array_offset =
1157 array_offset_type( padding() ,
1158 create_fad_array_layout<unsigned(Rank), unsigned(FadStaticDimension)>( internal_layout ) );
1159 const unsigned fad_dim =
1160 getFadDimension<unsigned(Rank)>( m_array_offset );
1161 if (unsigned(FadStaticDimension) == 0 && fad_dim == 0)
1162 Kokkos::abort("invalid fad dimension (0) supplied!");
1163 m_fad_size = fad_dim - 1 ;
1164 m_original_fad_size = m_fad_size.value ;
1165 m_fad_stride = 1;
1166 m_fad_index = 0;
1167
1168 const size_t alloc_size = m_array_offset.span() * sizeof(fad_value_type);
1169
1170 // Create shared memory tracking record with allocate memory from the memory space
1171 record_type * const record =
1172 record_type::allocate( ( (ViewCtorProp<void,memory_space> const &) prop ).value
1173 , ( (ViewCtorProp<void,std::string> const &) prop ).value
1174 , alloc_size );
1175
1176 // Only set the the pointer and initialize if the allocation is non-zero.
1177 // May be zero if one of the dimensions is zero.
1178 if ( alloc_size ) {
1179
1180 m_impl_handle = handle_type( reinterpret_cast< pointer_type >( record->data() ) );
1181
1182 if ( ctor_prop::initialize ) {
1183 // Assume destruction is only required when construction is requested.
1184 // The ViewValueFunctor has both value construction and destruction operators.
1185 if (execution_space_specified)
1186 record->m_destroy = functor_type( ( (ViewCtorProp<void,execution_space> const &) prop).value
1187 , (fad_value_type *) m_impl_handle
1188 , m_array_offset.span()
1189 , record->get_label()
1190 );
1191 else
1192 record->m_destroy = functor_type((fad_value_type *) m_impl_handle
1193 , m_array_offset.span()
1194 , record->get_label()
1195 );
1196
1197 // Construct values
1198 record->m_destroy.construct_shared_allocation();
1199 }
1200 }
1201
1202 return record ;
1203 }
1204
1205};
1206
1207} // namespace Impl
1208} // namespace Kokkos
1209
1210//----------------------------------------------------------------------------
1211
1212namespace Kokkos {
1213namespace Impl {
1214
1219template< class DstTraits , class SrcTraits >
1220class ViewMapping< DstTraits , SrcTraits ,
1221 typename std::enable_if<(
1222 Kokkos::Impl::MemorySpaceAccess
1223 < typename DstTraits::memory_space
1224 , typename SrcTraits::memory_space >::assignable
1225 &&
1226 // Destination view has FAD
1227 std::is_same< typename DstTraits::specialize
1228 , ViewSpecializeSacadoFadContiguous >::value
1229 &&
1230 // Source view has FAD
1231 std::is_same< typename SrcTraits::specialize
1232 , ViewSpecializeSacadoFadContiguous >::value
1233 )
1234 , typename DstTraits::specialize
1235 >::type >
1236{
1237public:
1238
1239 enum { is_assignable = true };
1240 enum { is_assignable_data_type = true };
1241
1242 typedef Kokkos::Impl::SharedAllocationTracker TrackType ;
1243 typedef ViewMapping< DstTraits , typename DstTraits::specialize > DstType ;
1244 typedef ViewMapping< SrcTraits , typename SrcTraits::specialize > SrcFadType ;
1245
1246 template< class DstType >
1247 KOKKOS_INLINE_FUNCTION static
1248 void assign( DstType & dst
1249 , const SrcFadType & src
1250 , const TrackType & )
1251 {
1252 static_assert(
1253 (
1254 std::is_same< typename DstTraits::array_layout
1255 , Kokkos::LayoutLeft >::value ||
1256 std::is_same< typename DstTraits::array_layout
1257 , Kokkos::LayoutRight >::value ||
1258 std::is_same< typename DstTraits::array_layout
1259 , Kokkos::LayoutStride >::value
1260 )
1261 &&
1262 (
1263 std::is_same< typename SrcTraits::array_layout
1264 , Kokkos::LayoutLeft >::value ||
1265 std::is_same< typename SrcTraits::array_layout
1266 , Kokkos::LayoutRight >::value ||
1267 std::is_same< typename SrcTraits::array_layout
1268 , Kokkos::LayoutStride >::value
1269 )
1270 , "View of FAD requires LayoutLeft, LayoutRight, or LayoutStride" );
1271
1272 static_assert(
1273 std::is_same< typename DstTraits::array_layout
1274 , typename SrcTraits::array_layout >::value ||
1275 std::is_same< typename DstTraits::array_layout
1276 , Kokkos::LayoutStride >::value ,
1277 "View assignment must have compatible layout" );
1278
1279 static_assert(
1280 std::is_same< typename DstTraits::value_type
1281 , typename SrcTraits::value_type >::value ||
1282 std::is_same< typename DstTraits::value_type
1283 , typename SrcTraits::const_value_type >::value ,
1284 "View assignment must have same value type or const = non-const" );
1285
1286 static_assert(
1287 ViewDimensionAssignable
1288 < typename DstType::offset_type::dimension_type
1289 , typename SrcFadType::offset_type::dimension_type >::value ,
1290 "View assignment must have compatible dimensions" );
1291
1292 static_assert(
1293 ViewDimensionAssignable
1294 < typename DstType::array_offset_type::dimension_type
1295 , typename SrcFadType::array_offset_type::dimension_type >::value ,
1296 "View assignment must have compatible dimensions" );
1297
1298 typedef typename DstType::offset_type dst_offset_type ;
1299 typedef typename DstType::array_offset_type dst_array_offset_type ;
1300
1301 dst.m_impl_handle = src.m_impl_handle ;
1302 dst.m_impl_offset = dst_offset_type( src.m_impl_offset );
1303 dst.m_array_offset = dst_array_offset_type( src.m_array_offset );
1304 dst.m_fad_size = src.m_fad_size.value ;
1305 dst.m_original_fad_size = src.m_original_fad_size ;
1306 dst.m_fad_stride = src.m_fad_stride ;
1307 dst.m_fad_index = src.m_fad_index ;
1308 }
1309};
1310
1315template< class DstTraits , class SrcTraits >
1316class ViewMapping< DstTraits , SrcTraits ,
1317 typename std::enable_if<(
1318 std::is_same< typename DstTraits::memory_space
1319 , typename SrcTraits::memory_space >::value
1320 &&
1321 // Destination view has FAD
1322 std::is_same< typename DstTraits::specialize
1323 , ViewSpecializeSacadoFad >::value
1324 &&
1325 // Source view has FAD contiguous
1326 std::is_same< typename SrcTraits::specialize
1327 , ViewSpecializeSacadoFadContiguous >::value
1328 &&
1329 // Destination view is LayoutStride
1330 std::is_same< typename DstTraits::array_layout
1331 , Kokkos::LayoutStride >::value
1332 )
1333 , typename DstTraits::specialize
1334 >::type >
1335{
1336public:
1337
1338 enum { is_assignable = true };
1339 enum { is_assignable_data_type = true };
1340
1341 typedef Kokkos::Impl::SharedAllocationTracker TrackType ;
1342 typedef ViewMapping< DstTraits , typename DstTraits::specialize > DstType ;
1343 typedef ViewMapping< SrcTraits , typename SrcTraits::specialize > SrcFadType ;
1344
1345 template< class DstType >
1346 KOKKOS_INLINE_FUNCTION static
1347 void assign( DstType & dst
1348 , const SrcFadType & src
1349 , const TrackType & )
1350 {
1351 static_assert(
1352 std::is_same< typename SrcTraits::array_layout
1353 , Kokkos::LayoutLeft >::value ||
1354 std::is_same< typename SrcTraits::array_layout
1355 , Kokkos::LayoutRight >::value ||
1356 std::is_same< typename SrcTraits::array_layout
1357 , Kokkos::LayoutStride >::value ,
1358 "View of FAD requires LayoutLeft, LayoutRight, or LayoutStride" );
1359
1360 static_assert(
1361 std::is_same< typename DstTraits::value_type
1362 , typename SrcTraits::value_type >::value ||
1363 std::is_same< typename DstTraits::value_type
1364 , typename SrcTraits::const_value_type >::value ,
1365 "View assignment must have same value type or const = non-const" );
1366
1367 static_assert(
1368 DstTraits::dimension::rank == SrcTraits::dimension::rank,
1369 "View assignment must have same rank" );
1370
1371 typedef typename DstType::array_offset_type dst_offset_type ;
1372
1373 dst.m_impl_handle = src.m_impl_handle ;
1374 dst.m_fad_size = src.m_fad_size.value ;
1375 dst.m_fad_stride = src.m_fad_stride ;
1376 dst.m_impl_offset = src.m_impl_offset;
1377
1378 size_t N[8], S[8];
1379 N[0] = src.m_array_offset.dimension_0();
1380 N[1] = src.m_array_offset.dimension_1();
1381 N[2] = src.m_array_offset.dimension_2();
1382 N[3] = src.m_array_offset.dimension_3();
1383 N[4] = src.m_array_offset.dimension_4();
1384 N[5] = src.m_array_offset.dimension_5();
1385 N[6] = src.m_array_offset.dimension_6();
1386 N[7] = src.m_array_offset.dimension_7();
1387 S[0] = src.m_array_offset.stride_0();
1388 S[1] = src.m_array_offset.stride_1();
1389 S[2] = src.m_array_offset.stride_2();
1390 S[3] = src.m_array_offset.stride_3();
1391 S[4] = src.m_array_offset.stride_4();
1392 S[5] = src.m_array_offset.stride_5();
1393 S[6] = src.m_array_offset.stride_6();
1394 S[7] = src.m_array_offset.stride_7();
1395
1396 // For LayoutLeft, we have to move the Sacado dimension from the first
1397 // to the last
1398 if (std::is_same< typename SrcTraits::array_layout
1399 , Kokkos::LayoutLeft >::value)
1400 {
1401 const size_t N_fad = N[0];
1402 const size_t S_fad = S[0];
1403 for (int i=0; i<7; ++i) {
1404 N[i] = N[i+1];
1405 S[i] = S[i+1];
1406 }
1407 N[DstTraits::dimension::rank] = N_fad;
1408 S[DstTraits::dimension::rank] = S_fad;
1409 }
1410 Kokkos::LayoutStride ls( N[0], S[0],
1411 N[1], S[1],
1412 N[2], S[2],
1413 N[3], S[3],
1414 N[4], S[4],
1415 N[5], S[5],
1416 N[6], S[6],
1417 N[7], S[7] );
1418 dst.m_array_offset = dst_offset_type(std::integral_constant<unsigned,0>(), ls);
1419 }
1420};
1421
1426template< class DstTraits , class SrcTraits >
1427class ViewMapping< DstTraits , SrcTraits ,
1428 typename std::enable_if<(
1429 std::is_same< typename DstTraits::memory_space
1430 , typename SrcTraits::memory_space >::value
1431 &&
1432 // Destination view has ordinary
1433 std::is_same< typename DstTraits::specialize , void >::value
1434 &&
1435 // Source view has FAD only
1436 std::is_same< typename SrcTraits::specialize
1437 , ViewSpecializeSacadoFadContiguous >::value
1438 )
1439 , typename DstTraits::specialize
1440 >::type >
1441{
1442public:
1443
1444 enum { is_assignable = true };
1445 enum { is_assignable_data_type = true };
1446
1447 typedef Kokkos::Impl::SharedAllocationTracker TrackType ;
1448 typedef ViewMapping< DstTraits , typename DstTraits::specialize > DstType ;
1449 typedef ViewMapping< SrcTraits , typename SrcTraits::specialize > SrcFadType ;
1450
1451
1452 // Helpers to assign, and generate if necessary, ViewOffset to the dst map
1453 // These are necessary to use Kokkos' deep_copy with nested fads
1454 template < class DstType, class SrcFadType, class Enable = void >
1455 struct AssignOffset;
1456
1457 template < class DstType, class SrcFadType >
1458 struct AssignOffset< DstType, SrcFadType, typename std::enable_if< ((int)DstType::offset_type::dimension_type::rank != (int)SrcFadType::array_offset_type::dimension_type::rank) >::type >
1459 {
1460 // ViewOffset's Dimensions Ranks do not match
1461 KOKKOS_INLINE_FUNCTION
1462 static void assign( DstType & dst, const SrcFadType & src )
1463 {
1464 typedef typename SrcTraits::value_type TraitsValueType;
1465
1468 )
1469 {
1470
1471 typedef typename DstType::offset_type::array_layout DstLayoutType;
1472 //typedef typename ViewArrayLayoutSelector<typename DstType::offset_type::array_layout>::type DstLayoutType;
1473 typedef typename SrcFadType::array_offset_type::dimension_type SrcViewDimension;
1474
1475 // This is the static dimension of the inner fad, missing from ViewDimension
1477
1478 static constexpr bool is_layout_left =
1479 std::is_same< DstLayoutType, Kokkos::LayoutLeft>::value;
1480
1481 typedef typename std::conditional< is_layout_left,
1482 typename SrcViewDimension:: template prepend< InnerStaticDim+1 >::type,
1483 typename SrcViewDimension:: template append < InnerStaticDim+1 >::type
1484 >::type SrcViewDimensionAppended;
1485
1486 typedef std::integral_constant< unsigned , 0 > padding ;
1487
1488 typedef ViewOffset< SrcViewDimensionAppended, DstLayoutType > TmpOffsetType;
1489
1490 auto src_layout = src.m_array_offset.layout();
1491
1492 if ( is_layout_left ) {
1493 auto prepend_layout = Kokkos::Impl::prependFadToLayout< DstLayoutType >::returnNewLayoutPlusFad(src_layout, InnerStaticDim+1);
1494 TmpOffsetType offset_tmp( padding(), prepend_layout );
1495 dst.m_impl_offset = offset_tmp;
1496 }
1497 else {
1498 TmpOffsetType offset_tmp( padding(), src_layout );
1499 dst.m_impl_offset = offset_tmp;
1500 }
1501 } else {
1502 Kokkos::abort("Sacado error: Applying AssignOffset for case with nested Fads, but without nested Fads - something went wrong");
1503 }
1504 }
1505 };
1506
1507 template < class DstType, class SrcFadType >
1508 struct AssignOffset< DstType, SrcFadType, typename std::enable_if< ((int)DstType::offset_type::dimension_type::rank == (int)SrcFadType::array_offset_type::dimension_type::rank) >::type >
1509 {
1510 KOKKOS_INLINE_FUNCTION
1511 static void assign( DstType & dst, const SrcFadType & src )
1512 {
1513 typedef typename DstType::offset_type dst_offset_type ;
1514 dst.m_impl_offset = dst_offset_type( src.m_array_offset );
1515 }
1516 };
1517
1518 template< class DstType >
1519 KOKKOS_INLINE_FUNCTION static
1520 void assign( DstType & dst
1521 , const SrcFadType & src
1522 , const TrackType &
1523 )
1524 {
1525 static_assert(
1526 (
1527 std::is_same< typename DstTraits::array_layout
1528 , Kokkos::LayoutLeft >::value ||
1529 std::is_same< typename DstTraits::array_layout
1530 , Kokkos::LayoutRight >::value ||
1531 std::is_same< typename DstTraits::array_layout
1532 , Kokkos::LayoutStride >::value
1533 )
1534 &&
1535 (
1536 std::is_same< typename SrcTraits::array_layout
1537 , Kokkos::LayoutLeft >::value ||
1538 std::is_same< typename SrcTraits::array_layout
1539 , Kokkos::LayoutRight >::value ||
1540 std::is_same< typename SrcTraits::array_layout
1541 , Kokkos::LayoutStride >::value
1542 )
1543 , "View of FAD requires LayoutLeft, LayoutRight, or LayoutStride" );
1544
1545 static_assert(
1546 std::is_same< typename DstTraits::array_layout
1547 , typename SrcTraits::array_layout >::value ||
1548 std::is_same< typename DstTraits::array_layout
1549 , Kokkos::LayoutStride >::value ,
1550 "View assignment must have compatible layout" );
1551
1552 if ( src.m_fad_index != 0 || src.m_fad_stride != 1 ) {
1553 Kokkos::abort("\n\n ****** Kokkos::View< Sacado::Fad ... > Cannot assign to array with partitioned view ******\n\n");
1554 }
1555
1556 AssignOffset< DstType, SrcFadType >::assign( dst, src );
1557 dst.m_impl_handle = reinterpret_cast< typename DstType::handle_type >(src.m_impl_handle) ;
1558 }
1559};
1560
1561} // namespace Impl
1562} // namespace Kokkos
1563
1564//----------------------------------------------------------------------------
1565
1566namespace Kokkos {
1567namespace Impl {
1568
1569// Rules for subview arguments and layouts matching
1570
1571template<class LayoutDest, class LayoutSrc, int RankDest, int RankSrc, int CurrentArg, class ... SubViewArgs>
1572struct SubviewLegalArgsCompileTime<Kokkos::LayoutContiguous<LayoutDest>,LayoutSrc,RankDest,RankSrc,CurrentArg,SubViewArgs...> {
1573 enum { value = SubviewLegalArgsCompileTime<LayoutDest,LayoutSrc,RankDest,RankSrc,CurrentArg,SubViewArgs...>::value };
1574};
1575
1576template<class LayoutDest, class LayoutSrc, int RankDest, int RankSrc, int CurrentArg, class ... SubViewArgs>
1577struct SubviewLegalArgsCompileTime<LayoutDest,Kokkos::LayoutContiguous<LayoutSrc>,RankDest,RankSrc,CurrentArg,SubViewArgs...> {
1578 enum { value = SubviewLegalArgsCompileTime<LayoutDest,LayoutSrc,RankDest,RankSrc,CurrentArg,SubViewArgs...>::value };
1579};
1580
1581template<class LayoutDest, class LayoutSrc, int RankDest, int RankSrc, int CurrentArg, class ... SubViewArgs>
1582struct SubviewLegalArgsCompileTime<Kokkos::LayoutContiguous<LayoutDest>,Kokkos::LayoutContiguous<LayoutSrc>,RankDest,RankSrc,CurrentArg,SubViewArgs...> {
1583 enum { value = SubviewLegalArgsCompileTime<LayoutDest,LayoutSrc,RankDest,RankSrc,CurrentArg,SubViewArgs...>::value };
1584};
1585
1586// Subview mapping
1587
1588template< class SrcTraits , class Arg0 , class ... Args >
1589struct ViewMapping
1590 < typename std::enable_if<(
1591 // Source view has FAD only
1592 std::is_same< typename SrcTraits::specialize
1593 , ViewSpecializeSacadoFadContiguous >::value
1594 &&
1595 (
1596 std::is_same< typename SrcTraits::array_layout
1597 , Kokkos::LayoutLeft >::value ||
1598 std::is_same< typename SrcTraits::array_layout
1599 , Kokkos::LayoutRight >::value ||
1600 std::is_same< typename SrcTraits::array_layout
1601 , Kokkos::LayoutStride >::value
1602 )
1603 && !Sacado::Fad::is_fad_partition<Arg0>::value
1604 )
1605 >::type
1606 , SrcTraits
1607 , Arg0, Args ... >
1608{
1609private:
1610
1611 static_assert( SrcTraits::rank == sizeof...(Args)+1 , "" );
1612
1613 enum
1614 { RZ = false
1615 , R0 = bool(is_integral_extent<0,Arg0,Args...>::value)
1616 , R1 = bool(is_integral_extent<1,Arg0,Args...>::value)
1617 , R2 = bool(is_integral_extent<2,Arg0,Args...>::value)
1618 , R3 = bool(is_integral_extent<3,Arg0,Args...>::value)
1619 , R4 = bool(is_integral_extent<4,Arg0,Args...>::value)
1620 , R5 = bool(is_integral_extent<5,Arg0,Args...>::value)
1621 , R6 = bool(is_integral_extent<6,Arg0,Args...>::value)
1622 };
1623
1624 // Public rank
1625 enum { rank = unsigned(R0) + unsigned(R1) + unsigned(R2) + unsigned(R3)
1626 + unsigned(R4) + unsigned(R5) + unsigned(R6) };
1627
1628 // Whether right-most non-FAD rank is a range.
1629 enum { R0_rev = ( 0 == SrcTraits::rank ? RZ : (
1630 1 == SrcTraits::rank ? R0 : (
1631 2 == SrcTraits::rank ? R1 : (
1632 3 == SrcTraits::rank ? R2 : (
1633 4 == SrcTraits::rank ? R3 : (
1634 5 == SrcTraits::rank ? R4 : (
1635 6 == SrcTraits::rank ? R5 : R6 ))))))) };
1636
1637 // Subview's layout
1638 typedef typename std::conditional<
1639 ( /* Same array layout IF */
1640 ( rank == 0 ) /* output rank zero */
1641 ||
1642 // OutputRank 1 or 2, InputLayout Left, Interval 0
1643 // because single stride one or second index has a stride.
1644 ( rank <= 2 && R0 && std::is_same< typename SrcTraits::array_layout , Kokkos::LayoutLeft >::value )
1645 ||
1646 // OutputRank 1 or 2, InputLayout Right, Interval [InputRank-1]
1647 // because single stride one or second index has a stride.
1648 ( rank <= 2 && R0_rev && std::is_same< typename SrcTraits::array_layout , Kokkos::LayoutRight >::value )
1650 >::type array_layout ;
1651
1652 typedef typename SrcTraits::value_type fad_type ;
1653
1654 typedef typename std::conditional< rank == 0 , fad_type ,
1655 typename std::conditional< rank == 1 , fad_type * ,
1656 typename std::conditional< rank == 2 , fad_type ** ,
1657 typename std::conditional< rank == 3 , fad_type *** ,
1658 typename std::conditional< rank == 4 , fad_type **** ,
1659 typename std::conditional< rank == 5 , fad_type ***** ,
1660 typename std::conditional< rank == 6 , fad_type ****** ,
1661 fad_type *******
1662 >::type >::type >::type >::type >::type >::type >::type
1663 data_type ;
1664
1665public:
1666
1667 typedef Kokkos::ViewTraits
1668 < data_type
1669 , array_layout
1670 , typename SrcTraits::device_type
1671 , typename SrcTraits::memory_traits > traits_type ;
1672
1673 typedef Kokkos::View
1674 < data_type
1675 , array_layout
1676 , typename SrcTraits::device_type
1677 , typename SrcTraits::memory_traits > type ;
1678
1679
1680 KOKKOS_INLINE_FUNCTION
1681 static void assign( ViewMapping< traits_type , typename traits_type::specialize > & dst
1682 , ViewMapping< SrcTraits , typename SrcTraits::specialize > const & src
1683 , Arg0 arg0 , Args ... args )
1684 {
1685 typedef ViewMapping< traits_type , typename traits_type::specialize > DstType ;
1686 typedef typename DstType::offset_type dst_offset_type ;
1687 typedef typename DstType::array_offset_type dst_array_offset_type ;
1688 typedef typename DstType::handle_type dst_handle_type ;
1689
1690 size_t offset;
1691 if (std::is_same< typename SrcTraits::array_layout, LayoutLeft >::value) {
1692 const SubviewExtents< SrcTraits::rank + 1 , rank + 1 >
1693 array_extents( src.m_array_offset.m_dim , Kokkos::ALL() , arg0 , args... );
1694 offset = src.m_array_offset( array_extents.domain_offset(0)
1695 , array_extents.domain_offset(1)
1696 , array_extents.domain_offset(2)
1697 , array_extents.domain_offset(3)
1698 , array_extents.domain_offset(4)
1699 , array_extents.domain_offset(5)
1700 , array_extents.domain_offset(6)
1701 , array_extents.domain_offset(7) );
1702 dst.m_array_offset = dst_array_offset_type( src.m_array_offset ,
1703 array_extents );
1704 }
1705 else {
1706 const SubviewExtents< SrcTraits::rank + 1 , rank + 1 >
1707 array_extents( src.m_array_offset.m_dim , arg0 , args... , Kokkos::ALL() );
1708 offset = src.m_array_offset( array_extents.domain_offset(0)
1709 , array_extents.domain_offset(1)
1710 , array_extents.domain_offset(2)
1711 , array_extents.domain_offset(3)
1712 , array_extents.domain_offset(4)
1713 , array_extents.domain_offset(5)
1714 , array_extents.domain_offset(6)
1715 , array_extents.domain_offset(7) );
1716 dst.m_array_offset = dst_array_offset_type( src.m_array_offset ,
1717 array_extents );
1718 }
1719
1720 const SubviewExtents< SrcTraits::rank , rank >
1721 extents( src.m_impl_offset.m_dim , arg0 , args... );
1722
1723 dst.m_impl_offset = dst_offset_type( src.m_impl_offset , extents );
1724 dst.m_impl_handle = dst_handle_type( src.m_impl_handle + offset );
1725 dst.m_fad_size = src.m_fad_size;
1726 dst.m_original_fad_size = src.m_original_fad_size;
1727 dst.m_fad_stride = src.m_fad_stride;
1728 dst.m_fad_index = src.m_fad_index;
1729 }
1730
1731};
1732
1733} // namespace Impl
1734} // namespace Kokkos
1735
1736//---------------------------------------------------------------------------
1737
1738namespace Kokkos {
1739namespace Impl {
1740
1741// Partition mapping
1742
1743template< class DataType, class ...P, unsigned Stride >
1744class ViewMapping<
1745 void,
1746 ViewTraits<DataType,P...> ,
1747 Sacado::Fad::Partition<Stride>
1748 >
1749{
1750public:
1751
1752 enum { is_assignable = true };
1753 enum { is_assignable_data_type = true };
1754
1755 typedef ViewTraits<DataType,P...> src_traits;
1756 typedef ViewMapping< src_traits , typename src_traits::specialize > src_type ;
1757
1758 typedef typename src_type::offset_type::dimension_type src_dimension;
1759 typedef typename src_traits::value_type fad_type;
1760 typedef typename Sacado::LocalScalarType<fad_type,Stride>::type strided_fad_type;
1761 typedef typename
1762 ViewDataType< strided_fad_type , src_dimension >::type strided_data_type;
1763 typedef ViewTraits<strided_data_type,P...> dst_traits;
1764 typedef View<strided_data_type,P...> type;
1765 typedef ViewMapping< dst_traits , typename dst_traits::specialize > dst_type ;
1766
1767 KOKKOS_INLINE_FUNCTION static
1768 void assign( dst_type & dst
1769 , const src_type & src
1770 , const Sacado::Fad::Partition<Stride> & part )
1771 {
1772 if ( Stride != part.stride && Stride != 0 ) {
1773 Kokkos::abort("\n\n ****** Kokkos::View< Sacado::Fad ... > Invalid size in partitioned view assignment ******\n\n");
1774 }
1775 if ( src.m_fad_stride != 1 ) {
1776 Kokkos::abort("\n\n ****** Kokkos::View< Sacado::Fad ... > Can't partition already partitioned view ******\n\n");
1777 }
1778
1779 dst.m_impl_handle = src.m_impl_handle ;
1780 dst.m_impl_offset = src.m_impl_offset ;
1781 dst.m_array_offset = src.m_array_offset ;
1782
1783 // Assuming the alignment was choosen correctly for the partitioning,
1784 // each partition should get the same size. This allows the use of SFad.
1785 dst.m_fad_size =
1786 (src.m_fad_size.value + part.stride-part.offset-1) / part.stride ;
1787
1788 dst.m_original_fad_size = src.m_original_fad_size ;
1789 dst.m_fad_stride = part.stride ;
1790 dst.m_fad_index = part.offset ;
1791 }
1792};
1793
1794} // namespace Impl
1795} // namespace Kokkos
1796
1797#endif // defined(HAVE_SACADO_VIEW_SPEC) && !defined(SACADO_DISABLE_FAD_VIEW_SPEC)
1798
1799#endif // defined(HAVE_SACADO_KOKKOSCORE)
1800
1801#endif /* #ifndef KOKKOS_EXPERIMENTAL_VIEW_SACADO_FAD_HPP */
View
expr expr expr bar false
expr expr dx(i)
#define T
Definition: Sacado_rad.hpp:573
#define D
Definition: Sacado_rad.hpp:577
const int N
const int fad_dim
int value
GeneralFad< DynamicStorage< T > > DFad
GeneralFad< StaticFixedStorage< T, Num > > SFad
GeneralFad< StaticStorage< T, Num > > SLFad
Base template specification for whether a type is a Fad type.
Base template specification for testing whether type is statically sized.
Base template specification for static size.
Get view type for any Fad type.