Sacado Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
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__)
47#include "Teuchos_Assert.hpp"
48
49// Example to demonstrate the use of Kokkos::ViewFactory for constructing
50// view's of Fad's without explicitly referencing the sacado dimension
51
52// An example function that takes two views of various ranks and scalar types
53// and produces a third allocated using the ViewFactory
54template <class View1, class View2>
55Kokkos::View< typename Kokkos::ViewFactory<View1,View2>::value_type*,
56 typename View1::device_type >
57my_func(const View1& v1, const View2& v2)
58{
60 typedef typename ViewFac::value_type value_type;
61 typedef typename View1::device_type device_type;
62 typedef typename View1::size_type size_type;
63 typedef Kokkos::View<value_type*,device_type> ViewTmp;
64
65 static_assert( unsigned(View1::rank) == 2, "" );
66 static_assert( unsigned(View2::rank) == 1, "" );
67
68 const size_type m = v1.extent(0);
69 const size_type n = v1.extent(1);
70 ViewTmp vtmp = ViewFac::template create_view<ViewTmp>(v1,v2,"tmp",m);
71
72 Kokkos::parallel_for(m, KOKKOS_LAMBDA (const size_type i) {
73 value_type t = 0;
74 for (size_type j=0; j<n; ++j)
75 t += v1(i,j)*v2(j);
76 vtmp(i) = t;
77 });
78
79 return vtmp;
80}
81
82#if defined(HAVE_SACADO_KOKKOSCONTAINERS)
83// An example function that takes two dynamic-rank views of various ranks and
84// scalar types and produces a third allocated using the ViewFactory
85template <class View1, class View2>
86Kokkos::DynRankView< typename Kokkos::ViewFactory<View1,View2>::value_type,
87 typename View1::device_type >
88my_func_dynamic(const View1& v1, const View2& v2)
89{
91 typedef typename ViewFac::value_type value_type;
92 typedef typename View1::device_type device_type;
93 typedef typename View1::size_type size_type;
94 typedef Kokkos::DynRankView<value_type,device_type> ViewTmp;
95
96 TEUCHOS_ASSERT( v1.rank() == 2 );
97 TEUCHOS_ASSERT( v2.rank() == 1 );
98
99 const size_type m = v1.extent(0);
100 const size_type n = v1.extent(1);
101 ViewTmp vtmp = ViewFac::template create_view<ViewTmp>(v1,v2,"tmp",m);
102
103 Kokkos::parallel_for(m, KOKKOS_LAMBDA (const size_type i) {
104 value_type t = 0;
105 for (size_type j=0; j<n; ++j)
106 t += v1(i,j)*v2(j);
107 vtmp(i) = t;
108 });
109
110 return vtmp;
111}
112#endif
113
114#endif
115
116int main(int argc, char* argv[]) {
117
118#if !defined(__CUDACC__)
120 typedef Kokkos::DefaultExecutionSpace execution_space;
121
122 Kokkos::initialize(argc, argv);
123
124 const size_t m = 10;
125 const size_t n = 2;
126 const size_t deriv_dim = 1;
127
128 // First use the static-rank view
129 {
130
131 // Calculation with double's
132 Kokkos::View<double**,execution_space> v1("v1",m,n);
133 Kokkos::View<double* ,execution_space> v2("v2",n);
134
135 Kokkos::deep_copy(v1, 2.0);
136 Kokkos::deep_copy(v2, 3.0);
137
138 auto v3 = my_func(v1,v2);
139
140 std::cout << "v3 = " << std::endl;
141 for (size_t i=0; i<m; ++i) {
142 std::cout << "\t" << v3(i) << std::endl;
143 }
144
145 // Calculation with Fad's (initialize each component of v2_fad with a
146 // Fad object with value == 3.0, one derivative == 1
147 Kokkos::View<FadType*,execution_space> v2_fad("v2_fad",n,deriv_dim+1);
148 Kokkos::deep_copy(v2_fad, FadType(deriv_dim, 0, 3.0));
149
150 auto v3_fad = my_func(v1,v2_fad);
151
152 std::cout << "v3_fad = " << std::endl;
153 for (size_t i=0; i<m; ++i) {
154 std::cout << "\t" << v3_fad(i) << std::endl;
155 }
156
157 }
158
159#if defined(HAVE_SACADO_KOKKOSCONTAINERS)
160 // Now use the dynamic-rank view
161 {
162
163 // Calculation with double's
164 Kokkos::DynRankView<double,execution_space> v1("v1",m,n);
165 Kokkos::DynRankView<double,execution_space> v2("v2",n);
166
167 Kokkos::deep_copy(v1, 2.0);
168 Kokkos::deep_copy(v2, 3.0);
169
170 auto v3 = my_func_dynamic(v1,v2);
171
172 std::cout << "v3 = " << std::endl;
173 for (size_t i=0; i<m; ++i) {
174 std::cout << "\t" << v3(i) << std::endl;
175 }
176
177 // Calculation with Fad's (initialize each component of v2_fad with a
178 // Fad object with value == 3.0, one derivative == 1
179 Kokkos::DynRankView<FadType,execution_space> v2_fad("v2_fad",n,deriv_dim+1);
180 Kokkos::deep_copy(v2_fad, FadType(deriv_dim, 0, 3.0));
181
182 auto v3_fad = my_func_dynamic(v1,v2_fad);
183
184 std::cout << "v3_fad = " << std::endl;
185 for (size_t i=0; i<m; ++i) {
186 std::cout << "\t" << v3_fad(i) << std::endl;
187 }
188
189 }
190#endif
191
192 Kokkos::finalize();
193#endif
194
195 return 0;
196}
int main()
Sacado::Fad::DFad< double > FadType
Fad specializations for Teuchos::BLAS wrappers.
Kokkos::View< typename Kokkos::ViewFactory< View1, View2 >::value_type *, typename View1::device_type > my_func(const View1 &v1, const View2 &v2)