Stokhos Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
Stokhos_Tpetra_UQ_PCE.hpp
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#ifndef STOKHOS_TPETRA_UQ_PCE_HPP
43#define STOKHOS_TPETRA_UQ_PCE_HPP
44
45// This header file should be included whenever compiling any Tpetra
46// code with Stokhos scalar types
47
48// MP includes and specializations
50
51// Kokkos includes
52#include "Tpetra_ConfigDefs.hpp"
53#include "Tpetra_MultiVector_fwd.hpp"
54#include "Tpetra_Vector_fwd.hpp"
55#include "Tpetra_Access.hpp"
56#include "Kokkos_Core.hpp"
57#include "KokkosCompat_ClassicNodeAPI_Wrapper.hpp"
58#include "KokkosCompat_View.hpp"
59#include "KokkosCompat_View_def.hpp"
60
61// Hack for mean-based prec-setup where we get the PCE size from the first
62// entry in the Teuchos::ArrayView. This is clearly quite dangerous and is
63// likely to bite us in the ass at some point!
64namespace Kokkos {
65 namespace Compat {
66 template <typename D, typename S>
67 Kokkos::View<Sacado::UQ::PCE<S>*,D>
68 getKokkosViewDeepCopy(const Teuchos::ArrayView< Sacado::UQ::PCE<S> >& a) {
69 typedef Sacado::UQ::PCE<S> T;
70 typedef typename std::conditional<
71 ::Kokkos::SpaceAccessibility< D, Kokkos::HostSpace>::accessible,
72 typename D::execution_space, Kokkos::HostSpace>::type
73 HostDevice;
74 typedef Kokkos::View<T*,D> view_type;
75 typedef Kokkos::View<T*,typename view_type::array_layout,HostDevice,Kokkos::MemoryUnmanaged> unmanaged_host_view_type;
76 if (a.size() == 0)
77 return view_type();
78 view_type v("", a.size(), a[0].size());
79 unmanaged_host_view_type hv(a.getRawPtr(), a.size(), a[0].size());
81 return v;
82 }
83
84 template <typename D, typename S>
85 Kokkos::View<const Sacado::UQ::PCE<S>*,D>
86 getKokkosViewDeepCopy(const Teuchos::ArrayView<const Sacado::UQ::PCE<S> >& a) {
87 typedef Sacado::UQ::PCE<S> T;
88 typedef typename std::conditional<
89 ::Kokkos::SpaceAccessibility< D, Kokkos::HostSpace>::accessible,
90 typename D::execution_space, Kokkos::HostSpace>::type
91 HostDevice;
92 typedef Kokkos::View<T*,D> view_type;
93 typedef Kokkos::View<const T*,typename view_type::array_layout,HostDevice,Kokkos::MemoryUnmanaged> unmanaged_host_view_type;
94 if (a.size() == 0)
95 return view_type();
96 view_type v("", a.size(), a[0].size());
97 unmanaged_host_view_type hv(a.getRawPtr(), a.size(), a[0].size());
99 return v;
100 }
101 }
102}
103
104// Kokkos-Linalg
107#include "Kokkos_MV_UQ_PCE.hpp"
115
116namespace Stokhos {
117
118// Traits for determining device type from node type
119template <typename Node>
121 // Prefer Serial execution space as the default, but if that's not
122 // available, use the Host memory space's default execution space.
123#if defined(KOKKOS_ENABLE_SERIAL)
124 typedef Kokkos::Serial type;
125#else
126 typedef Kokkos::HostSpace::execution_space type;
127#endif // defined(KOKKOS_ENABLE_SERIAL)
128};
129
130template <typename Device>
131struct DeviceForNode2< Kokkos::Compat::KokkosDeviceWrapperNode<Device> > {
132 typedef Device type;
133};
134
135}
136
137#include "Tpetra_Details_PackTraits.hpp"
138#include "Tpetra_Details_ScalarViewTraits.hpp"
139
140namespace Tpetra {
141namespace Details {
142
146template<class S>
147struct PackTraits<Sacado::UQ::PCE<S>> {
149
152 static const bool compileTimeSize = false;
153
154 using input_buffer_type = Kokkos::View<const char*, Kokkos::AnonymousSpace>;
155 using output_buffer_type = Kokkos::View<char*, Kokkos::AnonymousSpace>;
156 using input_array_type = Kokkos::View<const value_type*, Kokkos::AnonymousSpace>;
157 using output_array_type = Kokkos::View<value_type*, Kokkos::AnonymousSpace>;
158
159 using scalar_value_type = typename value_type::value_type;
160 using SPT = PackTraits<scalar_value_type>;
161 using scalar_input_array_type = typename SPT::input_array_type;
162 using scalar_output_array_type = typename SPT::output_array_type;
163
164 KOKKOS_INLINE_FUNCTION
165 static size_t numValuesPerScalar (const value_type& x) {
166 return x.size ();
167 }
168
169 KOKKOS_INLINE_FUNCTION
170 static Kokkos::pair<int, size_t>
171 packArray (char outBuf[],
172 const value_type inBuf[],
173 const size_t numEnt)
174 {
175 using return_type = Kokkos::pair<int, size_t>;
176 size_t numBytes = 0;
177 int errorCode = 0;
178
179 if (numEnt == 0) {
180 return return_type (errorCode, numBytes);
181 }
182 else {
183 // Check whether input array is contiguously allocated based on the size
184 // of the first entry. We can only pack contiguously allocated data
185 // since that is the only way we can guarrantee all of the PCE arrays
186 // are the same size and the buffer will allocated correctly.
187 const size_t scalar_size = numValuesPerScalar (inBuf[0]);
188 const scalar_value_type* last_coeff = inBuf[numEnt - 1].coeff ();
189 const scalar_value_type* last_coeff_expected =
190 inBuf[0].coeff () + (numEnt - 1)*scalar_size;
191 const bool is_contiguous = (last_coeff == last_coeff_expected);
192 if (! is_contiguous) {
193 // Cannot pack non-contiguous PCE array since buffer size calculation
194 // is likely wrong.
195 errorCode = 3;
196 return return_type (errorCode, numBytes);
197 }
198
199 // Check we are packing length-1 PCE arrays (mean-based preconditioner).
200 // We can technically pack length > 1, but the unpack assumes the
201 // output array is sized appropriately. Currently this is not the case
202 // in Tpetra::CrsMatrix::transferAndFillComplete() which allocates a
203 // local Teuchos::Array for the CSR values, which will only be length-1
204 // by default.
205 if (scalar_size != 1) {
206 // Cannot pack PCE array with pce_size > 1 since unpack array
207 // may not be allocated correctly.
208 errorCode = 4;
209 return return_type (errorCode, numBytes);
210 }
211
212 const size_t flat_numEnt = numEnt * scalar_size;
213 return SPT::packArray (outBuf, inBuf[0].coeff (), flat_numEnt);
214 }
215 }
216
217 KOKKOS_INLINE_FUNCTION
218 static Kokkos::pair<int, size_t>
220 const char inBuf[],
221 const size_t numEnt)
222 {
223 using return_type = Kokkos::pair<int, size_t>;
224 size_t numBytes = 0;
225 int errorCode = 0;
226
227 if (numEnt == 0) {
228 return return_type (errorCode, numBytes);
229 }
230 else {
231 // Check whether output array is contiguously allocated based on the size
232 // of the first entry. We have a simpler method to unpack in this case
233 const size_t scalar_size = numValuesPerScalar (outBuf[0]);
234 const scalar_value_type* last_coeff = outBuf[numEnt - 1].coeff ();
235 const scalar_value_type* last_coeff_expected =
236 outBuf[0].coeff () + (numEnt - 1) * scalar_size;
237 const bool is_contiguous = (last_coeff == last_coeff_expected);
238
239 if (is_contiguous) {
240 // Unpack all of the PCE coefficients for the whole array
241 const size_t flat_numEnt = numEnt * scalar_size;
242 return SPT::unpackArray (outBuf[0].coeff (), inBuf, flat_numEnt);
243 }
244 else {
245 // Unpack one entry at a time. This assumes each entry of outBuf
246 // is the correct size based on the packing. This is is only
247 // guarranteed to be true for pce_size == 1, hence the check in
248 // packArray().
249 size_t numBytesTotal = 0;
250 for (size_t i = 0; i < numEnt; ++i) {
251 const char* inBufVal = inBuf + numBytesTotal;
252 const size_t numBytes = unpackValue (outBuf[i], inBufVal);
253 numBytesTotal += numBytes;
254 }
255 return return_type (errorCode, numBytesTotal);
256 }
257 }
258 }
259
260 KOKKOS_INLINE_FUNCTION
261 static size_t
263 {
264 return inVal.size () * SPT::packValueCount (inVal.val ());
265 }
266
267 KOKKOS_INLINE_FUNCTION
268 static size_t
269 packValue (char outBuf[],
270 const value_type& inVal)
271 {
272 const size_t numBytes = packValueCount (inVal);
273 memcpy (outBuf, inVal.coeff (), numBytes);
274 return numBytes;
275 }
276
277 KOKKOS_INLINE_FUNCTION
278 static size_t
279 packValue (char outBuf[],
280 const size_t outBufIndex,
281 const value_type& inVal)
282 {
283 const size_t numBytes = packValueCount (inVal);
284 const size_t offset = outBufIndex * numBytes;
285 memcpy (outBuf + offset, inVal.coeff (), numBytes);
286 return numBytes;
287 }
288
289 KOKKOS_INLINE_FUNCTION
290 static size_t
291 unpackValue (value_type& outVal, const char inBuf[])
292 {
293 const size_t numBytes = packValueCount (outVal);
294 memcpy (outVal.coeff (), inBuf, numBytes);
295 return numBytes;
296 }
297}; // struct PackTraits
298
304template<typename S, typename D>
305struct ScalarViewTraits<Sacado::UQ::PCE<S>, D> {
307 using device_type = D;
308
309 static Kokkos::View<value_type*, device_type>
311 const size_t numEnt,
312 const std::string& label = "")
313 {
314 const size_t numVals = PackTraits<value_type>::numValuesPerScalar (x);
315 using view_type = Kokkos::View<value_type*, device_type>;
316 return view_type (label, numEnt, numVals);
317 }
318};
319
320} // namespace Details
321} // namespace Tpetra
322
323namespace Kokkos {
324 template <class S, class L, class G, class N>
325 size_t dimension_scalar(const Tpetra::MultiVector<S,L,G,N>& mv) {
326 if ( mv.need_sync_device() ) {
327 // NOTE (mfh 02 Apr 2019) This doesn't look right. Shouldn't I
328 // want the most recently updated View, which in this case would
329 // be the host View? However, this is what I found when I
330 // changed these lines not to call deprecated code, so I'm
331 // leaving it.
332 return dimension_scalar(mv.getLocalViewDevice(Tpetra::Access::ReadOnly));
333 }
334 return dimension_scalar(mv.getLocalViewHost(Tpetra::Access::ReadOnly));
335 }
336
337 template <class S, class L, class G, class N>
338 size_t dimension_scalar(const Tpetra::Vector<S,L,G,N>& v) {
339 if ( v.need_sync_device() ) {
340 // NOTE (mfh 02 Apr 2019) This doesn't look right. Shouldn't I
341 // want the most recently updated View, which in this case would
342 // be the host View? However, this is what I found when I
343 // changed these lines not to call deprecated code, so I'm
344 // leaving it.
345 return dimension_scalar(v.getLocalViewDevice(Tpetra::Access::ReadOnly));
346 }
347 return dimension_scalar(v.getLocalViewHost(Tpetra::Access::ReadOnly));
348 }
349}
350
351#endif // STOKHOS_TPETRA_UQ_PCE_HPP
Kokkos::View< Sacado::UQ::PCE< S > *, D > getKokkosViewDeepCopy(const Teuchos::ArrayView< Sacado::UQ::PCE< S > > &a)
KOKKOS_INLINE_FUNCTION constexpr std::enable_if< is_view_uq_pce< View< T, P... > >::value, unsigned >::type dimension_scalar(const View< T, P... > &view)
void deep_copy(const Stokhos::CrsMatrix< ValueType, DstDevice, Layout > &dst, const Stokhos::CrsMatrix< ValueType, SrcDevice, Layout > &src)
Top-level namespace for Stokhos classes and functions.
Kokkos::HostSpace::execution_space type
Kokkos::View< char *, Kokkos::AnonymousSpace > output_buffer_type
Kokkos::View< const value_type *, Kokkos::AnonymousSpace > input_array_type
static KOKKOS_INLINE_FUNCTION size_t packValue(char outBuf[], const size_t outBufIndex, const value_type &inVal)
static KOKKOS_INLINE_FUNCTION size_t packValueCount(const value_type &inVal)
static KOKKOS_INLINE_FUNCTION Kokkos::pair< int, size_t > packArray(char outBuf[], const value_type inBuf[], const size_t numEnt)
static KOKKOS_INLINE_FUNCTION size_t packValue(char outBuf[], const value_type &inVal)
static KOKKOS_INLINE_FUNCTION Kokkos::pair< int, size_t > unpackArray(value_type outBuf[], const char inBuf[], const size_t numEnt)
static KOKKOS_INLINE_FUNCTION size_t numValuesPerScalar(const value_type &x)
Kokkos::View< value_type *, Kokkos::AnonymousSpace > output_array_type
Kokkos::View< const char *, Kokkos::AnonymousSpace > input_buffer_type
static KOKKOS_INLINE_FUNCTION size_t unpackValue(value_type &outVal, const char inBuf[])
static Kokkos::View< value_type *, device_type > allocateArray(const value_type &x, const size_t numEnt, const std::string &label="")