Sacado Package Browser (Single Doxygen Collection)  Version of the Day
ViewFactoryTests.cpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Sacado Package
5 // Copyright (2006) Sandia Corporation
6 //
7 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8 // the U.S. Government retains certain rights in this software.
9 //
10 // This library is free software; you can redistribute it and/or modify
11 // it under the terms of the GNU Lesser General Public License as
12 // published by the Free Software Foundation; either version 2.1 of the
13 // License, or (at your option) any later version.
14 //
15 // This library is distributed in the hope that it will be useful, but
16 // WITHOUT ANY WARRANTY; without even the implied warranty of
17 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18 // Lesser General Public License for more details.
19 //
20 // You should have received a copy of the GNU Lesser General Public
21 // License along with this library; if not, write to the Free Software
22 // Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301
23 // USA
24 // Questions? Contact David M. Gay (dmgay@sandia.gov) or Eric T. Phipps
25 // (etphipp@sandia.gov).
26 //
27 // ***********************************************************************
28 // @HEADER
29 #include "Teuchos_UnitTestHarness.hpp"
30 #include "Teuchos_UnitTestRepository.hpp"
31 #include "Teuchos_GlobalMPISession.hpp"
32 
33 #include "Sacado.hpp"
34 
35 #include "Kokkos_ViewFactory.hpp"
36 
37 TEUCHOS_UNIT_TEST(view_factory, dyn_rank_views)
38 {
40  using Kokkos::View;
41  using Kokkos::DynRankView;
45  using Kokkos::dimension_scalar;
46  using Kokkos::view_alloc;
47  using Kokkos::WithoutInitializing;
48  const unsigned derivative_dim_plus_one = 7;
49 
50 
51  // Test constructing a View pod using Kokkos::view_alloc with deduce_value_type
52  {
53  // Typedef View
54  typedef View<double**, Kokkos::DefaultExecutionSpace> view_type;
55 
56  // Create two rank 2 Views that will be used for deducing types and Fad dims
57  view_type v1("v1", 10, 4);
58  view_type v2("v2", 10, 4);
59 
60  // Get common type of the Views
61  using CommonValueType = typename decltype( Kokkos::common_view_alloc_prop( v1, v2 ) )::value_type;
62  using ScalarArrayType = typename decltype( Kokkos::common_view_alloc_prop( v1, v2 ) )::scalar_array_type;
63  // Create an instance of this returned type to pass to ViewCtorProp via view_alloc function
64  auto cvt_for_ctorprop = Kokkos::common_view_alloc_prop(v1, v2);
65 
66  // Create a view with the common type and the max fad_dim of the views passed to deduce_value_type
67  typedef View< CommonValueType** > ViewCommonType;
68  ViewCommonType vct1( Kokkos::view_alloc("vct1", cvt_for_ctorprop), 10, 4 ); // fad_dim deduced and comes from the cvt_for_ctorprop
69 
70  TEST_EQUALITY(vct1.extent(0), v1.extent(0));
71  TEST_EQUALITY(vct1.extent(1), v1.extent(1));
72  TEST_EQUALITY(vct1.extent(2), v1.extent(2));
73  TEST_EQUALITY( Kokkos::dimension_scalar(vct1), 0);
75  bool check_eq_scalar_double = std::is_same < double, ScalarArrayType >::value;
76  TEST_EQUALITY(check_eq_kokkos_type, true);
77  TEST_EQUALITY(check_eq_scalar_double, true);
78  }
79  // Test constructing a View of Fad using Kokkos::view_alloc with deduce_value_type
80  {
81  // Typedef View
82  typedef View<FadType**, Kokkos::DefaultExecutionSpace> view_type;
83 
84  // Create two rank 2 Views that will be used for deducing types and Fad dims
85  view_type v1("v1", 10, 4, derivative_dim_plus_one );
86  view_type v2("v2", 10, 4, derivative_dim_plus_one );
87 
88  // Get common type of the Views
89  using CommonValueType = typename decltype( Kokkos::common_view_alloc_prop( v1, v2 ) )::value_type;
90  using ScalarArrayType = typename decltype( Kokkos::common_view_alloc_prop( v1, v2 ) )::scalar_array_type;
91  // Create an instance of this returned type to pass to ViewCtorProp via view_alloc function
92  auto cvt_for_ctorprop = Kokkos::common_view_alloc_prop(v1, v2);
93 
94  // Create a view with the common type and the max fad_dim of the views passed to deduce_value_type
95  typedef View< CommonValueType** > ViewCommonType;
96  ViewCommonType vct1( Kokkos::view_alloc("vct1", cvt_for_ctorprop), 10, 4 ); // fad_dim deduced and comes from the cvt_for_ctorprop
97 
98  TEST_EQUALITY(dimension_scalar(vct1), derivative_dim_plus_one);
99  TEST_EQUALITY(vct1.extent(0), v1.extent(0));
100  TEST_EQUALITY(vct1.extent(1), v1.extent(1));
101  TEST_EQUALITY(vct1.extent(2), v1.extent(2));
103  bool check_eq_fad_type = std::is_same < CommonValueType, FadType >::value;
104  bool check_eq_scalar_double = std::is_same < double, ScalarArrayType >::value;
105  TEST_EQUALITY(check_neq_kokkos_type, false);
106  TEST_EQUALITY(check_eq_fad_type, true);
107  TEST_EQUALITY(check_eq_scalar_double, true);
108  }
109  // Test constructing a View from mix of View and Viewof Fads using Kokkos::view_alloc with deduce_value_type
110  {
111  // Typedef View
112  typedef View<FadType**, Kokkos::DefaultExecutionSpace> view_of_fad_type;
113  typedef View<double**, Kokkos::DefaultExecutionSpace> view_of_pod_type;
114 
115  // Create two rank 2 Views that will be used for deducing types and Fad dims
116  view_of_fad_type v1("v1", 10, 4, derivative_dim_plus_one );
117  view_of_pod_type v2("v2", 10, 4);
118 
119  // Get common type of the Views
120  using CommonValueType = typename decltype( Kokkos::common_view_alloc_prop( v1, v2 ) )::value_type;
121  using ScalarArrayType = typename decltype( Kokkos::common_view_alloc_prop( v1, v2 ) )::scalar_array_type;
122  // Create an instance of this returned type to pass to ViewCtorProp via view_alloc function
123  auto cvt_for_ctorprop = Kokkos::common_view_alloc_prop(v1, v2);
124 
125  // Create a view with the common type and the max fad_dim of the views passed to deduce_value_type
126  typedef View< CommonValueType** > ViewCommonType;
127  ViewCommonType vct1( Kokkos::view_alloc("vct1", cvt_for_ctorprop), 10, 4 ); // fad_dim deduced and comes from the cvt_for_ctorprop
128 
129  TEST_EQUALITY(dimension_scalar(vct1), derivative_dim_plus_one);
130  TEST_EQUALITY(vct1.extent(0), v1.extent(0));
131  TEST_EQUALITY(vct1.extent(1), v1.extent(1));
132  TEST_EQUALITY(vct1.extent(2), v1.extent(2));
134  bool check_eq_fad_type = std::is_same < CommonValueType, FadType >::value;
135  bool check_eq_scalar_double = std::is_same < double, ScalarArrayType >::value;
136  TEST_EQUALITY(check_neq_kokkos_type, false);
137  TEST_EQUALITY(check_eq_fad_type, true);
138  TEST_EQUALITY(check_eq_scalar_double, true);
139  }
140  // Test constructing a DynRankView using Kokkos::view_alloc with deduce_value_type
141  {
142  // Typedef View
143  typedef DynRankView<FadType, Kokkos::DefaultExecutionSpace> view_type;
144 
145  // Create two rank 2 Views that will be used for deducing types and Fad dims
146  view_type v1("v1", 10, 4, derivative_dim_plus_one );
147  view_type v2("v2", 10, 4, derivative_dim_plus_one );
148 
149  // Get common type of the Views
150  using CommonValueType = typename decltype( Kokkos::common_view_alloc_prop( v1, v2 ) )::value_type;
151  using ScalarArrayType = typename decltype( Kokkos::common_view_alloc_prop( v1, v2 ) )::scalar_array_type;
152  // Create an instance of this returned type to pass to ViewCtorProp via view_alloc function
153  auto cvt_for_ctorprop = Kokkos::common_view_alloc_prop(v1, v2);
154 
155  // Create a view with the common type and the max fad_dim of the views passed to deduce_value_type
156  typedef DynRankView< CommonValueType > ViewCommonType;
157  ViewCommonType vct1( Kokkos::view_alloc("vct1", cvt_for_ctorprop), 10, 4 ); // fad_dim deduced and comes from the cvt_for_ctorprop
158 
159  TEST_EQUALITY(dimension_scalar(vct1), derivative_dim_plus_one);
160  TEST_EQUALITY(vct1.extent(0), v1.extent(0));
161  TEST_EQUALITY(vct1.extent(1), v1.extent(1));
162  TEST_EQUALITY(vct1.extent(2), v1.extent(2));
163  TEST_EQUALITY(Kokkos::rank(vct1), 2);
165  bool check_eq_fad_type = std::is_same < CommonValueType, FadType >::value;
166  bool check_eq_scalar_double = std::is_same < double, ScalarArrayType >::value;
167  TEST_EQUALITY(check_neq_kokkos_type, false);
168  TEST_EQUALITY(check_eq_fad_type, true);
169  TEST_EQUALITY(check_eq_scalar_double, true);
170  }
171  // Test constructing a DynRankView from mix of DynRankView and DynRankView of Fads using Kokkos::view_alloc with deduce_value_type
172  {
173  // Typedef View
174  typedef DynRankView<FadType, Kokkos::DefaultExecutionSpace> view_of_fad_type;
175  typedef DynRankView<double, Kokkos::DefaultExecutionSpace> view_of_pod_type;
176 
177  // Create two rank 2 Views that will be used for deducing types and Fad dims
178  view_of_fad_type v1("v1", 10, 4, derivative_dim_plus_one );
179  view_of_pod_type v2("v2", 10, 4);
180 
181  // Get common type of the Views
182  using CommonValueType = typename decltype( Kokkos::common_view_alloc_prop( v1, v2 ) )::value_type;
183  using ScalarArrayType = typename decltype( Kokkos::common_view_alloc_prop( v1, v2 ) )::scalar_array_type;
184  // Create an instance of this returned type to pass to ViewCtorProp via view_alloc function
185  auto cvt_for_ctorprop = Kokkos::common_view_alloc_prop(v1, v2);
186 
187  // Create a view with the common type and the max fad_dim of the views passed to deduce_value_type
188  typedef DynRankView< CommonValueType > ViewCommonType;
189  ViewCommonType vct1( Kokkos::view_alloc("vct1", cvt_for_ctorprop), 10, 4 ); // fad_dim deduced and comes from the cvt_for_ctorprop
190 
191  TEST_EQUALITY(dimension_scalar(vct1), derivative_dim_plus_one);
192  TEST_EQUALITY(vct1.extent(0), v1.extent(0));
193  TEST_EQUALITY(vct1.extent(1), v1.extent(1));
194  TEST_EQUALITY(vct1.extent(2), v1.extent(2));
195  TEST_EQUALITY(Kokkos::rank(vct1), 2);
197  bool check_eq_fad_type = std::is_same < CommonValueType, FadType >::value;
198  bool check_eq_scalar_double = std::is_same < double, ScalarArrayType >::value;
199  TEST_EQUALITY(check_neq_kokkos_type, false);
200  TEST_EQUALITY(check_eq_fad_type, true);
201  TEST_EQUALITY(check_eq_scalar_double, true);
202  }
203 
204 
205 
206 
207  // Test a DynRankView from a DynRankView
208  {
209  DynRankView<FadType> a("a",10,4,13,derivative_dim_plus_one);
210  TEST_EQUALITY(dimension_scalar(a),derivative_dim_plus_one);
211  TEST_EQUALITY(a.rank(),3);
212 
213  auto b = createDynRankView(a,"b",5,3,8);
214  TEST_EQUALITY(dimension_scalar(b),derivative_dim_plus_one);
215  TEST_EQUALITY(b.rank(),3);
216 
217  auto c = createDynRankView(a,view_alloc("c",WithoutInitializing),5,3,8);
218  TEST_EQUALITY(dimension_scalar(c),derivative_dim_plus_one);
219  TEST_EQUALITY(c.rank(),3);
220 
221  using d_type = Kokkos::DynRankView<FadType,Kokkos::LayoutRight>;
222  d_type d = createDynRankViewWithType<d_type>(a,"d",5,3,8);
223  TEST_EQUALITY(dimension_scalar(d),derivative_dim_plus_one);
224  TEST_EQUALITY(d.rank(),3);
225  }
226 
227  // Test a DynRankView from a View
228  {
229  View<FadType*> a("a",8,derivative_dim_plus_one);
230  TEST_EQUALITY(dimension_scalar(a),derivative_dim_plus_one);
231 
232  auto b = createDynRankView(a,"b",5,3,8);
233  TEST_EQUALITY(dimension_scalar(b),derivative_dim_plus_one);
234  TEST_EQUALITY(b.rank(),3);
235 
236  auto c = createDynRankView(a,view_alloc("c",WithoutInitializing),5,3,8);
237  TEST_EQUALITY(dimension_scalar(c),derivative_dim_plus_one);
238  TEST_EQUALITY(c.rank(),3);
239 
240  using d_type = Kokkos::DynRankView<FadType,Kokkos::LayoutRight>;
241  d_type d = createDynRankViewWithType<d_type>(a,"d",5,3,8);
242  TEST_EQUALITY(dimension_scalar(d),derivative_dim_plus_one);
243  TEST_EQUALITY(d.rank(),3);
244  }
245 
246  // Test a View from a View
247  {
248  View<FadType*> a("a",8,derivative_dim_plus_one);
249  TEST_EQUALITY(dimension_scalar(a),derivative_dim_plus_one);
250 
251  using b_type = Kokkos::View<FadType***>;
252  b_type b = createViewWithType<b_type>(a,"b",5,3,8);
253  TEST_EQUALITY(dimension_scalar(b),derivative_dim_plus_one);
254 
255  b_type c = createViewWithType<b_type>(a,view_alloc("c",WithoutInitializing),5,3,8);
256  TEST_EQUALITY(dimension_scalar(c),derivative_dim_plus_one);
257 
258  using d_type = Kokkos::View<FadType***,Kokkos::LayoutRight>;
259  d_type d = createViewWithType<d_type>(a,"d",5,3,8);
260  TEST_EQUALITY(dimension_scalar(d),derivative_dim_plus_one);
261  }
262 
263  // Test a View from a DynRankView
264  {
265  DynRankView<FadType> a("a",10,4,13,derivative_dim_plus_one);
266  TEST_EQUALITY(dimension_scalar(a),derivative_dim_plus_one);
267  TEST_EQUALITY(a.rank(),3);
268 
269  using b_type = Kokkos::View<FadType***>;
270  b_type b = createViewWithType<b_type>(a,"b",5,3,8);
271  TEST_EQUALITY(dimension_scalar(b),derivative_dim_plus_one);
272 
273  b_type c = createViewWithType<b_type>(a,view_alloc("c",WithoutInitializing),5,3,8);
274  TEST_EQUALITY(dimension_scalar(c),derivative_dim_plus_one);
275 
276  using d_type = Kokkos::View<FadType***,Kokkos::LayoutRight>;
277  d_type d = createViewWithType<d_type>(a,"d",5,3,8);
278  TEST_EQUALITY(dimension_scalar(d),derivative_dim_plus_one);
279  }
280 
281  // Test creation of a Fad DynRankView from a double DynRankView
282  {
283  DynRankView<double> a("a",10,4,13);
284  TEST_EQUALITY(dimension_scalar(a),0);
285  TEST_EQUALITY(a.rank(),3);
286 
287  using b_type = Kokkos::DynRankView<FadType,Kokkos::LayoutRight>;
288  b_type b = createDynRankViewWithType<b_type>(a,"b",5,3,8);
289  TEST_EQUALITY(dimension_scalar(b),1);
290  TEST_EQUALITY(b.rank(),3);
291  }
292 
293  // Test a double DynRankView from a double DynRankView
294  {
295  DynRankView<double> a("a",10,4,13);
296  TEST_EQUALITY(dimension_scalar(a),0);
297  TEST_EQUALITY(a.rank(),3);
298 
299  auto b = createDynRankView(a,"b",5,3,8);
300  TEST_EQUALITY(dimension_scalar(b),0);
301  TEST_EQUALITY(b.rank(),3);
302  }
303 
304  // Test double rank 0
305  {
306  DynRankView<double> a("a",10,4,13);
307  TEST_EQUALITY(dimension_scalar(a),0);
308  TEST_EQUALITY(a.rank(),3);
309 
310  auto b = createDynRankView(a,"b");
311  TEST_EQUALITY(dimension_scalar(b),0);
312  TEST_EQUALITY(b.rank(),0);
313  }
314 
315  // Test Fad rank 0
316  {
317  DynRankView<FadType> a("a",10,4,13,derivative_dim_plus_one);
318  TEST_EQUALITY(dimension_scalar(a),derivative_dim_plus_one);
319  TEST_EQUALITY(a.rank(),3);
320 
321  auto b = createDynRankView(a,"b");
322  TEST_EQUALITY(dimension_scalar(b),derivative_dim_plus_one);
323  TEST_EQUALITY(b.rank(),0);
324  }
325 
326  // Test unmanaged view of double
327  {
328  Kokkos::View<double*> a("a",5*3);
329  using b_type = Kokkos::View<double**,Kokkos::MemoryUnmanaged>;
330  b_type b = createViewWithType<b_type>(a,a.data(),5,3);
331  TEST_EQUALITY(b.extent(0),5);
332  TEST_EQUALITY(b.extent(1),3);
333  TEST_EQUALITY(dimension_scalar(b),0);
334  }
335 
336  // Test unmanaged view of Fad
337  {
338  Kokkos::View<FadType*> a("a",5*3,derivative_dim_plus_one);
339  using b_type = Kokkos::View<FadType**,Kokkos::MemoryUnmanaged>;
340  b_type b = createViewWithType<b_type>(a,a.data(),5,3);
341  TEST_EQUALITY(b.extent(0),5);
342  TEST_EQUALITY(b.extent(1),3);
343  TEST_EQUALITY(dimension_scalar(b),derivative_dim_plus_one);
344  }
345 
346  // Test LayoutStride view of double
347  {
348  Kokkos::DynRankView<double> a("a",10,13);
349  auto b = Kokkos::subview(a, std::make_pair(4,8), std::make_pair(5,11));
350  auto c = createDynRankView(b,"c",5,3);
351  using b_type = decltype(b);
352  using c_type = decltype(c);
353  using b_layout = typename b_type::array_layout;
354  using c_layout = typename c_type::array_layout;
355  using default_layout = typename b_type::device_type::execution_space::array_layout;
356  const bool is_b_layout_stride =
358  const bool is_c_default_layout =
360  TEST_EQUALITY(is_b_layout_stride,true);
361  TEST_EQUALITY(is_c_default_layout,true);
362  TEST_EQUALITY(c.rank(),2);
363  TEST_EQUALITY(c.extent(0),5);
364  TEST_EQUALITY(c.extent(1),3);
365  TEST_EQUALITY(dimension_scalar(b),0);
366  }
367 
368  // Test LayoutStride view of Fad
369  {
370  Kokkos::DynRankView<FadType> a("a",10,13,derivative_dim_plus_one);
371  auto b = Kokkos::subview(a, std::make_pair(4,8), std::make_pair(5,11));
372  auto c = createDynRankView(b,"c",5,3);
373  using b_type = decltype(b);
374  using c_type = decltype(c);
375  using b_layout = typename b_type::array_layout;
376  using c_layout = typename c_type::array_layout;
377  using default_layout = typename b_type::device_type::execution_space::array_layout;
378  const bool is_b_layout_stride =
380  const bool is_c_default_layout =
382  TEST_EQUALITY(is_b_layout_stride,true);
383  TEST_EQUALITY(is_c_default_layout,true);
384  TEST_EQUALITY(c.rank(),2);
385  TEST_EQUALITY(c.extent(0),5);
386  TEST_EQUALITY(c.extent(1),3);
387  TEST_EQUALITY(dimension_scalar(b),derivative_dim_plus_one);
388  }
389 
390 }
391 
392 int main( int argc, char* argv[] ) {
393  Teuchos::GlobalMPISession mpiSession(&argc, &argv);
394 
395  Kokkos::initialize();
396 
397  int res = Teuchos::UnitTestRepository::runUnitTestsFromMain(argc, argv);
398 
399  Kokkos::finalize();
400 
401  return res;
402 }
std::enable_if< is_view< InputViewType >::value||is_dyn_rank_view< InputViewType >::value, ResultViewType >::type createDynRankViewWithType(const InputViewType &a, const CtorProp &prop, const Dims... dims)
Wrapper to simplify use of Sacado ViewFactory.
std::enable_if< is_view< InputViewType >::value||is_dyn_rank_view< InputViewType >::value, typename Impl::ResultDynRankView< InputViewType >::type >::type createDynRankView(const InputViewType &a, const CtorProp &prop, const Dims... dims)
Wrapper to simplify use of Sacado ViewFactory.
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
int main(int argc, char *argv[])
std::enable_if< is_view< InputViewType >::value||is_dyn_rank_view< InputViewType >::value, ResultViewType >::type createViewWithType(const InputViewType &a, const CtorProp &prop, const Dims... dims)
Wrapper to simplify use of Sacado ViewFactory.
int value
TEUCHOS_UNIT_TEST(view_factory, dyn_rank_views)