Sacado Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Fad_KokkosTests.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Sacado Package
4 //
5 // Copyright 2006 NTESS and the Sacado contributors.
6 // SPDX-License-Identifier: LGPL-2.1-or-later
7 // *****************************************************************************
8 // @HEADER
9 
11 
12 #include "Sacado.hpp"
14 
15 template <typename T>
16 struct is_sfad {
17  static const bool value = false;
18 };
19 
20 template <typename T, int N>
21 struct is_sfad< Sacado::Fad::SFad<T,N> > {
22  static const bool value = true;
23 };
24 
25 template <typename T>
26 struct is_dfad {
27  static const bool value = false;
28 };
29 
30 template <typename T>
31 struct is_dfad< Sacado::Fad::DFad<T> > {
32  static const bool value = true;
33 };
34 
35 template <typename FadType1, typename FadType2>
36 bool checkFads(const FadType1& x, const FadType2& x2,
37  Teuchos::FancyOStream& out, double tol = 1.0e-15)
38 {
39  bool success = true;
40 
41  // Check sizes match
42  TEUCHOS_TEST_EQUALITY(x.size(), x2.size(), out, success);
43 
44  // Check values match
45  TEUCHOS_TEST_FLOATING_EQUALITY(x.val(), x2.val(), tol, out, success);
46 
47  // Check derivatives match
48  for (int i=0; i<x.size(); ++i)
49  TEUCHOS_TEST_FLOATING_EQUALITY(x.dx(i), x2.dx(i), tol, out, success);
50 
51  return success;
52 }
53 
54 template <typename fadtype, typename ordinal>
55 inline
56 fadtype generate_fad( const ordinal num_rows,
57  const ordinal num_cols,
58  const ordinal fad_size,
59  const ordinal row,
60  const ordinal col )
61 {
62  typedef typename fadtype::value_type scalar;
63  fadtype x(fad_size, scalar(0.0));
64 
65  const scalar x_row = 100.0 + scalar(num_rows) / scalar(row+1);
66  const scalar x_col = 10.0 + scalar(num_cols) / scalar(col+1);
67  x.val() = x_row + x_col;
68  for (ordinal i=0; i<fad_size; ++i) {
69  const scalar x_fad = 1.0 + scalar(fad_size) / scalar(i+1);
70  x.fastAccessDx(i) = x_row + x_col + x_fad;
71  }
72  return x;
73 }
74 
75 #ifndef GLOBAL_FAD_SIZE
76 #define GLOBAL_FAD_SIZE 5
77 #endif
78 const int global_num_rows = 11;
79 const int global_num_cols = 7;
81 
82 // Kernel to multiply two views
83 template <typename InputViewType1,
84  typename InputViewType2 = InputViewType1,
85  typename OutputViewType = InputViewType1>
87  typedef typename InputViewType1::execution_space execution_space;
88  typedef typename InputViewType1::size_type size_type;
89  typedef Kokkos::RangePolicy< execution_space> range_policy_type;
90  typedef Kokkos::TeamPolicy< execution_space> team_policy_type;
91  typedef typename team_policy_type::member_type team_handle;
92 
93  const InputViewType1 m_v1;
94  const InputViewType2 m_v2;
95  const OutputViewType m_v3;
96  const bool m_update;
97 
98  MultiplyKernel(const InputViewType1 v1,
99  const InputViewType2 v2,
100  const OutputViewType v3,
101  const bool update) :
102  m_v1(v1), m_v2(v2), m_v3(v3), m_update(update) {};
103 
104  // Multiply entries for row 'i' with a value
105  KOKKOS_INLINE_FUNCTION
106  void operator() (const size_type i) const {
107  if (m_update)
108  m_v3(i) += m_v1(i)*m_v2(i);
109  else
110  m_v3(i) = m_v1(i)*m_v2(i);
111  }
112 
113  KOKKOS_INLINE_FUNCTION
114  void operator()( const team_handle& team ) const
115  {
116  const size_type i = team.league_rank()*team.team_size() + team.team_rank();
117  if (i < m_v1.extent(0))
118  (*this)(i);
119  }
120 
121  // Kernel launch
122  static void apply(const InputViewType1 v1,
123  const InputViewType2 v2,
124  const OutputViewType v3,
125  const bool update = false) {
126  const size_type nrow = v1.extent(0);
127 
128 #if defined (KOKKOS_ENABLE_CUDA) && defined (SACADO_VIEW_CUDA_HIERARCHICAL)
129  const size_type stride = Kokkos::ViewScalarStride<InputViewType1>::stride;
130  const bool use_team =
134  ( stride > 1 );
135 #elif defined (KOKKOS_ENABLE_CUDA) && defined (SACADO_VIEW_CUDA_HIERARCHICAL_DFAD)
136  const size_type stride = team_policy_type::vector_length_max(); // 32
137  const bool use_team =
142 #elif defined (KOKKOS_ENABLE_HIP) && defined (SACADO_VIEW_CUDA_HIERARCHICAL)
143  const size_type stride = Kokkos::ViewScalarStride<InputViewType1>::stride;
144  const bool use_team =
148  ( stride > 1 );
149 #elif defined (KOKKOS_ENABLE_HIP) && defined (SACADO_VIEW_CUDA_HIERARCHICAL_DFAD)
150  const size_type stride = team_policy_type::vector_length_max(); // 64
151  const bool use_team =
156 #else
157  const size_type stride = 1;
158  const bool use_team = false;
159 #endif
160 
161  if (use_team) {
162  const size_type team_size = 256 / stride;
163  team_policy_type policy( (nrow+team_size-1)/team_size, team_size, stride );
164  Kokkos::parallel_for( policy, MultiplyKernel(v1,v2,v3,update) );
165  }
166  else {
167  range_policy_type policy( 0, nrow );
168  Kokkos::parallel_for( policy, MultiplyKernel(v1,v2,v3,update) );
169  }
170  }
171 };
172 
173 // Kernel to assign a constant to a view
174 template <typename ViewType>
176  typedef typename ViewType::execution_space execution_space;
177  typedef typename ViewType::size_type size_type;
178  typedef typename ViewType::value_type::value_type ScalarType;
179  typedef Kokkos::TeamPolicy< execution_space> team_policy_type;
180  typedef Kokkos::RangePolicy< execution_space> range_policy_type;
181  typedef typename team_policy_type::member_type team_handle;
182  static const size_type stride = Kokkos::ViewScalarStride<ViewType>::stride;
183 
184  const ViewType m_v;
186 
187  ScalarAssignKernel(const ViewType& v, const ScalarType& s) :
188  m_v(v), m_s(s) {};
189 
190  // Multiply entries for row 'i' with a value
191  KOKKOS_INLINE_FUNCTION
192  void operator() (const size_type i) const {
193  m_v(i) = m_s;
194  }
195 
196  KOKKOS_INLINE_FUNCTION
197  void operator()( const team_handle& team ) const
198  {
199  const size_type i = team.league_rank()*team.team_size() + team.team_rank();
200  if (i < m_v.extent(0))
201  (*this)(i);
202  }
203 
204  // Kernel launch
205  static void apply(const ViewType& v, const ScalarType& s) {
206  const size_type nrow = v.extent(0);
207 
208 #if defined (KOKKOS_ENABLE_CUDA) && defined (SACADO_VIEW_CUDA_HIERARCHICAL)
209  const bool use_team =
213  ( stride > 1 );
214 #elif defined (KOKKOS_ENABLE_CUDA) && defined (SACADO_VIEW_CUDA_HIERARCHICAL_DFAD)
215  const bool use_team =
220 #elif defined (KOKKOS_ENABLE_HIP) && defined (SACADO_VIEW_CUDA_HIERARCHICAL)
221  const bool use_team =
225  ( stride > 1 );
226 #elif defined (KOKKOS_ENABLE_HIP) && defined (SACADO_VIEW_CUDA_HIERARCHICAL_DFAD)
227  const bool use_team =
232 #else
233  const bool use_team = false;
234 #endif
235 
236  if (use_team) {
237  const size_type team_size = 256 / stride;
238  team_policy_type policy( (nrow+team_size-1)/team_size, team_size, stride );
239  Kokkos::parallel_for( policy, ScalarAssignKernel(v,s) );
240  }
241  else {
242  range_policy_type policy( 0, nrow );
243  Kokkos::parallel_for( policy, ScalarAssignKernel(v,s) );
244  }
245  }
246 };
247 
248 // Kernel to assign a constant to a view
249 template <typename ViewType, typename ScalarViewType>
251  typedef typename ViewType::execution_space execution_space;
252  typedef typename ViewType::size_type size_type;
253  typedef typename ViewType::value_type ValueType;
254  typedef Kokkos::TeamPolicy< execution_space> team_policy_type;
255  typedef Kokkos::RangePolicy< execution_space> range_policy_type;
256  typedef typename team_policy_type::member_type team_handle;
257  typedef typename Kokkos::ThreadLocalScalarType<ViewType>::type local_scalar_type;
258  static const size_type stride = Kokkos::ViewScalarStride<ViewType>::stride;
259 
260  const ViewType m_v;
261  const ScalarViewType m_s;
262 
263  ValueAssignKernel(const ViewType& v, const ScalarViewType& s) :
264  m_v(v), m_s(s) {};
265 
266  // Multiply entries for row 'i' with a value
267  KOKKOS_INLINE_FUNCTION
268  void operator() (const size_type i) const {
269  local_scalar_type s = Sacado::partition_scalar<stride>(m_s());
270  m_v(i) = s;
271  }
272 
273  KOKKOS_INLINE_FUNCTION
274  void operator()( const team_handle& team ) const
275  {
276  const size_type i = team.league_rank()*team.team_size() + team.team_rank();
277  if (i < m_v.extent(0))
278  (*this)(i);
279  }
280 
281  // Kernel launch
282  static void apply(const ViewType& v, const ScalarViewType& s) {
283  const size_type nrow = v.extent(0);
284 
285 #if defined (KOKKOS_ENABLE_CUDA) && defined (SACADO_VIEW_CUDA_HIERARCHICAL)
286  const bool use_team =
290  ( stride > 1 );
291 #elif defined (KOKKOS_ENABLE_CUDA) && defined (SACADO_VIEW_CUDA_HIERARCHICAL_DFAD)
292  const bool use_team =
297 #elif defined (KOKKOS_ENABLE_HIP) && defined (SACADO_VIEW_CUDA_HIERARCHICAL)
298  const bool use_team =
302  ( stride > 1 );
303 #elif defined (KOKKOS_ENABLE_HIP) && defined (SACADO_VIEW_CUDA_HIERARCHICAL_DFAD)
304  const bool use_team =
309 #else
310  const bool use_team = false;
311 #endif
312 
313  if (use_team) {
314  const size_type team_size = 256 / stride;
315  team_policy_type policy( (nrow+team_size-1)/team_size, team_size, stride );
316  Kokkos::parallel_for( policy, ValueAssignKernel(v,s) );
317  }
318  else {
319  range_policy_type policy( 0, nrow );
320  Kokkos::parallel_for( policy, ValueAssignKernel(v,s) );
321  }
322  }
323 };
324 
325 // Kernel to assign a column of a rank-2 to a rank-1
326 template <typename InputViewType,
327  typename OutputViewType,
328  typename Enabled = void>
330  typedef typename InputViewType::execution_space execution_space;
331  typedef typename InputViewType::size_type size_type;
332  typedef Kokkos::TeamPolicy< execution_space> team_policy_type;
333  typedef Kokkos::RangePolicy< execution_space> range_policy_type;
334  typedef typename team_policy_type::member_type team_handle;
335  static const size_type stride = Kokkos::ViewScalarStride<InputViewType>::stride;
336 
337  const InputViewType m_v1;
338  const OutputViewType m_v2;
340 
341  AssignRank2Rank1Kernel(const InputViewType v1,
342  const OutputViewType v2,
343  const size_type col) :
344  m_v1(v1), m_v2(v2), m_col(col) {
345  static_assert( unsigned(InputViewType::rank) == 2 ,
346  "Require rank-2 input view" );
347  static_assert( unsigned(OutputViewType::rank) == 1 ,
348  "Require rank-1 output view" );
349  };
350 
351  // Multiply entries for row 'i' with a value
352  KOKKOS_INLINE_FUNCTION
353  void operator() (const size_type i) const {
354  m_v2(i) = m_v1(i,m_col);
355  }
356 
357  KOKKOS_INLINE_FUNCTION
358  void operator()( const team_handle& team ) const
359  {
360  const size_type i = team.league_rank()*team.team_size() + team.team_rank();
361  if (i < m_v1.extent(0))
362  (*this)(i);
363  }
364 
365  // Kernel launch
366  static void apply(const InputViewType v1,
367  const OutputViewType v2,
368  const size_type col) {
369  const size_type nrow = v1.extent(0);
370 
371 #if defined (KOKKOS_ENABLE_CUDA) && defined (SACADO_VIEW_CUDA_HIERARCHICAL)
372  const bool use_team =
376  ( stride > 1 );
377 #elif defined (KOKKOS_ENABLE_CUDA) && defined (SACADO_VIEW_CUDA_HIERARCHICAL_DFAD)
378  const bool use_team =
383 #elif defined (KOKKOS_ENABLE_HIP) && defined (SACADO_VIEW_CUDA_HIERARCHICAL)
384  const bool use_team =
388  ( stride > 1 );
389 #elif defined (KOKKOS_ENABLE_HIP) && defined (SACADO_VIEW_CUDA_HIERARCHICAL_DFAD)
390  const bool use_team =
395 #else
396  const bool use_team = false;
397 #endif
398 
399  if (use_team) {
400  const size_type team_size = 256 / stride;
401  team_policy_type policy( (nrow+team_size-1)/team_size, team_size, stride );
402  Kokkos::parallel_for( policy, AssignRank2Rank1Kernel(v1,v2,col) );
403  }
404  else {
405  range_policy_type policy( 0, nrow );
406  Kokkos::parallel_for( policy, AssignRank2Rank1Kernel(v1,v2,col) );
407  }
408  }
409 };
410 
411 // Kernel to test atomic_add
412 template <typename ViewType, typename ScalarViewType>
414  typedef typename ViewType::execution_space execution_space;
415  typedef typename ViewType::size_type size_type;
416  typedef Kokkos::TeamPolicy< execution_space> team_policy_type;
417  typedef Kokkos::RangePolicy< execution_space> range_policy_type;
418  typedef typename team_policy_type::member_type team_handle;
419  typedef typename Kokkos::ThreadLocalScalarType<ViewType>::type local_scalar_type;
420  static const size_type stride = Kokkos::ViewScalarStride<ViewType>::stride;
421 
422  const ViewType m_v;
423  const ScalarViewType m_s;
424 
425  AtomicAddKernel(const ViewType& v, const ScalarViewType& s) :
426  m_v(v), m_s(s) {};
427 
428  // Multiply entries for row 'i' with a value
429  KOKKOS_INLINE_FUNCTION
430  void operator() (const size_type i) const {
431  local_scalar_type x = m_v(i);
432  Kokkos::atomic_add(&(m_s()), x);
433  }
434 
435  KOKKOS_INLINE_FUNCTION
436  void operator()( const team_handle& team ) const
437  {
438  const size_type i = team.league_rank()*team.team_size() + team.team_rank();
439  if (i < m_v.extent(0))
440  (*this)(i);
441  }
442 
443  // Kernel launch
444  static void apply(const ViewType& v, const ScalarViewType& s) {
445  const size_type nrow = v.extent(0);
446 
447 #if defined (KOKKOS_ENABLE_CUDA) && defined (SACADO_VIEW_CUDA_HIERARCHICAL)
448  const bool use_team =
452  ( stride > 1 );
453 #elif defined (KOKKOS_ENABLE_CUDA) && defined (SACADO_VIEW_CUDA_HIERARCHICAL_DFAD)
454  const bool use_team =
459 #elif defined (KOKKOS_ENABLE_HIP) && defined (SACADO_VIEW_CUDA_HIERARCHICAL)
460  const bool use_team =
464  ( stride > 1 );
465 #elif defined (KOKKOS_ENABLE_HIP) && defined (SACADO_VIEW_CUDA_HIERARCHICAL_DFAD)
466  const bool use_team =
471 #else
472  const bool use_team = false;
473 #endif
474 
475  if (use_team) {
476  const size_type team_size = 256 / stride;
477  team_policy_type policy( (nrow+team_size-1)/team_size, team_size, stride );
478  Kokkos::parallel_for( policy, AtomicAddKernel(v,s) );
479  }
480  else {
481  range_policy_type policy( 0, nrow );
482  Kokkos::parallel_for( policy, AtomicAddKernel(v,s) );
483  }
484  }
485 };
486 
488  Kokkos_View_Fad, Size, FadType, Layout, Device )
489 {
490  typedef Kokkos::View<FadType*,Layout,Device> ViewType;
491  typedef typename ViewType::size_type size_type;
492 
493  const size_type num_rows = global_num_rows;
494 
495  // Create and fill view
496  ViewType v;
497 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
498  v = ViewType("view", num_rows);
499 #else
500  const size_type fad_size = global_fad_size;
501  v = ViewType("view", num_rows, fad_size+1);
502 #endif
503  TEUCHOS_TEST_EQUALITY(v.size(), num_rows, out, success);
504 }
505 
507  Kokkos_View_Fad, DeepCopy, FadType, Layout, Device )
508 {
509  typedef Kokkos::View<FadType**,Layout,Device> ViewType;
510  typedef typename ViewType::size_type size_type;
511  typedef typename ViewType::host_mirror_type host_view_type;
512 
513  const size_type num_rows = global_num_rows;
514  const size_type num_cols = global_num_cols;
515  const size_type fad_size = global_fad_size;
516 
517  // Create and fill view
518  ViewType v;
519 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
520  v = ViewType ("view", num_rows, num_cols);
521 #else
522  v = ViewType ("view", num_rows, num_cols, fad_size+1);
523 #endif
524  host_view_type h_v = Kokkos::create_mirror_view(v);
525  for (size_type i=0; i<num_rows; ++i)
526  for (size_type j=0; j<num_cols; ++j)
527  h_v(i,j) = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
528  Kokkos::deep_copy(v, h_v);
529 
530  // Copy back
531  host_view_type h_v2 = Kokkos::create_mirror_view(v);
532  Kokkos::deep_copy(h_v2, v);
533 
534  // Check
535  success = true;
536  for (size_type i=0; i<num_rows; ++i) {
537  for (size_type j=0; j<num_cols; ++j) {
538  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
539  success = success && checkFads(f, h_v2(i,j), out);
540  }
541  }
542 }
543 
545  Kokkos_View_Fad, DeepCopy_ConstantScalar, FadType, Layout, Device )
546 {
547  typedef Kokkos::View<FadType**,Layout,Device> ViewType;
548  typedef typename ViewType::size_type size_type;
549  typedef typename ViewType::host_mirror_type host_view_type;
550  typedef typename FadType::value_type value_type;
551 
552  const size_type num_rows = global_num_rows;
553  const size_type num_cols = global_num_cols;
554 
555  // Create and fill view
556  ViewType v;
557 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
558  v = ViewType ("view", num_rows, num_cols);
559 #else
560  const size_type fad_size = global_fad_size;
561  v = ViewType ("view", num_rows, num_cols, fad_size+1);
562 #endif
563  typename ViewType::array_type va = v;
564  Kokkos::deep_copy( va, 1.0 );
565 
566  // Deep copy a constant scalar
567  value_type a = 2.3456;
568  Kokkos::deep_copy( v, a );
569 
570  // Copy to host
571  host_view_type hv = Kokkos::create_mirror_view(v);
572  Kokkos::deep_copy(hv, v);
573 
574  // Check
575  success = true;
576  for (size_type i=0; i<num_rows; ++i) {
577  for (size_type j=0; j<num_cols; ++j) {
578 #if defined(HAVE_SACADO_VIEW_SPEC) && !defined(SACADO_DISABLE_FAD_VIEW_SPEC)
579  FadType f = FadType(fad_size, a);
580 #else
581  FadType f = a;
582 #endif
583  success = success && checkFads(f, hv(i,j), out);
584  }
585  }
586 }
587 
589  Kokkos_View_Fad, DeepCopy_ConstantZero, FadType, Layout, Device )
590 {
591  typedef Kokkos::View<FadType**,Layout,Device> ViewType;
592  typedef typename ViewType::size_type size_type;
593  typedef typename ViewType::host_mirror_type host_view_type;
594  typedef typename FadType::value_type value_type;
595 
596  const size_type num_rows = global_num_rows;
597  const size_type num_cols = global_num_cols;
598 
599  // Create and fill view
600  ViewType v;
601 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
602  v = ViewType ("view", num_rows, num_cols);
603 #else
604  const size_type fad_size = global_fad_size;
605  v = ViewType ("view", num_rows, num_cols, fad_size+1);
606 #endif
607  typename ViewType::array_type va = v;
608  Kokkos::deep_copy( va, 1.0 );
609 
610  // Deep copy a constant scalar
611  value_type a = 0.0;
612  Kokkos::deep_copy( v, a );
613 
614  // Copy to host
615  host_view_type hv = Kokkos::create_mirror_view(v);
616  Kokkos::deep_copy(hv, v);
617 
618  // Check
619  success = true;
620  for (size_type i=0; i<num_rows; ++i) {
621  for (size_type j=0; j<num_cols; ++j) {
622 #if defined(HAVE_SACADO_VIEW_SPEC) && !defined(SACADO_DISABLE_FAD_VIEW_SPEC)
623  FadType f = FadType(fad_size, a);
624 #else
625  FadType f = a;
626 #endif
627  success = success && checkFads(f, hv(i,j), out);
628  }
629  }
630 }
631 
633  Kokkos_View_Fad, DeepCopy_ConstantFad, FadType, Layout, Device )
634 {
635  typedef Kokkos::View<FadType**,Layout,Device> ViewType;
636  typedef typename ViewType::size_type size_type;
637  typedef typename ViewType::host_mirror_type host_view_type;
638 
639  const size_type num_rows = global_num_rows;
640  const size_type num_cols = global_num_cols;
641 
642  // Create and fill view
643  ViewType v;
644 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
645  v = ViewType ("view", num_rows, num_cols);
646 #else
647  const size_type fad_size = global_fad_size;
648  v = ViewType ("view", num_rows, num_cols, fad_size+1);
649 #endif
650  typename ViewType::array_type va = v;
651  Kokkos::deep_copy( va, 1.0 );
652 
653  // Deep copy a constant scalar
654  FadType a = 2.3456;
655  Kokkos::deep_copy( v, a );
656 
657  // Copy to host
658  host_view_type hv = Kokkos::create_mirror_view(v);
659  Kokkos::deep_copy(hv, v);
660 
661  // Check
662  success = true;
663  for (size_type i=0; i<num_rows; ++i) {
664  for (size_type j=0; j<num_cols; ++j) {
665 #if defined(HAVE_SACADO_VIEW_SPEC) && !defined(SACADO_DISABLE_FAD_VIEW_SPEC)
666  FadType f = FadType(fad_size, a.val());
667 #else
668  FadType f = a;
669 #endif
670  success = success && checkFads(f, hv(i,j), out);
671  }
672  }
673 }
674 
676  Kokkos_View_Fad, DeepCopy_ConstantFadFull, FadType, Layout, Device )
677 {
678  typedef Kokkos::View<FadType**,Layout,Device> ViewType;
679  typedef typename ViewType::size_type size_type;
680  typedef typename ViewType::host_mirror_type host_view_type;
681 
682  const size_type num_rows = global_num_rows;
683  const size_type num_cols = global_num_cols;
684  const size_type fad_size = global_fad_size;
685 
686  // Create and fill view
687  ViewType v;
688 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
689  v = ViewType ("view", num_rows, num_cols);
690 #else
691  v = ViewType ("view", num_rows, num_cols, fad_size+1);
692 #endif
693  typename ViewType::array_type va = v;
694  Kokkos::deep_copy( va, 1.0 );
695 
696  // Deep copy a constant Fad
697  FadType a(fad_size, 2.3456);
698  for (size_type i=0; i<fad_size; ++i)
699  a.fastAccessDx(i) = 7.89 + (i+1);
700 
701  // Copy to host
702  host_view_type hv = Kokkos::create_mirror_view(v);
703  Kokkos::deep_copy(hv, a);
704 
705  // Check
706  success = true;
707  for (size_type i=0; i<num_rows; ++i) {
708  for (size_type j=0; j<num_cols; ++j) {
709  success = success && checkFads(a, hv(i,j), out);
710  }
711  }
712 }
713 
715  Kokkos_View_Fad, LocalDeepCopy, FadType, Layout, Device )
716 {
717  typedef Kokkos::View<FadType***,Layout,Device> ViewType;
718  typedef Kokkos::View<FadType,Layout,Device> ScalarViewType;
719  typedef typename ViewType::size_type size_type;
720  typedef typename ViewType::host_mirror_type host_view_type;
721  typedef typename ScalarViewType::host_mirror_type host_scalar_view_type;
722 
723  const size_type num_rows = global_num_rows;
724  const size_type num_cols = global_num_cols;
725  const size_type num_slices = 10;
726  const size_type fad_size = global_fad_size;
727 
728  // Create and fill view
729  ViewType v;
730 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
731  v = ViewType ("view", num_rows, num_cols, num_slices);
732 #else
733  v = ViewType ("view", num_rows, num_cols, num_slices, fad_size+1);
734 #endif
735  typename ViewType::array_type va = v;
736  Kokkos::deep_copy( va, 1.0 );
737 
738  // Deep copy a constant Fad to the device
739  // Can't deep_copy directly because that doesn't work with DFad
740  FadType a(fad_size, 2.3456);
741  for (size_type i=0; i<fad_size; ++i)
742  a.fastAccessDx(i) = 7.89 + (i+1);
743 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
744  ScalarViewType a_view("a");
745 #else
746  ScalarViewType a_view("a", fad_size+1);
747 #endif
748  host_scalar_view_type ha_view = Kokkos::create_mirror_view(a_view);
749  ha_view() = a;
750  Kokkos::deep_copy( a_view, ha_view );
751 
752  // Excersize local_deep_copy by setting each row of s to a
753  Kokkos::parallel_for(Kokkos::RangePolicy<Device>(0,num_rows),
754  KOKKOS_LAMBDA(const int i)
755  {
756  auto s = Kokkos::subview(v,i,Kokkos::ALL,Kokkos::ALL);
757  Kokkos::Experimental::local_deep_copy(s,a_view());
758  });
759 
760  // Copy back to host
761  host_view_type hv = Kokkos::create_mirror_view(v);
762  Kokkos::deep_copy(hv, a);
763 
764  // Check
765  success = true;
766  for (size_type i=0; i<num_rows; ++i) {
767  for (size_type j=0; j<num_cols; ++j) {
768  for (size_type k=0; k<num_slices; ++k) {
769  success = success && checkFads(a, hv(i,j,k), out);
770  }
771  }
772  }
773 }
774 
776  Kokkos_View_Fad, LocalDeepCopyTeam, FadType, Layout, Device )
777 {
778  typedef Kokkos::View<FadType***,Layout,Device> ViewType;
779  typedef Kokkos::View<FadType,Layout,Device> ScalarViewType;
780  typedef typename ViewType::size_type size_type;
781  typedef typename ViewType::host_mirror_type host_view_type;
782  typedef typename ScalarViewType::host_mirror_type host_scalar_view_type;
783 
784  const size_type num_rows = global_num_rows;
785  const size_type num_cols = global_num_cols;
786  const size_type num_slices = 10;
787  const size_type fad_size = global_fad_size;
788 
789  // Create and fill view
790  ViewType v;
791 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
792  v = ViewType ("view", num_rows, num_cols, num_slices);
793 #else
794  v = ViewType ("view", num_rows, num_cols, num_slices, fad_size+1);
795 #endif
796  typename ViewType::array_type va = v;
797  Kokkos::deep_copy( va, 1.0 );
798 
799  // Deep copy a constant Fad to the device
800  // Can't deep_copy directly because that doesn't work with DFad
801  FadType a(fad_size, 2.3456);
802  for (size_type i=0; i<fad_size; ++i)
803  a.fastAccessDx(i) = 7.89 + (i+1);
804 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
805  ScalarViewType a_view("a");
806 #else
807  ScalarViewType a_view("a", fad_size+1);
808 #endif
809  host_scalar_view_type ha_view = Kokkos::create_mirror_view(a_view);
810  ha_view() = a;
811  Kokkos::deep_copy( a_view, ha_view );
812 
813  // Excersize local_deep_copy by setting each row of s to a
814  typedef Kokkos::TeamPolicy<Device> Policy;
815  static const size_type stride = Kokkos::ViewScalarStride<ViewType>::stride;
816  Kokkos::parallel_for(Policy(num_rows,Kokkos::AUTO,stride),
817  KOKKOS_LAMBDA(const typename Policy::member_type& team)
818  {
819  int i = team.league_rank();
820  auto s = Kokkos::subview(v,i,Kokkos::ALL,Kokkos::ALL);
821  Kokkos::Experimental::local_deep_copy(team,s,a_view());
822  });
823 
824  // Copy back to host
825  host_view_type hv = Kokkos::create_mirror_view(v);
826  Kokkos::deep_copy(hv, a);
827 
828  // Check
829  success = true;
830  for (size_type i=0; i<num_rows; ++i) {
831  for (size_type j=0; j<num_cols; ++j) {
832  for (size_type k=0; k<num_slices; ++k) {
833  success = success && checkFads(a, hv(i,j,k), out);
834  }
835  }
836  }
837 }
838 
840  Kokkos_View_Fad, ScalarAssign, FadType, Layout, Device )
841 {
842  typedef Kokkos::View<FadType*,Layout,Device> ViewType;
843  typedef typename ViewType::size_type size_type;
844  typedef typename ViewType::host_mirror_type host_view_type;
845  typedef typename FadType::value_type value_type;
846 
847  const size_type num_rows = global_num_rows;
848 
849  // Create and fill view
850  ViewType v;
851 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
852  v = ViewType ("view", num_rows);
853 #else
854  const size_type fad_size = global_fad_size;
855  v = ViewType ("view", num_rows, fad_size+1);
856 #endif
857  typename ViewType::array_type va = v;
858  Kokkos::deep_copy( va, 1.0 );
859 
860  // Deep copy a constant scalar
861  value_type a = 2.3456;
863 
864  // Copy to host
865  host_view_type hv = Kokkos::create_mirror_view(v);
866  Kokkos::deep_copy(hv, v);
867 
868  // Check
869  success = true;
870  for (size_type i=0; i<num_rows; ++i) {
871 #if defined(HAVE_SACADO_VIEW_SPEC) && !defined(SACADO_DISABLE_FAD_VIEW_SPEC)
872  FadType f = FadType(fad_size, a);
873 #else
874  FadType f = a;
875 #endif
876  success = success && checkFads(f, hv(i), out);
877  }
878 }
879 
881  Kokkos_View_Fad, ValueAssign, FadType, Layout, Device )
882 {
883  typedef Kokkos::View<FadType*,Layout,Device> ViewType;
884  typedef Kokkos::View<FadType,Layout,Device> ScalarViewType;
885  typedef typename ViewType::size_type size_type;
886  typedef typename ViewType::host_mirror_type host_view_type;
887  typedef typename ScalarViewType::host_mirror_type host_scalar_view_type;
888 
889  const size_type num_rows = global_num_rows;
890  const size_type fad_size = global_fad_size;
891 
892  // Create and fill view
893  ViewType v;
894  ScalarViewType a;
895 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
896  v = ViewType ("view", num_rows);
897  a = ScalarViewType ("fad");
898 #else
899  v = ViewType ("view", num_rows, fad_size+1);
900  a = ScalarViewType ("fad", fad_size+1);
901 #endif
902  typename ViewType::array_type va = v;
903  Kokkos::deep_copy( va, 1.0 );
904 
905  // Deep copy a constant scalar
906  Kokkos::deep_copy(a, 2.3456);
907 
908  Kokkos::parallel_for(Kokkos::RangePolicy< Device>(0, fad_size), KOKKOS_LAMBDA(const int i) {
909  a().fastAccessDx(i) = 7.89 + i;
910  });
911  Kokkos::fence();
912 
914  Kokkos::fence();
915 
916  // Copy to host
917  host_view_type hv = Kokkos::create_mirror_view(v);
918  Kokkos::deep_copy(hv, v);
919 
920  host_scalar_view_type ha = Kokkos::create_mirror_view(a);
921  Kokkos::deep_copy(ha, a);
922 
923  // Check
924  success = true;
925  for (size_type i=0; i<num_rows; ++i) {
926  success = success && checkFads(ha(), hv(i), out);
927  }
928 }
929 
931  Kokkos_View_Fad, Resize, FadType, Layout, Device )
932 {
933  typedef Kokkos::View<FadType**,Layout,Device> ViewType;
934  typedef typename ViewType::size_type size_type;
935  typedef typename ViewType::host_mirror_type host_view_type;
936 
937  const size_type num_rows = global_num_rows;
938  const size_type num_cols = global_num_cols;
939  const size_type fad_size = global_fad_size;
940 
941  // Create and fill view
942  ViewType v;
943 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
944  v = ViewType ("view", num_rows, num_cols);
945 #else
946  v = ViewType ("view", num_rows, num_cols, fad_size+1);
947 #endif
948  host_view_type h_v = Kokkos::create_mirror_view(v);
949  for (size_type i=0; i<num_rows; ++i)
950  for (size_type j=0; j<num_cols; ++j)
951  h_v(i,j) = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
952  Kokkos::deep_copy(v, h_v);
953 
954  // Resize
955 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
956  Kokkos::resize(v, num_rows, num_cols+1);
957 #else
958  Kokkos::resize(v, num_rows, num_cols+1, fad_size+1);
959 #endif
960 
961  // Copy back
962  host_view_type h_v2 = Kokkos::create_mirror_view(v);
963  Kokkos::deep_copy(h_v2, v);
964 
965  // Check
966  success = true;
967  for (size_type i=0; i<num_rows; ++i) {
968  for (size_type j=0; j<num_cols; ++j) {
969  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
970  success = success && checkFads(f, h_v2(i,j), out);
971  }
972 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
973  FadType f = 0.0;
974 #else
975  FadType f(fad_size, 0.0);
976 #endif
977  success = success && checkFads(f, h_v2(i,num_cols), out);
978  }
979 }
980 
982  Kokkos_View_Fad, Multiply, FadType, Layout, Device )
983 {
984  typedef Kokkos::View<FadType*,Layout,Device> ViewType;
985  typedef typename ViewType::size_type size_type;
986  typedef typename ViewType::host_mirror_type host_view_type;
987 
988  const size_type num_rows = global_num_rows;
989  const size_type fad_size = global_fad_size;
990 
991  // Create and fill views
992  ViewType v1, v2;
993 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
994  v1 = ViewType ("view1", num_rows);
995  v2 = ViewType ("view2", num_rows);
996 #else
997  v1 = ViewType ("view1", num_rows, fad_size+1);
998  v2 = ViewType ("view2", num_rows, fad_size+1);
999 #endif
1000  host_view_type h_v1 = Kokkos::create_mirror_view(v1);
1001  host_view_type h_v2 = Kokkos::create_mirror_view(v2);
1002  for (size_type i=0; i<num_rows; ++i) {
1003  h_v1(i) = generate_fad<FadType>(
1004  num_rows, size_type(2), fad_size, i, size_type(0));
1005  h_v2(i) = generate_fad<FadType>(
1006  num_rows, size_type(2), fad_size, i, size_type(1));
1007  }
1008  Kokkos::deep_copy(v1, h_v1);
1009  Kokkos::deep_copy(v2, h_v2);
1010 
1011  // Launch kernel
1012  ViewType v3;
1013 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1014  v3 = ViewType ("view3", num_rows);
1015 #else
1016  v3 = ViewType ("view3", num_rows, fad_size+1);
1017 #endif
1019 
1020  // Copy back
1021  host_view_type h_v3 = Kokkos::create_mirror_view(v3);
1022  Kokkos::deep_copy(h_v3, v3);
1023 
1024  // Check
1025  success = true;
1026  for (size_type i=0; i<num_rows; ++i) {
1027  FadType f1 =
1028  generate_fad<FadType>(num_rows, size_type(2), fad_size, i, size_type(0));
1029  FadType f2 =
1030  generate_fad<FadType>(num_rows, size_type(2), fad_size, i, size_type(1));
1031  FadType f3 = f1*f2;
1032  success = success && checkFads(f3, h_v3(i), out);
1033  }
1034 }
1035 
1037  Kokkos_View_Fad, MultiplyUpdate, FadType, Layout, Device )
1038 {
1039  typedef Kokkos::View<FadType*,Layout,Device> ViewType;
1040  typedef typename ViewType::size_type size_type;
1041  typedef typename ViewType::host_mirror_type host_view_type;
1042 
1043  const size_type num_rows = global_num_rows;
1044  const size_type fad_size = global_fad_size;
1045 
1046  // Create and fill views
1047  ViewType v1, v2;
1048 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1049  v1 = ViewType ("view1", num_rows);
1050  v2 = ViewType ("view2", num_rows);
1051 #else
1052  v1 = ViewType ("view1", num_rows, fad_size+1);
1053  v2 = ViewType ("view2", num_rows, fad_size+1);
1054 #endif
1055  host_view_type h_v1 = Kokkos::create_mirror_view(v1);
1056  host_view_type h_v2 = Kokkos::create_mirror_view(v2);
1057  for (size_type i=0; i<num_rows; ++i) {
1058  h_v1(i) = generate_fad<FadType>(
1059  num_rows, size_type(2), fad_size, i, size_type(0));
1060  h_v2(i) = generate_fad<FadType>(
1061  num_rows, size_type(2), fad_size, i, size_type(1));
1062  }
1063  Kokkos::deep_copy(v1, h_v1);
1064  Kokkos::deep_copy(v2, h_v2);
1065 
1066  // Launch kernel
1067  ViewType v3;
1068 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1069  v3 = ViewType ("view3", num_rows);
1070 #else
1071  v3 = ViewType ("view3", num_rows, fad_size+1);
1072 #endif
1073  Kokkos::deep_copy(v3, 1.0);
1074  MultiplyKernel<ViewType>::apply(v1,v2,v3,true);
1075 
1076  // Copy back
1077  host_view_type h_v3 = Kokkos::create_mirror_view(v3);
1078  Kokkos::deep_copy(h_v3, v3);
1079 
1080  // Check
1081  success = true;
1082  for (size_type i=0; i<num_rows; ++i) {
1083  FadType f1 =
1084  generate_fad<FadType>(num_rows, size_type(2), fad_size, i, size_type(0));
1085  FadType f2 =
1086  generate_fad<FadType>(num_rows, size_type(2), fad_size, i, size_type(1));
1087  FadType f3 = 1.0 + f1*f2;
1088  success = success && checkFads(f3, h_v3(i), out);
1089  }
1090 }
1091 
1093  Kokkos_View_Fad, MultiplyConst, FadType, Layout, Device )
1094 {
1095  typedef Kokkos::View<const FadType*,Layout,Device,Kokkos::MemoryUnmanaged> ConstViewType;
1096  typedef Kokkos::View<FadType*,Layout,Device> ViewType;
1097  typedef typename ViewType::size_type size_type;
1098  typedef typename ViewType::host_mirror_type host_view_type;
1099 
1100  const size_type num_rows = global_num_rows;
1101  const size_type fad_size = global_fad_size;
1102 
1103  // Create and fill views
1104  ViewType v1, v2;
1105 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1106  v1 = ViewType ("view1", num_rows);
1107  v2 = ViewType ("view2", num_rows);
1108 #else
1109  v1 = ViewType ("view1", num_rows, fad_size+1);
1110  v2 = ViewType ("view2", num_rows, fad_size+1);
1111 #endif
1112  host_view_type h_v1 = Kokkos::create_mirror_view(v1);
1113  host_view_type h_v2 = Kokkos::create_mirror_view(v2);
1114  for (size_type i=0; i<num_rows; ++i) {
1115  h_v1(i) = generate_fad<FadType>(
1116  num_rows, size_type(2), fad_size, i, size_type(0));
1117  h_v2(i) = generate_fad<FadType>(
1118  num_rows, size_type(2), fad_size, i, size_type(1));
1119  }
1120  Kokkos::deep_copy(v1, h_v1);
1121  Kokkos::deep_copy(v2, h_v2);
1122 
1123  ConstViewType cv1 = v1;
1124 
1125  // Launch kernel
1126  ViewType v3;
1127 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1128  v3 = ViewType ("view3", num_rows);
1129 #else
1130  v3 = ViewType ("view3", num_rows, fad_size+1);
1131 #endif
1133 
1134  // Copy back
1135  host_view_type h_v3 = Kokkos::create_mirror_view(v3);
1136  Kokkos::deep_copy(h_v3, v3);
1137 
1138  // Check
1139  success = true;
1140  for (size_type i=0; i<num_rows; ++i) {
1141  FadType f1 =
1142  generate_fad<FadType>(num_rows, size_type(2), fad_size, i, size_type(0));
1143  FadType f2 =
1144  generate_fad<FadType>(num_rows, size_type(2), fad_size, i, size_type(1));
1145  FadType f3 = f1*f2;
1146  success = success && checkFads(f3, h_v3(i), out);
1147  }
1148 }
1149 
1151  Kokkos_View_Fad, MultiplyMixed, FadType, Layout, Device )
1152 {
1153  typedef Kokkos::View<FadType*,Layout,Device> ViewType;
1154  typedef typename ViewType::size_type size_type;
1155  typedef typename ViewType::host_mirror_type host_view_type;
1156 
1157  const size_type num_rows = 2;
1158  const size_type fad_size = global_fad_size;
1159 
1160  // Create and fill views -- do everything on the host for this test
1161  FadType f0 = generate_fad<FadType>(
1162  num_rows, size_type(2), fad_size, size_type(0), size_type(0));
1163  FadType f1 = generate_fad<FadType>(
1164  num_rows, size_type(2), fad_size, size_type(1), size_type(0));
1165  host_view_type h_v;
1166 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1167  h_v = host_view_type ("view1", num_rows);
1168 #else
1169  h_v = host_view_type ("view1", num_rows, fad_size+1);
1170 #endif
1171  h_v(0) = f0;
1172  h_v(1) = f1;
1173 
1174  FadType f2 = f0 * h_v(1);
1175 
1176  // Check
1177  FadType f3 = f0 * f1;
1178  success = checkFads(f3, f2, out);
1179 }
1180 
1182  Kokkos_View_Fad, AtomicAdd, FadType, Layout, Device )
1183 {
1184  typedef Kokkos::View<FadType*,Layout,Device> ViewType;
1185  typedef Kokkos::View<FadType,Layout,Device> ScalarViewType;
1186  typedef typename ViewType::size_type size_type;
1187  typedef typename ScalarViewType::host_mirror_type host_scalar_view_type;
1188 
1189  const size_type num_rows = global_num_rows;
1190  const size_type fad_size = global_fad_size;
1191 
1192  // Create and fill view
1193  ViewType v;
1194 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1195  v = ViewType ("view", num_rows);
1196 #else
1197  v = ViewType ("view", num_rows, fad_size+1);
1198 #endif
1199  Kokkos::deep_copy(v, 2.3456);
1200 
1201  Kokkos::parallel_for(Kokkos::RangePolicy<Device>(0, num_rows), KOKKOS_LAMBDA(const size_type i) {
1202  for (size_type j = 0; j < fad_size; ++j)
1203  v(i).fastAccessDx(j) = 7.89 + j;
1204  });
1205 
1206  // Create scalar view
1207  ScalarViewType s;
1208 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1209  s = ScalarViewType ("scalar view");
1210 #else
1211  s = ScalarViewType ("scalar view", fad_size+1);
1212 #endif
1213 
1214  // Call atomic_add kernel, which adds up entries in v
1216 
1217  // Copy to host
1218  host_scalar_view_type hs = Kokkos::create_mirror_view(s);
1219  Kokkos::deep_copy(hs, s);
1220 
1221  // Check
1222  auto hv = Kokkos::create_mirror_view(v);
1223  Kokkos::deep_copy(hv, v);
1224 
1225  FadType b = num_rows*hv(0);
1226  success = checkFads(b, hs(), out);
1227 }
1228 
1230  Kokkos_View_Fad, Rank8, FadType, Layout, Device )
1231 {
1232  typedef Kokkos::View<FadType*******,Layout,Device> ViewType;
1233  typedef typename ViewType::size_type size_type;
1234  typedef typename ViewType::host_mirror_type host_view_type;
1235 
1236  const size_type fad_size = global_fad_size;
1237 
1238  // Create and fill view
1239  ViewType v;
1240 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1241  v = ViewType ("view", 100, 1, 2, 3, 4, 5, 6);
1242 #else
1243  v = ViewType ("view", 100, 1, 2, 3, 4, 5, 6, fad_size+1);
1244 #endif
1245  host_view_type h_v = Kokkos::create_mirror_view(v);
1246  typename host_view_type::array_type h_a = h_v;
1247  Kokkos::deep_copy(h_a, 1.0);
1248 
1249  FadType f1 = FadType(fad_size, 2.0);
1250  h_v(99,0,1,2,3,4,5) = f1;
1251  FadType f2 = h_v(99,0,1,2,3,4,5);
1252 
1253  // Check
1254  success = checkFads(f1, f2, out);
1255 }
1256 
1258  Kokkos_View_Fad, Roger, FadType, Layout, Device )
1259 {
1260  Kokkos::View<FadType*,Layout,Device> a;
1261  Kokkos::View<FadType**,Layout,Device> b;
1262  Kokkos::View<FadType***,Layout,Device> c;
1263  Kokkos::View<FadType****,Layout,Device> d;
1264  Kokkos::View<FadType*****,Layout,Device> e;
1265  Kokkos::View<FadType******,Layout,Device> f;
1266  Kokkos::View<FadType*******,Layout,Device> g;
1267 
1268 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1269  a = Kokkos::View<FadType*,Layout,Device>("a",4);
1270  b = Kokkos::View<FadType**,Layout,Device> ("b",4,4);
1271  c = Kokkos::View<FadType***,Layout,Device> ("c",4,4,4);
1272  d = Kokkos::View<FadType****,Layout,Device> ("d",4,4,4,4);
1273  e = Kokkos::View<FadType*****,Layout,Device> ("e",4,4,4,4,4);
1274  f = Kokkos::View<FadType******,Layout,Device> ("f",4,4,4,4,4,4);
1275  g = Kokkos::View<FadType*******,Layout,Device> ("g",4,4,4,4,4,4,4);
1276 #else
1277  const unsigned fad_size = global_fad_size;
1278  a = Kokkos::View<FadType*,Layout,Device>("a",4,fad_size+1);
1279  b = Kokkos::View<FadType**,Layout,Device> ("b",4,4,fad_size+1);
1280  c = Kokkos::View<FadType***,Layout,Device> ("c",4,4,4,fad_size+1);
1281  d = Kokkos::View<FadType****,Layout,Device> ("d",4,4,4,4,fad_size+1);
1282  e = Kokkos::View<FadType*****,Layout,Device> ("e",4,4,4,4,4,fad_size+1);
1283  f = Kokkos::View<FadType******,Layout,Device> ("f",4,4,4,4,4,4,fad_size+1);
1284  g = Kokkos::View<FadType*******,Layout,Device> ("g",4,4,4,4,4,4,4,fad_size+1);
1285 #endif
1286 
1287  typedef typename Device::memory_space memory_space;
1288  const bool is_accessible =
1289  Kokkos::Impl::MemorySpaceAccess<Kokkos::HostSpace,
1290  memory_space>::accessible;
1291  if (is_accessible) {
1292  a(0) = FadType(1.0);
1293  f(0,0,0,0,0,0) = FadType(1.0);
1294  g(0,0,0,0,0,0,0) = FadType(1.0);
1295  }
1296 
1297  // Check
1298  success = true;
1299 }
1300 
1302  Kokkos_View_Fad, AssignDifferentStrides, FadType, Layout, Device )
1303 {
1304  typedef Kokkos::View<FadType**,Layout,Device> ViewType1;
1305  typedef Kokkos::View<FadType*,Layout,Device> ViewType2;
1306  typedef typename ViewType1::size_type size_type;
1307  typedef typename ViewType1::host_mirror_type host_view_type1;
1308  typedef typename ViewType2::host_mirror_type host_view_type2;
1309 
1310  const size_type num_rows = global_num_rows;
1311  const size_type num_cols = global_num_cols;
1312  const size_type fad_size = global_fad_size;
1313 
1314  // Create and fill views
1315  ViewType1 v1;
1316 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1317  v1 = ViewType1 ("view1", num_rows, num_cols);
1318 #else
1319  v1 = ViewType1 ("view1", num_rows, num_cols, fad_size+1);
1320 #endif
1321  host_view_type1 h_v1 = Kokkos::create_mirror_view(v1);
1322  for (size_type i=0; i<num_rows; ++i) {
1323  for (size_type j=0; j<num_cols; ++j) {
1324  h_v1(i,j) = generate_fad<FadType>(
1325  num_rows, num_cols, fad_size, i, j);
1326  }
1327  }
1328  Kokkos::deep_copy(v1, h_v1);
1329 
1330  // Launch kernel
1331  ViewType2 v2;
1332 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1333  v2 = ViewType2 ("view2", num_rows);
1334 #else
1335  v2 = ViewType2 ("view2", num_rows, fad_size+1);
1336 #endif
1338 
1339  // Copy back
1340  host_view_type2 h_v2 = Kokkos::create_mirror_view(v2);
1341  Kokkos::deep_copy(h_v2, v2);
1342 
1343  // Check
1344  success = true;
1345  for (size_type i=0; i<num_rows; ++i) {
1346  FadType f =
1347  generate_fad<FadType>(num_rows, num_cols, fad_size, i, size_type(1));
1348  success = success && checkFads(f, h_v2(i), out);
1349  }
1350 }
1351 
1353  Kokkos_View_Fad, ScalarValue, FadType, Layout, Device )
1354 {
1355  typedef typename Sacado::ScalarType<FadType>::type ScalarType;
1356  typedef Kokkos::View<FadType,Layout,Device> ViewType1;
1357  typedef Kokkos::View<ScalarType,Layout,Device> ViewType2;
1358  typedef typename ViewType1::size_type size_type;
1359  typedef typename ViewType1::host_mirror_type host_view_type1;
1360  typedef typename ViewType2::host_mirror_type host_view_type2;
1361 
1362  const int fad_size = global_fad_size;
1363 
1364  // Create and fill views
1365  ViewType1 v1;
1366 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1367  v1 = ViewType1 ("view1");
1368 #else
1369  v1 = ViewType1 ("view1", fad_size+1);
1370 #endif
1371  host_view_type1 h_v1 = Kokkos::create_mirror_view(v1);
1372  h_v1() = generate_fad<FadType>(1, 1, fad_size, 0, 0);
1373  Kokkos::deep_copy(v1, h_v1);
1374 
1375  // Launch kernel
1376  ViewType2 v2 = ViewType2 ("view2");
1377  Kokkos::parallel_for(Kokkos::RangePolicy<Device>(0,1),
1378  KOKKOS_LAMBDA(const size_type i)
1379  {
1380  v2() = Sacado::scalarValue(v1());
1381  });
1382 
1383  // Copy back
1384  host_view_type2 h_v2 = Kokkos::create_mirror_view(v2);
1385  Kokkos::deep_copy(h_v2, v2);
1386 
1387  // Check
1388  success = true;
1389  TEUCHOS_TEST_EQUALITY(h_v1().val(), h_v2(), out, success);
1390 }
1391 
1392 #if defined(HAVE_SACADO_KOKKOS) && defined(HAVE_SACADO_VIEW_SPEC) && !defined(SACADO_DISABLE_FAD_VIEW_SPEC)
1393 
1395  Kokkos_View_Fad, DynRankDimensionScalar, FadType, Layout, Device )
1396 {
1397  typedef Kokkos::DynRankView<double,Layout,Device> DoubleViewType;
1398  typedef Kokkos::DynRankView<FadType,Layout,Device> FadViewType;
1399  typedef typename FadViewType::size_type size_type;
1400 
1401  const size_type num_rows = global_num_rows;
1402  const size_type fad_size = global_fad_size;
1403 
1404  // Create views
1405  DoubleViewType v1("view1", num_rows);
1406  FadViewType v2 ("view2", num_rows, fad_size+1);
1407 
1408  // Check dimension scalar works
1409  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(v1), 0, out, success);
1410  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(v2), fad_size+1, out, success);
1411 }
1412 
1414  Kokkos_View_Fad, DynRankAssignStatic0, FadType, Layout, Device )
1415 {
1416  typedef Kokkos::View<FadType,Layout,Device> StaticViewType;
1417  typedef Kokkos::DynRankView<FadType,Layout,Device> DynamicViewType;
1418  typedef typename StaticViewType::size_type size_type;
1419 
1420  const size_type num_rows = global_num_rows;
1421  const size_type num_cols = global_num_cols;
1422  const size_type fad_size = global_fad_size;
1423 
1424  // Create and fill views
1425  StaticViewType v1("view", fad_size+1);
1426  auto h_v1 = Kokkos::create_mirror_view(v1);
1427  h_v1() = generate_fad<FadType>(num_rows, num_cols, fad_size, size_type(0), size_type(0));
1428  Kokkos::deep_copy(v1, h_v1);
1429 
1430  // Assign static to dynamic
1431  DynamicViewType v2 = v1;
1432 
1433  // Copy back
1434  auto h_v2 = Kokkos::create_mirror_view(v2);
1435  Kokkos::deep_copy(h_v2, v2);
1436 
1437  // Check dimensions are correct
1438  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(v2), fad_size+1, out, success);
1439  TEUCHOS_TEST_EQUALITY(v2.stride(0), v1.stride(0), out, success);
1440 
1441  // Check values
1442  FadType f =
1443  generate_fad<FadType>(num_rows, num_cols, fad_size, size_type(0), size_type(0));
1444  success = success && checkFads(f, h_v2(), out);
1445 }
1446 
1448  Kokkos_View_Fad, DynRankAssignStatic1, FadType, Layout, Device )
1449 {
1450  typedef Kokkos::View<FadType*,Layout,Device> StaticViewType;
1451  typedef Kokkos::DynRankView<FadType,Layout,Device> DynamicViewType;
1452  typedef typename StaticViewType::size_type size_type;
1453 
1454  const size_type num_rows = global_num_rows;
1455  const size_type num_cols = global_num_cols;
1456  const size_type fad_size = global_fad_size;
1457 
1458  // Create and fill views
1459  StaticViewType v1("view", num_rows, fad_size+1);
1460  auto h_v1 = Kokkos::create_mirror_view(v1);
1461  for (size_type i=0; i<num_rows; ++i)
1462  h_v1(i) =
1463  generate_fad<FadType>(num_rows, num_cols, fad_size, i, size_type(0));
1464  Kokkos::deep_copy(v1, h_v1);
1465 
1466  // Assign static to dynamic
1467  DynamicViewType v2 = v1;
1468 
1469  // Copy back
1470  auto h_v2 = Kokkos::create_mirror_view(v2);
1471  Kokkos::deep_copy(h_v2, v2);
1472 
1473  // Check dimensions are correct
1474  TEUCHOS_TEST_EQUALITY(v2.extent(0), num_rows, out, success);
1475  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(v2), fad_size+1, out, success);
1476  TEUCHOS_TEST_EQUALITY(v2.stride(0), v1.stride(0), out, success);
1477  TEUCHOS_TEST_EQUALITY(v2.stride(1), v1.stride(1), out, success);
1478 
1479  // Check values
1480  for (size_type i=0; i<num_rows; ++i) {
1481  FadType f =
1482  generate_fad<FadType>(num_rows, num_cols, fad_size, i, size_type(0));
1483  success = success && checkFads(f, h_v2(i), out);
1484  }
1485 }
1486 
1488  Kokkos_View_Fad, DynRankAssignStatic2, FadType, Layout, Device )
1489 {
1490  typedef Kokkos::View<FadType**,Layout,Device> StaticViewType;
1491  typedef Kokkos::DynRankView<FadType,Layout,Device> DynamicViewType;
1492  typedef typename StaticViewType::size_type size_type;
1493 
1494  const size_type num_rows = global_num_rows;
1495  const size_type num_cols = global_num_cols;
1496  const size_type fad_size = global_fad_size;
1497 
1498  // Create and fill views
1499  StaticViewType v1("view", num_rows, num_cols, fad_size+1);
1500  auto h_v1 = Kokkos::create_mirror_view(v1);
1501  for (size_type i=0; i<num_rows; ++i)
1502  for (size_type j=0; j<num_cols; ++j)
1503  h_v1(i,j) = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
1504  Kokkos::deep_copy(v1, h_v1);
1505 
1506  // Assign static to dynamic
1507  DynamicViewType v2 = v1;
1508 
1509  // Copy back
1510  auto h_v2 = Kokkos::create_mirror_view(v2);
1511  Kokkos::deep_copy(h_v2, v2);
1512 
1513  // Check dimensions are correct
1514  TEUCHOS_TEST_EQUALITY(v2.extent(0), num_rows, out, success);
1515  TEUCHOS_TEST_EQUALITY(v2.extent(1), num_cols, out, success);
1516  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(v2), fad_size+1, out, success);
1517  TEUCHOS_TEST_EQUALITY(v2.stride(0), v1.stride(0), out, success);
1518  TEUCHOS_TEST_EQUALITY(v2.stride(1), v1.stride(1), out, success);
1519  TEUCHOS_TEST_EQUALITY(v2.stride(2), v1.stride(2), out, success);
1520 
1521  // Check values
1522  for (size_type i=0; i<num_rows; ++i) {
1523  for (size_type j=0; j<num_cols; ++j) {
1524  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
1525  success = success && checkFads(f, h_v2(i,j), out);
1526  }
1527  }
1528 }
1529 
1531  Kokkos_View_Fad, DynRankMultiply, FadType, Layout, Device )
1532 {
1533  typedef Kokkos::DynRankView<FadType,Layout,Device> ViewType;
1534  typedef typename ViewType::size_type size_type;
1535  typedef typename ViewType::host_mirror_type host_view_type;
1536 
1537  const size_type num_rows = global_num_rows;
1538  const size_type fad_size = global_fad_size;
1539 
1540  // Create and fill views
1541  ViewType v1("view1", num_rows, fad_size+1);
1542  ViewType v2("view2", num_rows, fad_size+1);
1543  host_view_type h_v1 = Kokkos::create_mirror_view(v1);
1544  host_view_type h_v2 = Kokkos::create_mirror_view(v2);
1545  for (size_type i=0; i<num_rows; ++i) {
1546  h_v1(i) = generate_fad<FadType>(
1547  num_rows, size_type(2), fad_size, i, size_type(0));
1548  h_v2(i) = generate_fad<FadType>(
1549  num_rows, size_type(2), fad_size, i, size_type(1));
1550  }
1551  Kokkos::deep_copy(v1, h_v1);
1552  Kokkos::deep_copy(v2, h_v2);
1553 
1554  // Launch kernel
1555  ViewType v3("view3", num_rows, fad_size+1);
1557 
1558  // Copy back
1559  host_view_type h_v3 = Kokkos::create_mirror_view(v3);
1560  Kokkos::deep_copy(h_v3, v3);
1561 
1562  // Check
1563  success = true;
1564  TEUCHOS_TEST_EQUALITY(v3.rank(), 1, out, success);
1565  for (size_type i=0; i<num_rows; ++i) {
1566  FadType f1 =
1567  generate_fad<FadType>(num_rows, size_type(2), fad_size, i, size_type(0));
1568  FadType f2 =
1569  generate_fad<FadType>(num_rows, size_type(2), fad_size, i, size_type(1));
1570  FadType f3 = f1*f2;
1571  success = success && checkFads(f3, h_v3(i), out);
1572  }
1573 }
1574 
1576  Kokkos_View_Fad, SubdynrankviewCol, FadType, Layout, Device )
1577 {
1578  typedef Kokkos::DynRankView<FadType,Layout,Device> ViewType;
1579  typedef typename ViewType::size_type size_type;
1580  typedef typename ViewType::host_mirror_type host_view_type;
1581 
1582  const size_type num_rows = global_num_rows;
1583  const size_type num_cols = global_num_cols;
1584  const size_type fad_size = global_fad_size;
1585 
1586  // Create and fill view
1587  ViewType v("view", num_rows, num_cols, fad_size+1);
1588  host_view_type h_v = Kokkos::create_mirror_view(v);
1589  for (size_type i=0; i<num_rows; ++i) {
1590  for (size_type j=0; j<num_cols; ++j) {
1591  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
1592  h_v(i,j) = f;
1593  }
1594  }
1595  Kokkos::deep_copy(v, h_v);
1596 
1597  // Create subview of first column
1598  size_type col = 1;
1599  auto s = Kokkos::subdynrankview(v, Kokkos::ALL(), col);
1600 
1601  // Copy back
1602  typedef decltype(s) SubviewType;
1603  typedef typename SubviewType::host_mirror_type HostSubviewType;
1604 
1605  // Note: don't create h_s through create_mirror_view and deep_copy
1606  // since Kokkos doesn't support deep_copy of non-contiguous views
1607  //HostSubviewType h_s = Kokkos::create_mirror_view(s);
1608  //Kokkos::deep_copy(h_s, s);
1609  HostSubviewType h_s = Kokkos::subdynrankview(h_v, Kokkos::ALL(), col);
1610 
1611  // Check
1612  success = true;
1613  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(s), fad_size+1, out, success);
1614  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(h_s), fad_size+1, out, success);
1615  TEUCHOS_TEST_EQUALITY(h_s.extent(0), num_rows, out, success);
1616  TEUCHOS_TEST_EQUALITY(h_s.extent(1), 1, out, success);
1617  TEUCHOS_TEST_EQUALITY(h_s.extent(7), 1, out, success);
1618 
1619  for (size_type i=0; i<num_rows; ++i) {
1620  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, col);
1621  success = success && checkFads(f, h_s(i), out);
1622  }
1623 }
1624 
1626  Kokkos_View_Fad, SubdynrankviewRow, FadType, Layout, Device )
1627 {
1628  typedef Kokkos::DynRankView<FadType,Layout,Device> ViewType;
1629  typedef typename ViewType::size_type size_type;
1630  typedef typename ViewType::host_mirror_type host_view_type;
1631 
1632  const size_type num_rows = global_num_rows;
1633  const size_type num_cols = global_num_cols;
1634  const size_type num_planes = 9;
1635  const size_type fad_size = global_fad_size;
1636 
1637  // Create and fill view
1638  ViewType v("view", num_rows, num_cols, num_planes, fad_size+1);
1639  host_view_type h_v = Kokkos::create_mirror_view(v);
1640  for (size_type i=0; i<num_rows; ++i) {
1641  for (size_type j=0; j<num_cols; ++j) {
1642  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
1643  for (size_type k=0; k<num_planes; ++k) {
1644  h_v(i,j,k) = (k+1)*f;
1645  }
1646  }
1647  }
1648  Kokkos::deep_copy(v, h_v);
1649 
1650  // Create subview of first column
1651  size_type row = 2;
1652  auto s = Kokkos::subdynrankview(v, row, Kokkos::ALL(), Kokkos::ALL());
1653 
1654  // Copy back
1655  typedef decltype(s) SubviewType;
1656  typedef typename SubviewType::host_mirror_type HostSubviewType;
1657 
1658  // Note: don't create h_s through create_mirror_view and deep_copy
1659  // since Kokkos doesn't support deep_copy of non-contiguous views
1660  //HostSubviewType h_s = Kokkos::create_mirror_view(s);
1661  //Kokkos::deep_copy(h_s, s);
1662  HostSubviewType h_s =
1663  Kokkos::subdynrankview(h_v, row, Kokkos::ALL(), Kokkos::ALL());
1664 
1665  // Check
1666  success = true;
1667  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(s), fad_size+1, out, success);
1668  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(h_s), fad_size+1, out, success);
1669  TEUCHOS_TEST_EQUALITY(h_s.extent(0), num_cols, out, success);
1670  TEUCHOS_TEST_EQUALITY(h_s.extent(1), num_planes, out, success);
1671  TEUCHOS_TEST_EQUALITY(h_s.extent(2), 1, out, success);
1672  TEUCHOS_TEST_EQUALITY(h_s.extent(7), 1, out, success);
1673 
1674  for (size_type j=0; j<num_cols; ++j) {
1675  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, row, j);
1676  for (size_type k=0; k<num_planes; ++k) {
1677  FadType g = (k+1)*f;
1678  success = success && checkFads(g, h_s(j,k), out);
1679  }
1680  }
1681 }
1682 
1684  Kokkos_View_Fad, SubdynrankviewScalar, FadType, Layout, Device )
1685 {
1686  typedef Kokkos::DynRankView<FadType,Layout,Device> ViewType;
1687  typedef typename ViewType::size_type size_type;
1688  typedef typename ViewType::host_mirror_type host_view_type;
1689 
1690  const size_type num_rows = global_num_rows;
1691  const size_type num_cols = global_num_cols;
1692  const size_type fad_size = global_fad_size;
1693 
1694  // Create and fill view
1695  ViewType v("view", num_rows, num_cols, fad_size+1);
1696  host_view_type h_v = Kokkos::create_mirror_view(v);
1697  for (size_type i=0; i<num_rows; ++i) {
1698  for (size_type j=0; j<num_cols; ++j) {
1699  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
1700  h_v(i,j) = f;
1701  }
1702  }
1703  Kokkos::deep_copy(v, h_v);
1704 
1705  // Create subview of first column
1706  size_type row = 3;
1707  size_type col = 1;
1708  auto s = Kokkos::subdynrankview(v, row, col);
1709 
1710  // Copy back
1711  typedef decltype(s) SubviewType;
1712  typedef typename SubviewType::host_mirror_type HostSubviewType;
1713 
1714  // Note: don't create h_s through create_mirror_view and deep_copy
1715  // since Kokkos doesn't support deep_copy of non-contiguous views
1716  //HostSubviewType h_s = Kokkos::create_mirror_view(s);
1717  //Kokkos::deep_copy(h_s, s);
1718  HostSubviewType h_s = Kokkos::subdynrankview(h_v, row, col);
1719 
1720  // Check
1721  success = true;
1722  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(s), fad_size+1, out, success);
1723  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(h_s), fad_size+1, out, success);
1724  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, row, col);
1725  success = success && checkFads(f, h_s(), out);
1726 }
1727 
1728 #else
1729 
1731  Kokkos_View_Fad, DynRankDimensionScalar, FadType, Layout, Device ) {}
1733  Kokkos_View_Fad, DynRankAssignStatic0, FadType, Layout, Device ) {}
1735  Kokkos_View_Fad, DynRankAssignStatic1, FadType, Layout, Device ) {}
1737  Kokkos_View_Fad, DynRankAssignStatic2, FadType, Layout, Device ) {}
1739  Kokkos_View_Fad, DynRankMultiply, FadType, Layout, Device ) {}
1741  Kokkos_View_Fad, SubdynrankviewCol, FadType, Layout, Device ) {}
1743  Kokkos_View_Fad, SubdynrankviewRow, FadType, Layout, Device ) {}
1745  Kokkos_View_Fad, SubdynrankviewScalar, FadType, Layout, Device ) {}
1746 
1747 #endif
1748 
1750  Kokkos_View_Fad, Subview, FadType, Layout, Device )
1751 {
1752  typedef Kokkos::View<FadType**,Layout,Device> ViewType;
1753  typedef typename ViewType::size_type size_type;
1754  typedef typename ViewType::host_mirror_type host_view_type;
1755 
1756  const size_type num_rows = global_num_rows;
1757  const size_type num_cols = global_num_cols;
1758  const size_type fad_size = global_fad_size;
1759 
1760  // Create and fill view
1761  ViewType v;
1762 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1763  v = ViewType ("view", num_rows, num_cols);
1764 #else
1765  v = ViewType ("view", num_rows, num_cols, fad_size+1);
1766 #endif
1767  host_view_type h_v = Kokkos::create_mirror_view(v);
1768  for (size_type i=0; i<num_rows; ++i) {
1769  for (size_type j=0; j<num_cols; ++j) {
1770  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
1771  h_v(i,j) = f;
1772  }
1773  }
1774  Kokkos::deep_copy(v, h_v);
1775 
1776  // Create subview of first column
1777  size_type col = 1;
1778  auto s = Kokkos::subview(v, Kokkos::ALL(), col);
1779 
1780  // Copy back
1781  typedef decltype(s) SubviewType;
1782  typedef typename SubviewType::host_mirror_type HostSubviewType;
1783 
1784  // Note: don't create h_s through create_mirror_view and deep_copy
1785  // since Kokkos doesn't support deep_copy of non-contiguous views
1786  //HostSubviewType h_s = Kokkos::create_mirror_view(s);
1787  //Kokkos::deep_copy(h_s, s);
1788  HostSubviewType h_s = Kokkos::subview(h_v, Kokkos::ALL(), col);
1789 
1790  // Check
1791  success = true;
1792 #if defined(HAVE_SACADO_VIEW_SPEC) && !defined(SACADO_DISABLE_FAD_VIEW_SPEC)
1793  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(s), fad_size+1, out, success);
1794  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(h_s), fad_size+1, out, success);
1795 #endif
1796  for (size_type i=0; i<num_rows; ++i) {
1797  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, col);
1798  success = success && checkFads(f, h_s(i), out);
1799  }
1800 }
1801 
1803  Kokkos_View_Fad, Subview2, FadType, Layout, Device )
1804 {
1805  typedef Kokkos::View<FadType***,Layout,Device> ViewType;
1806  typedef typename ViewType::host_mirror_type host_view_type;
1807 
1808  // Test various subview operations to check the resulting indexing is correct.
1809  // We only need to run these tests on the host because the indexing does
1810  // not depend on the device.
1811 
1812  const int num_cell = 5;
1813  const int num_qp = 4;
1814  const int num_dim = 3;
1815  const int num_deriv = 2;
1816 
1817  // Create and fill view
1818  host_view_type v;
1819 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1820  v = host_view_type ("view", num_cell, num_qp, num_dim);
1821 #else
1822  v = host_view_type ("view", num_cell, num_qp, num_dim, num_deriv+1);
1823 #endif
1824  for (int cell=0; cell < num_cell; ++cell) {
1825  for (int qp=0; qp < num_qp; ++qp) {
1826  for (int dim = 0; dim < num_dim; ++dim) {
1827  v(cell,qp,dim).val() = 100.*cell + 10.*qp + 1.*dim;
1828  for (int deriv = 0; deriv < num_deriv; ++deriv) {
1829  v(cell,qp,dim).fastAccessDx(deriv) = v(cell,qp,dim).val() + (1.0*deriv)/10.;
1830  }
1831  }
1832  }
1833  }
1834 
1835  success = true;
1836 
1837  out << "checking subview(v,ALL,*,*)..." << std::endl;
1838  for (int qp=0; qp < num_qp; ++qp) {
1839  for (int dim=0; dim < num_dim; ++dim) {
1840  auto v_tmp = subview(v,Kokkos::ALL(),qp,dim);
1841  for (int cell=0; cell < num_cell; ++cell) {
1842  out << "\tChecking (" << cell << "," << qp << "," << dim << ")" << std::endl;
1843  success = success && checkFads(v(cell,qp,dim), v_tmp(cell), out);
1844  }
1845  }
1846  }
1847 
1848  out << "checking subview(v,*,ALL,*)..." << std::endl;
1849  for (int cell=0; cell < num_cell; ++cell) {
1850  for (int dim=0; dim < num_dim; ++dim) {
1851  auto v_tmp = subview(v,cell,Kokkos::ALL(),dim);
1852  for (int qp=0; qp < num_qp; ++qp) {
1853  out << "\tChecking (" << cell << "," << qp << "," << dim << ")" << std::endl;
1854  success = success && checkFads(v(cell,qp,dim), v_tmp(qp), out);
1855  }
1856  }
1857  }
1858 
1859  out << "checking subview(v,*,*,ALL)..." << std::endl;
1860  for (int cell=0; cell < num_cell; ++cell) {
1861  for (int qp=0; qp < num_qp; ++qp) {
1862  auto v_tmp = subview(v,cell,qp,Kokkos::ALL());
1863  for (int dim=0; dim < num_dim; ++dim) {
1864  out << "\tChecking (" << cell << "," << qp << "," << dim << ")" << std::endl;
1865  success = success && checkFads(v(cell,qp,dim), v_tmp(dim), out);
1866  }
1867  }
1868  }
1869 
1870  out << "checking subview(v,ALL,ALL,*)..." << std::endl;
1871  for (int dim=0; dim < num_dim; ++dim) {
1872  auto v_tmp = subview(v,Kokkos::ALL(),Kokkos::ALL(),dim);
1873  for (int cell=0; cell < num_cell; ++cell) {
1874  for (int qp=0; qp < num_qp; ++qp) {
1875  out << "\tChecking (" << cell << "," << qp << "," << dim << ")" << std::endl;
1876  success = success && checkFads(v(cell,qp,dim), v_tmp(cell,qp), out);
1877  }
1878  }
1879  }
1880 
1881  out << "checking subview(v,*,ALL,ALL)..." << std::endl;
1882  for (int cell=0; cell < num_cell; ++cell) {
1883  auto v_tmp = subview(v,cell,Kokkos::ALL(),Kokkos::ALL());
1884  for (int qp=0; qp < num_qp; ++qp) {
1885  for (int dim=0; dim < num_dim; ++dim) {
1886  out << "\tChecking (" << cell << "," << qp << "," << dim << ")" << std::endl;
1887  success = success && checkFads(v(cell,qp,dim), v_tmp(qp,dim), out);
1888  }
1889  }
1890  }
1891 
1892  out << "checking subview(v,ALL,*,ALL)..." << std::endl;
1893  for (int qp=0; qp < num_qp; ++qp) {
1894  auto v_tmp = subview(v,Kokkos::ALL(),qp,Kokkos::ALL());
1895  for (int cell=0; cell < num_cell; ++cell) {
1896  for (int dim=0; dim < num_dim; ++dim) {
1897  out << "\tChecking (" << cell << "," << qp << "," << dim << ")" << std::endl;
1898  success = success && checkFads(v(cell,qp,dim), v_tmp(cell,dim), out);
1899  }
1900  }
1901  }
1902 
1903  out << "checking subview(v,range,range,range)..." << std::endl;
1904  auto v_tmp = subview(v,std::make_pair(1,5),std::make_pair(1,4),std::make_pair(1,3));
1905  for (int cell=1; cell < num_cell; ++cell) {
1906  for (int qp=1; qp < num_qp; ++qp) {
1907  for (int dim=1; dim < num_dim; ++dim) {
1908  out << "\tChecking (" << cell << "," << qp << "," << dim << ")" << std::endl;
1909  success = success && checkFads(v(cell,qp,dim), v_tmp(cell-1,qp-1,dim-1), out);
1910  }
1911  }
1912  }
1913 }
1914 
1915 #ifdef SACADO_NEW_FAD_DESIGN_IS_DEFAULT
1917  Kokkos_View_Fad, ConstViewAssign, FadType, Layout, Device )
1918 {
1919  typedef Kokkos::View<FadType*,Layout,Device> ViewType;
1920  typedef Kokkos::View<const FadType,Layout,Device> ConstViewType;
1921  typedef typename ViewType::size_type size_type;
1922  typedef typename ViewType::host_mirror_type host_view_type;
1923  typedef typename ViewType::execution_space exec_space;
1924 
1925  const size_type num_rows = global_num_rows;
1926  const size_type fad_size = global_fad_size;
1927 
1928  // Create and fill view
1929 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1930  ViewType v1("view1", num_rows);
1931 #else
1932  ViewType v1("view1", num_rows, fad_size+1);
1933 #endif
1934  host_view_type h_v1 = Kokkos::create_mirror_view(v1);
1935  for (size_type i=0; i<num_rows; ++i) {
1936  FadType f = generate_fad<FadType>(num_rows, size_type(1), fad_size, i,
1937  size_type(0));
1938  h_v1(i) = f;
1939  }
1940  Kokkos::deep_copy(v1, h_v1);
1941 
1942 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
1943  ViewType v2("view2", num_rows);
1944 #else
1945  ViewType v2("view2", num_rows, fad_size+1);
1946 #endif
1947 
1948  static const size_type stride = Kokkos::ViewScalarStride<ViewType>::stride;
1949 #if defined (KOKKOS_ENABLE_CUDA) && defined (SACADO_VIEW_CUDA_HIERARCHICAL)
1950  const bool use_team =
1954  ( stride > 1 );
1955 #elif defined (KOKKOS_ENABLE_CUDA) && defined (SACADO_VIEW_CUDA_HIERARCHICAL_DFAD)
1956  const bool use_team =
1961 #elif defined (KOKKOS_ENABLE_HIP) && defined (SACADO_VIEW_CUDA_HIERARCHICAL)
1962  const bool use_team =
1966  ( stride > 1 );
1967 #elif defined (KOKKOS_ENABLE_HIP) && defined (SACADO_VIEW_CUDA_HIERARCHICAL_DFAD)
1968  const bool use_team =
1973 #else
1974  const bool use_team = false;
1975 #endif
1976 
1977  if (use_team) {
1978  typedef Kokkos::TeamPolicy<exec_space> team_policy;
1979  Kokkos::parallel_for(team_policy(num_rows, 1, stride),
1980  KOKKOS_LAMBDA(typename team_policy::member_type team)
1981  {
1982  const int i = team.league_rank();
1983  typename ConstViewType::reference_type x = v1(i);
1984  v2(i) = x;
1985  });
1986  }
1987  else {
1988  Kokkos::parallel_for(Kokkos::RangePolicy<exec_space>(0,num_rows),
1989  KOKKOS_LAMBDA(const int i)
1990  {
1991  typename ConstViewType::reference_type x = v1(i);
1992  v2(i) = x;
1993  });
1994  }
1995 
1996  // Copy back
1997  host_view_type h_v2 = Kokkos::create_mirror_view(v2);
1998  Kokkos::deep_copy(h_v2, v2);
1999 
2000  // Check
2001  success = true;
2002  for (size_type i=0; i<num_rows; ++i) {
2003  FadType f = generate_fad<FadType>(num_rows, size_type(1), fad_size, i,
2004  size_type(0));
2005  success = success && checkFads(f, h_v2(i), out);
2006  }
2007 }
2008 #else
2010  Kokkos_View_Fad, ConstViewAssign, FadType, Layout, Device ) {}
2011 #endif
2012 
2013 // Tests that require view spec
2014 
2015 #if defined(HAVE_SACADO_VIEW_SPEC) && !defined(SACADO_DISABLE_FAD_VIEW_SPEC)
2017  Kokkos_View_Fad, ShmemSize, FadType, Layout, Device )
2018 {
2019  typedef Kokkos::View<FadType**,Layout,Device> ViewType;
2020  typedef typename FadType::value_type value_type;
2021  typedef typename ViewType::size_type size_type;
2022 
2023  const size_type num_rows = global_num_rows;
2024  const size_type num_cols = global_num_cols;
2025  const size_type fad_size = global_fad_size;
2026 
2027  // Compute shared memory size for View
2028  const size_type shmem_size =
2029  ViewType::shmem_size(num_rows, num_cols, fad_size+1);
2030 
2031  // Check
2032  const size_type align = 8;
2033  const size_type mask = align - 1;
2034  ViewType v;
2035 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
2036  v = ViewType ("view", num_rows, num_cols);
2037 #else
2038  v = ViewType ("view", num_rows, num_cols, fad_size+1);
2039 #endif
2040 #if defined(SACADO_HAS_NEW_KOKKOS_VIEW_IMPL) || (defined(SACADO_DISABLE_FAD_VIEW_SPEC) && !defined(KOKKOS_ENABLE_IMPL_VIEW_LEGACY))
2041  size_t scratch_value_alignment =
2042  Kokkos::max({sizeof(value_type),
2043  alignof(value_type),
2044  static_cast<size_t>(
2045  ViewType::execution_space::scratch_memory_space::ALIGN)});
2046  const size_type shmem_size_expected =
2047  sizeof(value_type) * global_num_rows * global_num_cols * (fad_size+1) + scratch_value_alignment;
2048 #else
2049  const size_type shmem_size_expected =
2050  (( sizeof(value_type) * global_num_rows * global_num_cols * (fad_size+1) + mask ) & ~mask) + sizeof(typename ViewType::traits::value_type);
2051 #endif
2052  TEUCHOS_TEST_EQUALITY(shmem_size, shmem_size_expected, out, success);
2053 }
2054 
2056  Kokkos_View_Fad, Unmanaged, FadType, Layout, Device )
2057 {
2058  // For LayoutContiguous or LayoutNatural, strip out the layout they are templated on
2059  typedef typename Kokkos::inner_layout<Layout>::type TestLayout;
2060 
2061  typedef typename FadType::value_type scalar_type;
2062  typedef Kokkos::View<scalar_type***,TestLayout,Device> ViewType;
2063  typedef Kokkos::View<FadType**,TestLayout,Device,Kokkos::MemoryUnmanaged> FadViewType;
2064  typedef typename ViewType::size_type size_type;
2065  typedef typename ViewType::host_mirror_type host_view_type;
2066  typedef typename FadViewType::host_mirror_type fad_host_view_type;
2067 
2068  const size_type num_rows = global_num_rows;
2069  const size_type num_cols = global_num_cols;
2070  const size_type fad_size = global_fad_size;
2071 
2072  // Create and fill view
2073  ViewType v;
2074  host_view_type h_v;
2077  v = ViewType ("view", fad_size+1, num_rows, num_cols);
2078  h_v = Kokkos::create_mirror_view(v);
2079  for (size_type i=0; i<num_rows; ++i) {
2080  for (size_type j=0; j<num_cols; ++j) {
2081  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2082  for (size_type k=0; k<fad_size; k++)
2083  h_v(k,i,j) = f.dx(k);
2084  h_v(fad_size,i,j) = f.val();
2085  }
2086  }
2087  }
2088  else {
2089  v = ViewType ("view", num_rows, num_cols, fad_size+1);
2090  h_v = Kokkos::create_mirror_view(v);
2091  for (size_type i=0; i<num_rows; ++i) {
2092  for (size_type j=0; j<num_cols; ++j) {
2093  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2094  for (size_type k=0; k<fad_size; k++)
2095  h_v(i,j,k) = f.dx(k);
2096  h_v(i,j,fad_size) = f.val();
2097  }
2098  }
2099  }
2100  Kokkos::deep_copy(v, h_v);
2101 
2102  // Create unmanaged view
2103  FadViewType v_fad;//( v.data(), num_rows, num_cols, fad_size+1);
2104  fad_host_view_type h_v_fad;
2105 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
2106  v_fad = FadViewType ( v.data(), num_rows, num_cols);
2107  h_v_fad = fad_host_view_type ("host_view_fad", num_rows, num_cols);
2108 #else
2109  v_fad = FadViewType ( v.data(), num_rows, num_cols, fad_size+1);
2110  h_v_fad = fad_host_view_type ("host_view_fad", num_rows, num_cols, fad_size+1);
2111 #endif
2112 
2113  // Copy back -- can't use create_mirror_view() because v_fad is unmanaged
2114  Kokkos::deep_copy(h_v_fad, v_fad);
2115 
2116  // Check
2117  success = true;
2118  for (size_type i=0; i<num_rows; ++i) {
2119  for (size_type j=0; j<num_cols; ++j) {
2120  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2121  success = success && checkFads(f, h_v_fad(i,j), out);
2122  }
2123  }
2124 }
2125 
2127  Kokkos_View_Fad, Unmanaged2, FadType, Layout, Device )
2128 {
2129  // For LayoutContiguous or LayoutNatural, strip out the layout they are templated on
2130  typedef typename Kokkos::inner_layout<Layout>::type TestLayout;
2131 
2132  typedef typename FadType::value_type scalar_type;
2133  typedef Kokkos::View<scalar_type***,TestLayout,Device> ViewType;
2134  typedef Kokkos::View<FadType**,TestLayout,Device> FadViewType;
2135  typedef typename ViewType::size_type size_type;
2136  typedef typename ViewType::host_mirror_type host_view_type;
2137  typedef typename FadViewType::host_mirror_type fad_host_view_type;
2138 
2139  const size_type num_rows = global_num_rows;
2140  const size_type num_cols = global_num_cols;
2141  const size_type fad_size = global_fad_size;
2142 
2143  // Create and fill view
2144  ViewType v;
2145  host_view_type h_v;
2148  v = ViewType ("view", fad_size+1, num_rows, num_cols);
2149  h_v = Kokkos::create_mirror_view(v);
2150  for (size_type i=0; i<num_rows; ++i) {
2151  for (size_type j=0; j<num_cols; ++j) {
2152  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2153  for (size_type k=0; k<fad_size; k++)
2154  h_v(k,i,j) = f.dx(k);
2155  h_v(fad_size,i,j) = f.val();
2156  }
2157  }
2158  }
2159  else {
2160  v = ViewType ("view", num_rows, num_cols, fad_size+1);
2161  h_v = Kokkos::create_mirror_view(v);
2162  for (size_type i=0; i<num_rows; ++i) {
2163  for (size_type j=0; j<num_cols; ++j) {
2164  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2165  for (size_type k=0; k<fad_size; k++)
2166  h_v(i,j,k) = f.dx(k);
2167  h_v(i,j,fad_size) = f.val();
2168  }
2169  }
2170  }
2171  Kokkos::deep_copy(v, h_v);
2172 
2173  // Create unmanaged view
2174  FadViewType v_fad;//( v.data(), num_rows, num_cols, fad_size+1);
2175  fad_host_view_type h_v_fad;
2176 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
2177  v_fad = FadViewType ( v.data(), num_rows, num_cols);
2178  h_v_fad = fad_host_view_type ("host_view_fad", num_rows, num_cols);
2179 #else
2180  v_fad = FadViewType ( v.data(), num_rows, num_cols, fad_size+1);
2181  h_v_fad = fad_host_view_type ("host_view_fad", num_rows, num_cols, fad_size+1);
2182 #endif
2183 
2184  // Copy back -- can't use create_mirror_view() because v_fad is unmanaged
2185  Kokkos::deep_copy(h_v_fad, v_fad);
2186 
2187  // Check
2188  success = true;
2189  for (size_type i=0; i<num_rows; ++i) {
2190  for (size_type j=0; j<num_cols; ++j) {
2191  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2192  success = success && checkFads(f, h_v_fad(i,j), out);
2193  }
2194  }
2195 }
2196 
2198  Kokkos_View_Fad, UnmanagedConst, FadType, Layout, Device )
2199 {
2200  // For LayoutContiguous or LayoutNatural, strip out the layout they are templated on
2201  typedef typename Kokkos::inner_layout<Layout>::type TestLayout;
2202 
2203  typedef typename FadType::value_type scalar_type;
2204  typedef Kokkos::View<scalar_type***,TestLayout,Device> ViewType;
2205  typedef Kokkos::View<const scalar_type***,TestLayout,Device> ConstViewType;
2206  typedef Kokkos::View<FadType**,TestLayout,Device,Kokkos::MemoryUnmanaged> FadViewType;
2207  typedef Kokkos::View<const FadType**,TestLayout,Device,Kokkos::MemoryUnmanaged> ConstFadViewType;
2208  typedef typename ViewType::size_type size_type;
2209  typedef typename ViewType::host_mirror_type host_view_type;
2210  typedef typename FadViewType::host_mirror_type fad_host_view_type;
2211 
2212  const size_type num_rows = global_num_rows;
2213  const size_type num_cols = global_num_cols;
2214  const size_type fad_size = global_fad_size;
2215 
2216  // Create and fill view
2217  ViewType v;
2218  host_view_type h_v;
2221  v = ViewType ("view", fad_size+1, num_rows, num_cols);
2222  h_v = Kokkos::create_mirror_view(v);
2223  for (size_type i=0; i<num_rows; ++i) {
2224  for (size_type j=0; j<num_cols; ++j) {
2225  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2226  for (size_type k=0; k<fad_size; k++)
2227  h_v(k,i,j) = f.dx(k);
2228  h_v(fad_size,i,j) = f.val();
2229  }
2230  }
2231  }
2232  else {
2233  v = ViewType ("view", num_rows, num_cols, fad_size+1);
2234  h_v = Kokkos::create_mirror_view(v);
2235  for (size_type i=0; i<num_rows; ++i) {
2236  for (size_type j=0; j<num_cols; ++j) {
2237  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2238  for (size_type k=0; k<fad_size; k++)
2239  h_v(i,j,k) = f.dx(k);
2240  h_v(i,j,fad_size) = f.val();
2241  }
2242  }
2243  }
2244  Kokkos::deep_copy(v, h_v);
2245  ConstViewType v_const = v;
2246 
2247  // Create unmanaged view
2248 
2249  ConstFadViewType v_fad;//( v.data(), num_rows, num_cols, fad_size+1);
2250  fad_host_view_type h_v_fad;
2251 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
2252  v_fad = ConstFadViewType ( v_const.data(), num_rows, num_cols);
2253  h_v_fad = fad_host_view_type ("host_view_fad", num_rows, num_cols);
2254 #else
2255  v_fad = ConstFadViewType ( v_const.data(), num_rows, num_cols, fad_size+1);
2256  h_v_fad = fad_host_view_type ("host_view_fad", num_rows, num_cols, fad_size+1);
2257 #endif
2258 
2259  // Copy back -- can't use create_mirror_view() because v_fad is unmanaged
2260  Kokkos::deep_copy(h_v_fad, v_fad);
2261 
2262  // Check
2263  success = true;
2264  for (size_type i=0; i<num_rows; ++i) {
2265  for (size_type j=0; j<num_cols; ++j) {
2266  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2267  success = success && checkFads(f, h_v_fad(i,j), out);
2268  }
2269  }
2270 }
2271 
2273  Kokkos_View_Fad, UnmanagedConst2, FadType, Layout, Device )
2274 {
2275  // For LayoutContiguous or LayoutNatural, strip out the layout they are templated on
2276  typedef typename Kokkos::inner_layout<Layout>::type TestLayout;
2277  typedef typename FadType::value_type scalar_type;
2278  typedef Kokkos::View<scalar_type***,TestLayout,Device> ViewType;
2279  typedef Kokkos::View<const scalar_type***,TestLayout,Device> ConstViewType;
2280  typedef Kokkos::View<FadType**,TestLayout,Device> FadViewType;
2281  typedef Kokkos::View<const FadType**,TestLayout,Device> ConstFadViewType;
2282  typedef typename ViewType::size_type size_type;
2283  typedef typename ViewType::host_mirror_type host_view_type;
2284  typedef typename FadViewType::host_mirror_type fad_host_view_type;
2285 
2286  const size_type num_rows = global_num_rows;
2287  const size_type num_cols = global_num_cols;
2288  const size_type fad_size = global_fad_size;
2289 
2290  // Create and fill view
2291  ViewType v;
2292  host_view_type h_v;
2295  v = ViewType ("view", fad_size+1, num_rows, num_cols);
2296  h_v = Kokkos::create_mirror_view(v);
2297  for (size_type i=0; i<num_rows; ++i) {
2298  for (size_type j=0; j<num_cols; ++j) {
2299  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2300  for (size_type k=0; k<fad_size; k++)
2301  h_v(k,i,j) = f.dx(k);
2302  h_v(fad_size,i,j) = f.val();
2303  }
2304  }
2305  }
2306  else {
2307  v = ViewType ("view", num_rows, num_cols, fad_size+1);
2308  h_v = Kokkos::create_mirror_view(v);
2309  for (size_type i=0; i<num_rows; ++i) {
2310  for (size_type j=0; j<num_cols; ++j) {
2311  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2312  for (size_type k=0; k<fad_size; k++)
2313  h_v(i,j,k) = f.dx(k);
2314  h_v(i,j,fad_size) = f.val();
2315  }
2316  }
2317  }
2318  Kokkos::deep_copy(v, h_v);
2319  ConstViewType v_const = v;
2320 
2321  // Create unmanaged view
2322  ConstFadViewType v_fad;
2323  fad_host_view_type h_v_fad;
2324 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
2325  v_fad = ConstFadViewType (v_const.data(), num_rows, num_cols);
2326  h_v_fad = fad_host_view_type ("host_view_fad", num_rows, num_cols);
2327 #else
2328  v_fad = ConstFadViewType (v_const.data(), num_rows, num_cols, fad_size+1);
2329  h_v_fad = fad_host_view_type ("host_view_fad", num_rows, num_cols, fad_size+1);
2330 #endif
2331 
2332  // Copy back -- can't use create_mirror_view() because v_fad is unmanaged
2333  Kokkos::deep_copy(h_v_fad, v_fad);
2334 
2335  // Check
2336  success = true;
2337  for (size_type i=0; i<num_rows; ++i) {
2338  for (size_type j=0; j<num_cols; ++j) {
2339  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2340  success = success && checkFads(f, h_v_fad(i,j), out);
2341  }
2342  }
2343 }
2344 
2345 // This test checks we can allocate a view
2346 // with SFad without specifying the fad size in the constructor
2348  Kokkos_View_Fad, SFadNoSizeArg, FadType, Layout, Device )
2349 {
2350  typedef Kokkos::View<FadType**,Layout,Device> ViewType;
2351  typedef typename ViewType::size_type size_type;
2352  typedef typename ViewType::host_mirror_type host_view_type;
2353 
2354  const size_type num_rows = global_num_rows;
2355  const size_type num_cols = global_num_cols;
2356  const size_type fad_size = global_fad_size;
2357 
2358  // Create and fill view
2359  ViewType v("view", num_rows, num_cols);
2360  host_view_type h_v = Kokkos::create_mirror_view(v);
2361  for (size_type i=0; i<num_rows; ++i) {
2362  for (size_type j=0; j<num_cols; ++j) {
2363  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2364  h_v(i,j) = f;
2365  }
2366  }
2367  Kokkos::deep_copy(v, h_v);
2368 
2369  // Copy back
2370  Kokkos::deep_copy(h_v, v);
2371 
2372  // Check
2373  success = true;
2374  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(v), fad_size+1, out, success);
2375  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(h_v), fad_size+1, out, success);
2376  for (size_type i=0; i<num_rows; ++i) {
2377  for (size_type j=0; j<num_cols; ++j) {
2378  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2379  success = success && checkFads(f, h_v(i,j), out);
2380  }
2381  }
2382 }
2383 
2384 #ifndef SACADO_HAS_NEW_KOKKOS_VIEW_IMPL
2386  Kokkos_View_Fad, Partition, FadType, Layout, Device )
2387 {
2388 #if !defined(SACADO_VIEW_CUDA_HIERARCHICAL) && !defined(SACADO_VIEW_CUDA_HIERARCHICAL_DFAD)
2389  typedef Kokkos::View<FadType**,Layout,Device> ViewType;
2390  typedef typename ViewType::size_type size_type;
2391  typedef typename ViewType::host_mirror_type host_view_type;
2392 
2393  const size_type num_rows = global_num_rows;
2394  const size_type num_cols = global_num_cols;
2395  const size_type fad_size = global_fad_size;
2396 
2397  // Create and fill view
2398  ViewType v;
2399 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
2400  v = ViewType ("view", num_rows, num_cols);
2401 #else
2402  v = ViewType ("view", num_rows, num_cols, fad_size+1);
2403 #endif
2404  host_view_type h_v = Kokkos::create_mirror_view(v);
2405 
2406  for (size_type i=0; i<num_rows; ++i) {
2407  for (size_type j=0; j<num_cols; ++j) {
2408  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2409  h_v(i,j) = f;
2410  }
2411  }
2412  Kokkos::deep_copy(v, h_v);
2413 
2414  // Copy back
2415  Kokkos::deep_copy(h_v, v);
2416 
2417  // Partition derivative array of h_v into 2, first one starting at index 0,
2418  // the second at 1
2419  const size_type stride = 2;
2420  auto h_v1 = Kokkos::partition<2>(h_v, 0, stride);
2421  auto h_v2 = Kokkos::partition<2>(h_v, 1, stride);
2422 
2423  // Check
2424  const size_type fad_size_1 = (fad_size + stride - 0 - 1) / stride;
2425  const size_type fad_size_2 = (fad_size + stride - 1 - 1) / stride;
2426  success = true;
2427  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(h_v1), fad_size_1+1, out, success);
2428  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(h_v2), fad_size_2+1, out, success);
2429  for (size_type i=0; i<num_rows; ++i) {
2430  for (size_type j=0; j<num_cols; ++j) {
2431  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2432  Sacado::Fad::DFad<double> f1( fad_size_1, f.val() );
2433  Sacado::Fad::DFad<double> f2( fad_size_2, f.val() );
2434  for (unsigned int k=0; k<fad_size_1; ++k)
2435  if (2*k < fad_size) f1.fastAccessDx(k) = f.dx(2*k);
2436  for (unsigned int k=0; k<fad_size_2; ++k)
2437  if (2*k+1 < fad_size) f2.fastAccessDx(k) = f.dx(2*k+1);
2438  success = success && checkFads(f1, h_v1(i,j), out);
2439  success = success && checkFads(f2, h_v2(i,j), out);
2440  }
2441  }
2442 #endif
2443 }
2444 #endif // SACADO_HAS_NEW_KOKKOS_VIEW_IMPL
2445 
2447  Kokkos_View_Fad, AssignLayoutContiguousToLayoutStride, FadType, Layout, Device )
2448 {
2449  typedef Kokkos::View<FadType**,Kokkos::LayoutContiguous<Layout>,Device> ContViewType;
2450  typedef Kokkos::View<FadType**,Kokkos::LayoutStride,Device> StrideViewType;
2451  typedef typename ContViewType::size_type size_type;
2452  typedef typename ContViewType::host_mirror_type cont_host_view_type;
2453  typedef typename StrideViewType::host_mirror_type stride_host_view_type;
2454 
2455  const size_type num_rows = global_num_rows;
2456  const size_type num_cols = global_num_cols;
2457  const size_type fad_size = global_fad_size;
2458 
2459  // Create and fill view
2460  ContViewType v;
2461 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
2462  v = ContViewType ("view", num_rows, num_cols);
2463 #else
2464  v = ContViewType ("view", num_rows, num_cols, fad_size+1);
2465 #endif
2466  cont_host_view_type h_v = Kokkos::create_mirror_view(v);
2467 
2468  for (size_type i=0; i<num_rows; ++i) {
2469  for (size_type j=0; j<num_cols; ++j) {
2470  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2471  h_v(i,j) = f;
2472  }
2473  }
2474  Kokkos::deep_copy(v, h_v);
2475 
2476  // Assign to LayoutStride view
2477  StrideViewType vs = v;
2478 
2479  // Copy back
2480  // Note: don't create h_vs through create_mirror_view and deep_copy
2481  // since Kokkos doesn't support deep_copy of non-contiguous views
2482  //stride_host_view_type h_vs = Kokkos::create_mirror_view(vs);
2483  //Kokkos::deep_copy(h_vs, vs);
2484  stride_host_view_type h_vs = h_v;
2485 
2486  // Check
2487  success = true;
2488  TEUCHOS_TEST_EQUALITY(h_vs.extent(0), num_rows, out, success);
2489  TEUCHOS_TEST_EQUALITY(h_vs.extent(1), num_cols, out, success);
2490  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(h_vs), fad_size+1, out, success);
2491  for (size_type i=0; i<num_rows; ++i) {
2492  for (size_type j=0; j<num_cols; ++j) {
2493  FadType f = generate_fad<FadType>(num_rows, num_cols, fad_size, i, j);
2494  success = success && checkFads(f, h_vs(i,j), out);
2495  }
2496  }
2497 }
2498 
2500  Kokkos_View_Fad, CommonViewAllocMixedSpec, FadType, Layout, Device )
2501 {
2502  typedef Kokkos::View<FadType**,Kokkos::LayoutContiguous<Layout>,Device> ContViewType;
2503  typedef Kokkos::View<FadType**,Layout,Device> ViewType;
2504  typedef typename ContViewType::size_type size_type;
2505 
2506  const size_type num_rows = global_num_rows;
2507  const size_type num_cols = global_num_cols;
2508  const size_type fad_size = global_fad_size;
2509 
2510  // Create contiguous view
2511  ContViewType v1;
2512 #if defined (SACADO_DISABLE_FAD_VIEW_SPEC)
2513  v1 = ContViewType ("view", num_rows, num_cols);
2514 #else
2515  v1 = ContViewType ("view", num_rows, num_cols, fad_size+1);
2516 #endif
2517 
2518  // Create non-contiguous view using commen_view_alloc_prop
2519  auto cprop = Kokkos::common_view_alloc_prop(v1);
2520  ViewType v2(Kokkos::view_alloc("v2",cprop), num_rows, num_cols);
2521 
2522  // Check dimensions are correct for v2
2523  success = true;
2524  TEUCHOS_TEST_EQUALITY(v2.extent(0), num_rows, out, success);
2525  TEUCHOS_TEST_EQUALITY(v2.extent(1), num_cols, out, success);
2526  TEUCHOS_TEST_EQUALITY(Kokkos::dimension_scalar(v2), fad_size+1, out, success);
2527 }
2528 
2529 #else
2530 
2532  Kokkos_View_Fad, ShmemSize, FadType, Layout, Device )
2533 {
2534  typedef Kokkos::View<FadType**,Layout,Device> ViewType;
2535  typedef typename ViewType::size_type size_type;
2536 
2537  const size_type num_rows = global_num_rows;
2538  const size_type num_cols = global_num_cols;
2539 
2540  // Compute shared memory size for View
2541  const size_type shmem_size =
2542  ViewType::shmem_size(num_rows, num_cols);
2543 
2544  // Check
2545 #if !defined(KOKKOS_ENABLE_IMPL_VIEW_LEGACY)
2546  size_t scratch_value_alignment =
2547  Kokkos::max({sizeof(FadType),
2548  alignof(FadType),
2549  static_cast<size_t>(
2550  ViewType::execution_space::scratch_memory_space::ALIGN)});
2551  const size_type shmem_size_expected =
2552  sizeof(FadType) * global_num_rows * global_num_cols + scratch_value_alignment;
2553 #else
2554  static const size_type align = 8;
2555  static const size_type mask = align - 1;
2556  const size_type shmem_size_expected =
2557  (( sizeof(FadType) * global_num_rows * global_num_cols + mask ) & ~mask) + sizeof(typename ViewType::traits::value_type);
2558 #endif
2559  TEUCHOS_TEST_EQUALITY(shmem_size, shmem_size_expected, out, success);
2560 }
2561 
2563  Kokkos_View_Fad, Unmanaged, FadType, Layout, Device ) {}
2564 
2566  Kokkos_View_Fad, Unmanaged2, FadType, Layout, Device ) {}
2567 
2569  Kokkos_View_Fad, UnmanagedConst, FadType, Layout, Device ) {}
2570 
2572  Kokkos_View_Fad, UnmanagedConst2, FadType, Layout, Device ) {}
2573 
2575  Kokkos_View_Fad, SFadNoSizeArg, FadType, Layout, Device ) {}
2576 
2578  Kokkos_View_Fad, Partition, FadType, Layout, Device ) {}
2579 
2581  Kokkos_View_Fad, AssignLayoutContiguousToLayoutStride, FadType, Layout, Device ) {}
2582 
2584  Kokkos_View_Fad, CommonViewAllocMixedSpec, FadType, Layout, Device ) {}
2585 
2586 #endif
2587 
2588 #define VIEW_FAD_TESTS_FLD( F, L, D ) \
2589  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, Size, F, L, D ) \
2590  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, DeepCopy, F, L, D ) \
2591  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, DeepCopy_ConstantScalar, F, L, D ) \
2592  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, DeepCopy_ConstantZero, F, L, D ) \
2593  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, DeepCopy_ConstantFad, F, L, D ) \
2594  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, DeepCopy_ConstantFadFull, F, L, D ) \
2595  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, LocalDeepCopy, F, L, D ) \
2596  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, LocalDeepCopyTeam, F, L, D ) \
2597  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, ScalarAssign, F, L, D ) \
2598  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, ValueAssign, F, L, D ) \
2599  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, Resize, F, L, D ) \
2600  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, Unmanaged, F, L, D ) \
2601  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, Unmanaged2, F, L, D ) \
2602  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, UnmanagedConst, F, L, D ) \
2603  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, UnmanagedConst2, F, L, D ) \
2604  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, Multiply, F, L, D ) \
2605  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, MultiplyUpdate, F, L, D ) \
2606  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, MultiplyConst, F, L, D ) \
2607  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, MultiplyMixed, F, L, D ) \
2608  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, Rank8, F, L, D ) \
2609  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, Roger, F, L, D ) \
2610  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, AtomicAdd, F, L, D ) \
2611  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, AssignDifferentStrides, F, L, D ) \
2612  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, ScalarValue, F, L, D ) \
2613  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, DynRankDimensionScalar, F, L, D ) \
2614  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, DynRankAssignStatic0, F, L, D ) \
2615  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, DynRankAssignStatic1, F, L, D ) \
2616  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, DynRankAssignStatic2, F, L, D ) \
2617  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, DynRankMultiply, F, L, D ) \
2618  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, SubdynrankviewCol, F, L, D ) \
2619  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, SubdynrankviewRow, F, L, D ) \
2620  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, SubdynrankviewScalar, F, L, D ) \
2621  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, Subview, F, L, D ) \
2622  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, Subview2, F, L, D ) \
2623  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, ShmemSize, F, L, D ) \
2624  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, ConstViewAssign, F, L, D )
2625 
2626 #define VIEW_FAD_TESTS_SFLD( F, L, D ) \
2627  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, SFadNoSizeArg, F, L, D )
2628 
2629 #define VIEW_FAD_TESTS_FDI( F, D ) \
2630  using Kokkos::LayoutLeft; \
2631  using Kokkos::LayoutRight; \
2632  VIEW_FAD_TESTS_FLD( F, LayoutLeft, D ) \
2633  VIEW_FAD_TESTS_FLD( F, LayoutRight, D ) \
2634  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, AssignLayoutContiguousToLayoutStride, F, LayoutLeft, D ) \
2635  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, AssignLayoutContiguousToLayoutStride, F, LayoutRight, D ) \
2636  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, CommonViewAllocMixedSpec, F, LayoutLeft, D ) \
2637  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, CommonViewAllocMixedSpec, F, LayoutRight, D )
2638 
2639 #define VIEW_FAD_TESTS_SFDI( F, D ) \
2640  using Kokkos::LayoutLeft; \
2641  using Kokkos::LayoutRight; \
2642  VIEW_FAD_TESTS_SFLD( F, LayoutLeft, D ) \
2643  VIEW_FAD_TESTS_SFLD( F, LayoutRight, D )
2644 
2645 #if defined(HAVE_SACADO_VIEW_SPEC) && !defined(SACADO_DISABLE_FAD_VIEW_SPEC)
2648 
2649 #ifndef SACADO_HAS_NEW_KOKKOS_VIEW_IMPL
2650 #define VIEW_FAD_TESTS_FDC( F, D ) \
2651  VIEW_FAD_TESTS_FLD( F, LeftContiguous, D ) \
2652  VIEW_FAD_TESTS_FLD( F, RightContiguous, D ) \
2653  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, Partition, F, LeftContiguous, D ) \
2654  TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( Kokkos_View_Fad, Partition, F, RightContiguous, D )
2655 #else
2656 #define VIEW_FAD_TESTS_FDC( F, D ) \
2657  VIEW_FAD_TESTS_FLD( F, LeftContiguous, D ) \
2658  VIEW_FAD_TESTS_FLD( F, RightContiguous, D )
2659 #endif
2660 
2661 #define VIEW_FAD_TESTS_SFDC( F, D ) \
2662  VIEW_FAD_TESTS_SFLD( F, LeftContiguous, D ) \
2663  VIEW_FAD_TESTS_SFLD( F, RightContiguous, D )
2664 #else
2665 #define VIEW_FAD_TESTS_FDC( F, D ) /* */
2666 #define VIEW_FAD_TESTS_SFDC( F, D ) /* */
2667 #endif
2668 
2669 #define VIEW_FAD_TESTS_FD( F, D ) \
2670  VIEW_FAD_TESTS_FDI( F, D ) \
2671  VIEW_FAD_TESTS_FDC( F, D )
2672 
2673 #define VIEW_FAD_TESTS_SFD( F, D ) \
2674  VIEW_FAD_TESTS_SFDI( F, D ) \
2675  VIEW_FAD_TESTS_SFDC( F, D )
2676 
2677 // We've unified the implementation for the different Fad variants, so
2678 // there is no reason to test ELRFad, CacheFad, and ELRCacheFad.
2682 
2683 /*
2684 typedef Sacado::ELRFad::DFad<double> ELRDFadType;
2685 typedef Sacado::ELRFad::SLFad<double,2*global_fad_size> ELRSLFadType;
2686 typedef Sacado::ELRFad::SFad<double,global_fad_size> ELRSFadType;
2687 
2688 typedef Sacado::CacheFad::DFad<double> CacheDFadType;
2689 typedef Sacado::CacheFad::SLFad<double,2*global_fad_size> CacheSLFadType;
2690 typedef Sacado::CacheFad::SFad<double,global_fad_size> CacheSFadType;
2691 
2692 typedef Sacado::ELRCacheFad::DFad<double> ELRCacheDFadType;
2693 typedef Sacado::ELRCacheFad::SLFad<double,2*global_fad_size> ELRCacheSLFadType;
2694 typedef Sacado::ELRCacheFad::SFad<double,global_fad_size> ELRCacheSFadType;
2695 */
2696 
2697 // We can't use DFad unless we use the View specialization
2698 #if defined(HAVE_SACADO_VIEW_SPEC) && !defined(SACADO_DISABLE_FAD_VIEW_SPEC) && SACADO_TEST_DFAD
2699 #define VIEW_FAD_TESTS_D( D ) \
2700  VIEW_FAD_TESTS_FD( SFadType, D ) \
2701  VIEW_FAD_TESTS_FD( SLFadType, D ) \
2702  VIEW_FAD_TESTS_FD( DFadType, D ) \
2703  VIEW_FAD_TESTS_SFD( SFadType, D )
2704 
2705 #if 0
2706  VIEW_FAD_TESTS_FD( ELRSFadType, D ) \
2707  VIEW_FAD_TESTS_FD( ELRSLFadType, D ) \
2708  VIEW_FAD_TESTS_FD( ELRDFadType, D ) \
2709  VIEW_FAD_TESTS_FD( CacheSFadType, D ) \
2710  VIEW_FAD_TESTS_FD( CacheSLFadType, D ) \
2711  VIEW_FAD_TESTS_FD( CacheDFadType, D ) \
2712  VIEW_FAD_TESTS_FD( ELRCacheSFadType, D ) \
2713  VIEW_FAD_TESTS_FD( ELRCacheSLFadType, D ) \
2714  VIEW_FAD_TESTS_FD( ELRCacheDFadType, D ) \
2715  VIEW_FAD_TESTS_SFD( SFadType, D ) \
2716  VIEW_FAD_TESTS_SFD( ELRSFadType, D ) \
2717  VIEW_FAD_TESTS_SFD( CacheSFadType, D ) \
2718  VIEW_FAD_TESTS_SFD( ELRCacheSFadType, D )
2719 #endif
2720 
2721 #else
2722 
2723 #define VIEW_FAD_TESTS_D( D ) \
2724  VIEW_FAD_TESTS_FD( SFadType, D ) \
2725  VIEW_FAD_TESTS_FD( SLFadType, D ) \
2726  VIEW_FAD_TESTS_SFD( SFadType, D )
2727 
2728 #if 0
2729  VIEW_FAD_TESTS_FD( ELRSFadType, D ) \
2730  VIEW_FAD_TESTS_FD( ELRSLFadType, D ) \
2731  VIEW_FAD_TESTS_FD( CacheSFadType, D ) \
2732  VIEW_FAD_TESTS_FD( CacheSLFadType, D ) \
2733  VIEW_FAD_TESTS_FD( ELRCacheSFadType, D ) \
2734  VIEW_FAD_TESTS_FD( ELRCacheSLFadType, D ) \
2735  VIEW_FAD_TESTS_SFD( SFadType, D ) \
2736  VIEW_FAD_TESTS_SFD( ELRSFadType, D ) \
2737  VIEW_FAD_TESTS_SFD( CacheSFadType, D ) \
2738  VIEW_FAD_TESTS_SFD( ELRCacheSFadType, D )
2739 #endif
2740 
2741 #endif
SACADO_INLINE_FUNCTION ScalarType< T >::type scalarValue(const T &x)
A simple template function for invoking ScalarValue&lt;&gt;
InputViewType::size_type size_type
team_policy_type::member_type team_handle
Kokkos::TeamPolicy< execution_space > team_policy_type
static const size_type stride
ViewType::value_type::value_type ScalarType
ViewType::execution_space execution_space
Kokkos::TeamPolicy< execution_space > team_policy_type
static const size_type stride
ViewType::execution_space execution_space
KOKKOS_INLINE_FUNCTION void operator()(const team_handle &team) const
void f()
ViewType::size_type size_type
const ViewType m_v
Kokkos::LayoutContiguous< Kokkos::LayoutRight > RightContiguous
Kokkos::RangePolicy< execution_space > range_policy_type
ViewType::value_type ValueType
Kokkos::LayoutContiguous< Kokkos::LayoutLeft > LeftContiguous
static const size_type stride
#define TEUCHOS_TEST_FLOATING_EQUALITY(v1, v2, tol, out, success)
KOKKOS_INLINE_FUNCTION void operator()(const team_handle &team) const
static const bool value
const int global_fad_size
KOKKOS_INLINE_FUNCTION void operator()(const size_type i) const
KOKKOS_INLINE_FUNCTION auto subview(const View< D, Kokkos::LayoutContiguous< LayoutSrc, StrideSrc >, P...> &src, Args...args)
const ScalarViewType m_s
#define VIEW_FAD_TESTS_SFD(F, D)
KOKKOS_INLINE_FUNCTION void operator()(const size_type i) const
Sacado::Fad::DFad< double > FadType
KOKKOS_INLINE_FUNCTION void operator()(const size_type i) const
#define VIEW_FAD_TESTS_FD(F, D)
expr true
const InputViewType1 m_v1
static const size_type stride
ViewType::size_type size_type
InputViewType::execution_space execution_space
static void apply(const InputViewType v1, const OutputViewType v2, const size_type col)
team_policy_type::member_type team_handle
bool checkFads(const FadType1 &x, const FadType2 &x2, Teuchos::FancyOStream &out, double tol=1.0e-15)
const ScalarViewType m_s
const int global_num_rows
Kokkos::RangePolicy< execution_space > range_policy_type
KOKKOS_INLINE_FUNCTION void operator()(const size_type i) const
ViewType::execution_space execution_space
const ScalarType m_s
const ViewType m_v
scalar generate_fad(const size_t n0, const size_t n1, const size_t n2, const size_t n3, const int fad_size, const size_t i0, const size_t i1, const size_t i2, const size_t i3, const int i_fad)
Sacado::Fad::SFad< double, fad_dim > SFadType
const OutputViewType m_v2
expr val()
TEUCHOS_UNIT_TEST_TEMPLATE_3_DECL(Kokkos_View_FadFad, DeepCopy, FadFadType, Layout, Device)
AssignRank2Rank1Kernel(const InputViewType v1, const OutputViewType v2, const size_type col)
expr expr1 expr1 expr1 c expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr1 c expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr1 c *expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr1 c expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr2 expr1 expr2 expr1 expr1 expr1 c
team_policy_type::member_type team_handle
const int global_num_cols
KOKKOS_INLINE_FUNCTION void operator()(const team_handle &team) const
ValueAssignKernel(const ViewType &v, const ScalarViewType &s)
ViewType::size_type size_type
GeneralFad< DynamicStorage< T > > DFad
team_policy_type::member_type team_handle
const bool m_update
static void apply(const InputViewType1 v1, const InputViewType2 v2, const OutputViewType v3, const bool update=false)
InputViewType1::execution_space execution_space
#define GLOBAL_FAD_SIZE
int value
Kokkos::TeamPolicy< execution_space > team_policy_type
static void apply(const ViewType &v, const ScalarViewType &s)
Sacado::Fad::SLFad< double, fad_dim > SLFadType
void g()
Sacado::Fad::DFad< double > DFadType
static void apply(const ViewType &v, const ScalarViewType &s)
Kokkos::TeamPolicy< execution_space > team_policy_type
ScalarAssignKernel(const ViewType &v, const ScalarType &s)
KOKKOS_INLINE_FUNCTION void operator()(const team_handle &team) const
static const bool value
AtomicAddKernel(const ViewType &v, const ScalarViewType &s)
KOKKOS_INLINE_FUNCTION auto subdynrankview(const DynRankView< T, LayoutContiguous< LayoutSrc, StrideSrc >, DRVArgs...> &drv, SubArg0 arg0=SubArg0{}, SubArg1 arg1=SubArg1{}, SubArg2 arg2=SubArg2{}, SubArg3 arg3=SubArg3{}, SubArg4 arg4=SubArg4{}, SubArg5 arg5=SubArg5{}, SubArg6 arg6=SubArg6{})
expr expr expr fastAccessDx(i)) FAD_UNARYOP_MACRO(exp
#define TEUCHOS_TEST_EQUALITY(v1, v2, out, success)
KOKKOS_INLINE_FUNCTION void operator()(const team_handle &team) const
team_policy_type::member_type team_handle
Kokkos::ThreadLocalScalarType< ViewType >::type local_scalar_type
const double tol
SimpleFad< ValueT > max(const SimpleFad< ValueT > &a, const SimpleFad< ValueT > &b)
static void apply(const ViewType &v, const ScalarType &s)
Kokkos::ThreadLocalScalarType< ViewType >::type local_scalar_type
MultiplyKernel(const InputViewType1 v1, const InputViewType2 v2, const OutputViewType v3, const bool update)
const InputViewType m_v1
InputViewType1::size_type size_type
Kokkos::RangePolicy< execution_space > range_policy_type
Kokkos::RangePolicy< execution_space > range_policy_type
KOKKOS_INLINE_FUNCTION void operator()(const size_type i) const
GeneralFad< StaticFixedStorage< T, Num > > SFad
const InputViewType2 m_v2
Kokkos::RangePolicy< execution_space > range_policy_type
Kokkos::TeamPolicy< execution_space > team_policy_type
const OutputViewType m_v3