Sacado Package Browser (Single Doxygen Collection)  Version of the Day
view_factory_example.cpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Stokhos Package
5 // Copyright (2009) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Eric T. Phipps (etphipp@sandia.gov).
38 //
39 // ***********************************************************************
40 // @HEADER
41 
42 #include "Sacado.hpp"
43 
44 #if !defined(__CUDACC__)
46 #include "Kokkos_ViewFactory.hpp"
47 
48 // Example to demonstrate the use of Kokkos::ViewFactory for constructing
49 // view's of Fad's without explicitly referencing the sacado dimension
50 
51 // An example function that takes two views of various ranks and scalar types
52 // and produces a third allocated using the ViewFactory
53 template <class View1, class View2>
54 Kokkos::View< typename Kokkos::ViewFactory<View1,View2>::value_type*,
55  typename View1::device_type >
56 my_func(const View1& v1, const View2& v2)
57 {
58  typedef Kokkos::ViewFactory<View1,View2> ViewFac;
59  typedef typename ViewFac::value_type value_type;
60  typedef typename View1::device_type device_type;
61  typedef typename View1::size_type size_type;
62  typedef Kokkos::View<value_type*,device_type> ViewTmp;
63 
64  static_assert( unsigned(View1::rank) == 2, "" );
65  static_assert( unsigned(View2::rank) == 1, "" );
66 
67  const size_type m = v1.extent(0);
68  const size_type n = v1.extent(1);
69  ViewTmp vtmp = ViewFac::template create_view<ViewTmp>(v1,v2,"tmp",m);
70 
71  Kokkos::parallel_for(m, KOKKOS_LAMBDA (const size_type i) {
72  value_type t = 0;
73  for (size_type j=0; j<n; ++j)
74  t += v1(i,j)*v2(j);
75  vtmp(i) = t;
76  });
77 
78  return vtmp;
79 }
80 
81 #if defined(HAVE_SACADO_KOKKOSCONTAINERS)
82 // An example function that takes two dynamic-rank views of various ranks and
83 // scalar types and produces a third allocated using the ViewFactory
84 template <class View1, class View2>
85 Kokkos::DynRankView< typename Kokkos::ViewFactory<View1,View2>::value_type,
86  typename View1::device_type >
87 my_func_dynamic(const View1& v1, const View2& v2)
88 {
89  typedef Kokkos::ViewFactory<View1,View2> ViewFac;
90  typedef typename ViewFac::value_type value_type;
91  typedef typename View1::device_type device_type;
92  typedef typename View1::size_type size_type;
93  typedef Kokkos::DynRankView<value_type,device_type> ViewTmp;
94 
95  TEUCHOS_ASSERT( v1.rank() == 2 );
96  TEUCHOS_ASSERT( v2.rank() == 1 );
97 
98  const size_type m = v1.extent(0);
99  const size_type n = v1.extent(1);
100  ViewTmp vtmp = ViewFac::template create_view<ViewTmp>(v1,v2,"tmp",m);
101 
102  Kokkos::parallel_for(m, KOKKOS_LAMBDA (const size_type i) {
103  value_type t = 0;
104  for (size_type j=0; j<n; ++j)
105  t += v1(i,j)*v2(j);
106  vtmp(i) = t;
107  });
108 
109  return vtmp;
110 }
111 #endif
112 
113 #endif
114 
115 int main(int argc, char* argv[]) {
116 
117 #if !defined(__CUDACC__)
119  typedef Kokkos::DefaultExecutionSpace execution_space;
120 
121  Kokkos::initialize(argc, argv);
122 
123  const size_t m = 10;
124  const size_t n = 2;
125  const size_t deriv_dim = 1;
126 
127  // First use the static-rank view
128  {
129 
130  // Calculation with double's
131  Kokkos::View<double**,execution_space> v1("v1",m,n);
132  Kokkos::View<double* ,execution_space> v2("v2",n);
133 
134  Kokkos::deep_copy(v1, 2.0);
135  Kokkos::deep_copy(v2, 3.0);
136 
137  auto v3 = my_func(v1,v2);
138 
139  std::cout << "v3 = " << std::endl;
140  for (size_t i=0; i<m; ++i) {
141  std::cout << "\t" << v3(i) << std::endl;
142  }
143 
144  // Calculation with Fad's (initialize each component of v2_fad with a
145  // Fad object with value == 3.0, one derivative == 1
146  Kokkos::View<FadType*,execution_space> v2_fad("v2_fad",n,deriv_dim+1);
147  Kokkos::deep_copy(v2_fad, FadType(deriv_dim, 0, 3.0));
148 
149  auto v3_fad = my_func(v1,v2_fad);
150 
151  std::cout << "v3_fad = " << std::endl;
152  for (size_t i=0; i<m; ++i) {
153  std::cout << "\t" << v3_fad(i) << std::endl;
154  }
155 
156  }
157 
158 #if defined(HAVE_SACADO_KOKKOSCONTAINERS)
159  // Now use the dynamic-rank view
160  {
161 
162  // Calculation with double's
163  Kokkos::DynRankView<double,execution_space> v1("v1",m,n);
164  Kokkos::DynRankView<double,execution_space> v2("v2",n);
165 
166  Kokkos::deep_copy(v1, 2.0);
167  Kokkos::deep_copy(v2, 3.0);
168 
169  auto v3 = my_func_dynamic(v1,v2);
170 
171  std::cout << "v3 = " << std::endl;
172  for (size_t i=0; i<m; ++i) {
173  std::cout << "\t" << v3(i) << std::endl;
174  }
175 
176  // Calculation with Fad's (initialize each component of v2_fad with a
177  // Fad object with value == 3.0, one derivative == 1
178  Kokkos::DynRankView<FadType,execution_space> v2_fad("v2_fad",n,deriv_dim+1);
179  Kokkos::deep_copy(v2_fad, FadType(deriv_dim, 0, 3.0));
180 
181  auto v3_fad = my_func_dynamic(v1,v2_fad);
182 
183  std::cout << "v3_fad = " << std::endl;
184  for (size_t i=0; i<m; ++i) {
185  std::cout << "\t" << v3_fad(i) << std::endl;
186  }
187 
188  }
189 #endif
190 
191  Kokkos::finalize();
192 #endif
193 
194  return 0;
195 }
Kokkos::View< typename Kokkos::ViewFactory< View1, View2 >::value_type *, typename View1::device_type > my_func(const View1 &v1, const View2 &v2)
Sacado::Fad::DFad< double > FadType
int main(int argc, char *argv[])
#define TEUCHOS_ASSERT(assertion_test)
int n