Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Kokkos_MV_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 KOKKOS_MV_UQ_PCE_HPP
43 #define KOKKOS_MV_UQ_PCE_HPP
44 
45 #include "Sacado_UQ_PCE.hpp"
46 #include "Kokkos_View_UQ_PCE.hpp"
47 #include "Kokkos_InnerProductSpaceTraits_UQ_PCE.hpp"
48 #include "Kokkos_Blas1_UQ_PCE.hpp"
49 /*
50 //----------------------------------------------------------------------------
51 // Specializations of Kokkos Vector/MultiVector math functions
52 //----------------------------------------------------------------------------
53 
54 namespace Kokkos {
55 
56 // Rank-1 vector add with Sacado::UQ::PCE scalar type, constant a, b
57 template <typename RS, typename RL, typename RD, typename RM,
58  typename XS, typename XL, typename XD, typename XM,
59  typename YS, typename YL, typename YD, typename YM>
60 Kokkos::View< Sacado::UQ::PCE<RS>*, RL, RD, RM>
61 V_Add( const Kokkos::View< Sacado::UQ::PCE<RS>*, RL, RD, RM >& r,
62  const typename Sacado::UQ::PCE<XS>::value_type& av,
63  const Kokkos::View< Sacado::UQ::PCE<XS>*, XL, XD, XM >& x,
64  const typename Sacado::UQ::PCE<XS>::value_type& bv,
65  const Kokkos::View< Sacado::UQ::PCE<YS>*, YL, YD, YM >& y,
66  int n = -1)
67 {
68  typedef Kokkos::View< Sacado::UQ::PCE<RS>*, RL, RD, RM > RVector;
69  typedef Kokkos::View< Sacado::UQ::PCE<XS>*, XL, XD, XM > XVector;
70  typedef Kokkos::View< Sacado::UQ::PCE<YS>*, YL, YD, YM > YVector;
71 
72  typename RVector::flat_array_type r_flat = r;
73  typename XVector::flat_array_type x_flat = x;
74  typename YVector::flat_array_type y_flat = y;
75  if (n != -1) n = n * r.sacado_size();
76 
77  V_Add( r_flat, av, x_flat, bv, y_flat, n );
78 
79  return r;
80 }
81 
82 // Rank-1 vector add with Sacado::UQ::PCE scalar type, non-constant a, b
83 template <typename RS, typename RL, typename RD, typename RM,
84  typename XS, typename XL, typename XD, typename XM,
85  typename YS, typename YL, typename YD, typename YM>
86 Kokkos::View< Sacado::UQ::PCE<RS>*, RL, RD, RM>
87 V_Add( const Kokkos::View< Sacado::UQ::PCE<RS>*, RL, RD, RM >& r,
88  const Sacado::UQ::PCE<XS>& av,
89  const Kokkos::View< Sacado::UQ::PCE<XS>*, XL, XD, XM >& x,
90  const Sacado::UQ::PCE<XS>& bv,
91  const Kokkos::View< Sacado::UQ::PCE<YS>*, YL, YD, YM >& y,
92  int n = -1)
93 {
94  if (Sacado::is_constant(av) && Sacado::is_constant(bv)) {
95  return V_Add( r, av.fastAccessCoeff(0), x, bv.fastAccessCoeff(0), y, n );
96  }
97  else {
98  Impl::raise_error("V_Add not implemented for non-constant a or b");
99  }
100  return r;
101 }
102 
103 // Rank-2 vector add with Sacado::UQ::PCE scalar type, constant a, b
104 template <typename RS, typename RL, typename RD, typename RM,
105  typename XS, typename XL, typename XD, typename XM,
106  typename YS, typename YL, typename YD, typename YM>
107 Kokkos::View< Sacado::UQ::PCE<RS>**, RL, RD, RM>
108 MV_Add( const Kokkos::View< Sacado::UQ::PCE<RS>**, RL, RD, RM >& r,
109  const typename Sacado::UQ::PCE<XS>::value_type& av,
110  const Kokkos::View< Sacado::UQ::PCE<XS>**, XL, XD, XM >& x,
111  const typename Sacado::UQ::PCE<XS>::value_type& bv,
112  const Kokkos::View< Sacado::UQ::PCE<YS>**, YL, YD, YM >& y,
113  int n = -1)
114 {
115  typedef Kokkos::View< Sacado::UQ::PCE<RS>**, RL, RD, RM > RVector;
116  typedef Kokkos::View< Sacado::UQ::PCE<XS>**, XL, XD, XM > XVector;
117  typedef Kokkos::View< Sacado::UQ::PCE<YS>**, YL, YD, YM > YVector;
118 
119  typename RVector::flat_array_type r_flat = r;
120  typename XVector::flat_array_type x_flat = x;
121  typename YVector::flat_array_type y_flat = y;
122  if (n != -1) n = n * r.sacado_size();
123 
124  MV_Add( r_flat, av, x_flat, bv, y_flat, n );
125 
126  return r;
127 }
128 
129 // Rank-2 vector add with Sacado::UQ::PCE scalar type, non-constant a, b
130 template <typename RS, typename RL, typename RD, typename RM,
131  typename XS, typename XL, typename XD, typename XM,
132  typename YS, typename YL, typename YD, typename YM>
133 Kokkos::View< Sacado::UQ::PCE<RS>**, RL, RD, RM>
134 MV_Add( const Kokkos::View< Sacado::UQ::PCE<RS>**, RL, RD, RM >& r,
135  const Sacado::UQ::PCE<XS>& av,
136  const Kokkos::View< Sacado::UQ::PCE<XS>**, XL, XD, XM >& x,
137  const Sacado::UQ::PCE<XS>& bv,
138  const Kokkos::View< Sacado::UQ::PCE<YS>**, YL, YD, YM >& y,
139  int n = -1)
140 {
141  if (Sacado::is_constant(av) && Sacado::is_constant(bv)) {
142  return MV_Add( r, av.fastAccessCoeff(0), x, bv.fastAccessCoeff(0), y, n );
143  }
144  else {
145  Impl::raise_error("MV_Add not implemented for non-constant a or b");
146  }
147  return r;
148 }
149 
150 // Rank-1 dot product
151 template <typename XS, typename XL, typename XD, typename XM,
152  typename YS, typename YL, typename YD, typename YM>
153 typename Details::InnerProductSpaceTraits< Sacado::UQ::PCE<XS> >::dot_type
154 V_Dot( const Kokkos::View< Sacado::UQ::PCE<XS>*, XL, XD, XM >& x,
155  const Kokkos::View< Sacado::UQ::PCE<YS>*, YL, YD, YM >& y,
156  int n = -1 )
157 {
158  typedef Kokkos::View< Sacado::UQ::PCE<XS>*, XL, XD, XM > XVector;
159  typedef Kokkos::View< Sacado::UQ::PCE<YS>*, YL, YD, YM > YVector;
160 
161  typename XVector::flat_array_type x_flat = x;
162  typename YVector::flat_array_type y_flat = y;
163  if (n != -1) n = n * x.sacado_size();
164 
165  return V_Dot( x_flat, y_flat, n );
166 }
167 
168 // Rank-2 dot product
169 template <typename rVector,
170  typename XS, typename XL, typename XD, typename XM,
171  typename YS, typename YL, typename YD, typename YM>
172 void
173 MV_Dot( const rVector& r,
174  const Kokkos::View< Sacado::UQ::PCE<XS>**, XL, XD, XM >& x,
175  const Kokkos::View< Sacado::UQ::PCE<YS>**, YL, YD, YM >& y,
176  int n = -1 )
177 {
178  typedef Kokkos::View< Sacado::UQ::PCE<XS>**, XL, XD, XM > XVector;
179  typedef Kokkos::View< Sacado::UQ::PCE<YS>**, YL, YD, YM > YVector;
180 
181  typename XVector::flat_array_type x_flat = x;
182  typename YVector::flat_array_type y_flat = y;
183  if (n != -1) n = n * x.sacado_size();
184 
185  MV_Dot( r, x_flat, y_flat, n );
186 }
187 
188 template<class VT1, class VT2, class VT3>
189 struct MV_ElementWiseMultiplyFunctor;
190 
191 template <typename CT, typename CD, typename CM,
192  typename AT, typename AD, typename AM,
193  typename BT, typename BD, typename BM>
194 struct MV_ElementWiseMultiplyFunctor<
195  View< CT,LayoutLeft,CD,CM,Impl::ViewPCEContiguous >,
196  View< AT,LayoutLeft,AD,AM,Impl::ViewPCEContiguous >,
197  View< BT,LayoutLeft,BD,BM,Impl::ViewPCEContiguous > >
198 {
199  typedef View< CT,LayoutLeft,CD,CM,Impl::ViewPCEContiguous > CVector;
200  typedef View< AT,LayoutLeft,AD,AM,Impl::ViewPCEContiguous > AVector;
201  typedef View< BT,LayoutLeft,BD,BM,Impl::ViewPCEContiguous > BVector;
202 
203  typedef typename CVector::array_type CArray;
204  typedef typename AVector::array_type AArray;
205  typedef typename BVector::array_type BArray;
206 
207  typedef typename CArray::execution_space execution_space;
208  typedef typename CArray::size_type size_type;
209 
210  typename CArray::const_value_type m_c;
211  CArray m_C;
212  typename AArray::const_value_type m_ab;
213  typename AArray::const_type m_A ;
214  typename BArray::const_type m_B ;
215  const size_type m_n;
216  const size_type m_n_pce;
217 
218  MV_ElementWiseMultiplyFunctor(typename CVector::const_value_type c,
219  CVector C,
220  typename AVector::const_value_type ab,
221  typename AVector::const_type A,
222  typename BVector::const_type B,
223  const size_type n) :
224  m_c(c.fastAccessCoeff(0)),
225  m_C(C),
226  m_ab(ab.fastAccessCoeff(0)),
227  m_A(A),
228  m_B(B),
229  m_n(n),
230  m_n_pce(C.sacado_size()) {}
231 
232  KOKKOS_INLINE_FUNCTION
233  void operator()( const size_type i) const
234  {
235  // Currently specialized for use case where A is degree-0
236  typename AArray::const_value_type Ai = m_A(0,i);
237  for (size_type k=0; k<m_n; ++k) {
238 #ifdef KOKKOS_ENABLE_PRAGMA_IVDEP
239 #pragma ivdep
240 #endif
241  for (size_type l=0; l<m_n_pce; ++l)
242  m_C(l,i,k) = m_c*m_C(l,i,k) + m_ab*Ai*m_B(l,i,k);
243  }
244  }
245 };
246 
247 template<class VT1, class VT2, class VT3>
248 struct V_ElementWiseMultiplyFunctor;
249 
250 template <typename CT, typename CD, typename CM,
251  typename AT, typename AD, typename AM,
252  typename BT, typename BD, typename BM>
253 struct V_ElementWiseMultiplyFunctor<
254  View< CT,LayoutLeft,CD,CM,Impl::ViewPCEContiguous >,
255  View< AT,LayoutLeft,AD,AM,Impl::ViewPCEContiguous >,
256  View< BT,LayoutLeft,BD,BM,Impl::ViewPCEContiguous > >
257 {
258  typedef View< CT,LayoutLeft,CD,CM,Impl::ViewPCEContiguous > CVector;
259  typedef View< AT,LayoutLeft,AD,AM,Impl::ViewPCEContiguous > AVector;
260  typedef View< BT,LayoutLeft,BD,BM,Impl::ViewPCEContiguous > BVector;
261 
262  typedef typename CVector::array_type CArray;
263  typedef typename AVector::array_type AArray;
264  typedef typename BVector::array_type BArray;
265 
266  typedef typename CArray::execution_space execution_space;
267  typedef typename CArray::size_type size_type;
268 
269  typename CArray::const_value_type m_c;
270  CArray m_C;
271  typename AArray::const_value_type m_ab;
272  typename AArray::const_type m_A ;
273  typename BArray::const_type m_B ;
274  const size_type m_n_pce;
275 
276  V_ElementWiseMultiplyFunctor(typename CVector::const_value_type c,
277  CVector C,
278  typename AVector::const_value_type ab,
279  typename AVector::const_type A,
280  typename BVector::const_type B) :
281  m_c(c.fastAccessCoeff(0)),
282  m_C(C),
283  m_ab(ab.fastAccessCoeff(0)),
284  m_A(A),
285  m_B(B),
286  m_n_pce(C.sacado_size()) {}
287 
288  KOKKOS_INLINE_FUNCTION
289  void operator()( const size_type i) const
290  {
291  // Currently specialized for use case where A is degree-0
292  typename AArray::const_value_type Ai = m_A(0,i);
293 #ifdef KOKKOS_ENABLE_PRAGMA_IVDEP
294 #pragma ivdep
295 #endif
296  for (size_type l=0; l<m_n_pce; ++l)
297  m_C(l,i) = m_c*m_C(l,i) + m_ab*Ai*m_B(l,i);
298  }
299 };
300 
301 // Rank-1 vector multiply with Sacado::UQ::PCE scalar type, constant c, ab
302 template <typename CS, typename CL, typename CD, typename CM,
303  typename AS, typename AL, typename AD, typename AM,
304  typename BS, typename BL, typename BD, typename BM>
305 Kokkos::View< Sacado::UQ::PCE<CS>*, CL, CD, CM>
306 V_ElementWiseMultiply(
307  const typename Sacado::UQ::PCE<CS>::value_type& c,
308  const Kokkos::View< Sacado::UQ::PCE<CS>*, CL, CD, CM >& C,
309  const typename Sacado::UQ::PCE<AS>::value_type& ab,
310  const Kokkos::View< Sacado::UQ::PCE<AS>*, AL, AD, AM >& A,
311  const Kokkos::View< Sacado::UQ::PCE<BS>*, BL, BD, BM >& B )
312 {
313  typedef View< Sacado::UQ::PCE<CS>*, CL, CD, CM > CVector;
314  typedef View< Sacado::UQ::PCE<AS>*, AL, AD, AM > AVector;
315  typedef View< Sacado::UQ::PCE<BS>*, BL, BD, BM > BVector;
316 
317  typedef View< typename CVector::data_type, typename CVector::array_layout, typename CVector::execution_space, typename CVector::memory_traits > CView;
318  typedef View< typename AVector::data_type, typename AVector::array_layout, typename AVector::execution_space, typename AVector::memory_traits > AView;
319  typedef View< typename BVector::data_type, typename BVector::array_layout, typename BVector::execution_space, typename BVector::memory_traits > BView;
320 
321  V_ElementWiseMultiplyFunctor<CView,AView,BView> op(c,C,ab,A,B) ;
322  Kokkos::parallel_for( C.extent(0) , op );
323  return C;
324 }
325 
326 // Rank-1 vector multiply with Sacado::UQ::PCE scalar type, non-constant c, ab
327 template <typename CS, typename CL, typename CD, typename CM,
328  typename AS, typename AL, typename AD, typename AM,
329  typename BS, typename BL, typename BD, typename BM>
330 Kokkos::View< Sacado::UQ::PCE<CS>*, CL, CD, CM>
331 V_ElementWiseMultiply(
332  const Sacado::UQ::PCE<CS>& c,
333  const Kokkos::View< Sacado::UQ::PCE<CS>*, CL, CD, CM >& C,
334  const Sacado::UQ::PCE<AS>& ab,
335  const Kokkos::View< Sacado::UQ::PCE<AS>*, AL, AD, AM >& A,
336  const Kokkos::View< Sacado::UQ::PCE<BS>*, BL, BD, BM >& B )
337 {
338  if (Sacado::is_constant(c) && Sacado::is_constant(ab)) {
339  return V_ElementWiseMultiply( c.fastAccessCoeff(0), C,
340  ab.fastAccessCoeff(0), A, B );
341  }
342  else {
343  Impl::raise_error(
344  "V_ElementWiseMultiply not implemented for non-constant c or ab");
345  }
346  return C;
347 }
348 
349 // Rank-2 vector multiply with Sacado::UQ::PCE scalar type, constant c, ab
350 template <typename CS, typename CL, typename CD, typename CM,
351  typename AS, typename AL, typename AD, typename AM,
352  typename BS, typename BL, typename BD, typename BM>
353 Kokkos::View< Sacado::UQ::PCE<CS>**, CL, CD, CM>
354 MV_ElementWiseMultiply(
355  const typename Sacado::UQ::PCE<CS>::value_type& c,
356  const Kokkos::View< Sacado::UQ::PCE<CS>**, CL, CD, CM >& C,
357  const typename Sacado::UQ::PCE<AS>::value_type& ab,
358  const Kokkos::View< Sacado::UQ::PCE<AS>*, AL, AD, AM >& A,
359  const Kokkos::View< Sacado::UQ::PCE<BS>**, BL, BD, BM >& B )
360 {
361  typedef View< Sacado::UQ::PCE<CS>**, CL, CD, CM > CVector;
362  typedef View< Sacado::UQ::PCE<AS>*, AL, AD, AM > AVector;
363  typedef View< Sacado::UQ::PCE<BS>**, BL, BD, BM > BVector;
364 
365  typedef View< typename CVector::data_type, typename CVector::array_layout, typename CVector::execution_space, typename CVector::memory_traits > CView;
366  typedef View< typename AVector::data_type, typename AVector::array_layout, typename AVector::execution_space, typename AVector::memory_traits > AView;
367  typedef View< typename BVector::data_type, typename BVector::array_layout, typename BVector::execution_space, typename BVector::memory_traits > BView;
368 
369  MV_ElementWiseMultiplyFunctor<CView,AView,BView> op(c,C,ab,A,B,C.extent(1)) ;
370  Kokkos::parallel_for( C.extent(0) , op );
371  return C;
372 }
373 
374 // Rank-2 vector multiply with Sacado::UQ::PCE scalar type, non-constant c, ab
375 template <typename CS, typename CL, typename CD, typename CM,
376  typename AS, typename AL, typename AD, typename AM,
377  typename BS, typename BL, typename BD, typename BM>
378 Kokkos::View< Sacado::UQ::PCE<CS>**, CL, CD, CM>
379 MV_ElementWiseMultiply(
380  const Sacado::UQ::PCE<CS>& c,
381  const Kokkos::View< Sacado::UQ::PCE<CS>**, CL, CD, CM >& C,
382  const Sacado::UQ::PCE<AS>& ab,
383  const Kokkos::View< Sacado::UQ::PCE<AS>*, AL, AD, AM >& A,
384  const Kokkos::View< Sacado::UQ::PCE<BS>**, BL, BD, BM >& B )
385 {
386  if (Sacado::is_constant(c) && Sacado::is_constant(ab)) {
387  return MV_ElementWiseMultiply( c.fastAccessCoeff(0), C,
388  ab.fastAccessCoeff(0), A, B );
389  }
390  else {
391  Impl::raise_error(
392  "MV_ElementWiseMultiply not implemented for non-constant c or ab");
393  }
394  return C;
395 }
396 
397 // Rank-1 vector scale with Sacado::UQ::PCE scalar type, constant a, b
398 template <typename RS, typename RL, typename RD, typename RM,
399  typename XS, typename XL, typename XD, typename XM>
400 Kokkos::View< Sacado::UQ::PCE<RS>*, RL, RD, RM>
401 V_MulScalar( const Kokkos::View< Sacado::UQ::PCE<RS>*, RL, RD, RM >& r,
402  const typename Sacado::UQ::PCE<XS>::value_type& a,
403  const Kokkos::View< Sacado::UQ::PCE<XS>*, XL, XD, XM >& x )
404 {
405  typedef Kokkos::View< Sacado::UQ::PCE<RS>*, RL, RD, RM > RVector;
406  typedef Kokkos::View< Sacado::UQ::PCE<XS>*, XL, XD, XM > XVector;
407 
408  typename RVector::flat_array_type r_flat = r;
409  typename XVector::flat_array_type x_flat = x;
410 
411  V_MulScalar( r_flat, a, x_flat );
412 
413  return r;
414 }
415 
416 // Rank-1 vector scale with Sacado::UQ::PCE scalar type, non-constant a, b
417 template <typename RS, typename RL, typename RD, typename RM,
418  typename XS, typename XL, typename XD, typename XM>
419 Kokkos::View< Sacado::UQ::PCE<RS>*, RL, RD, RM>
420 V_MulScalar( const Kokkos::View< Sacado::UQ::PCE<RS>*, RL, RD, RM >& r,
421  const Sacado::UQ::PCE<XS>& a,
422  const Kokkos::View< Sacado::UQ::PCE<XS>*, XL, XD, XM >& x )
423 {
424  if (Sacado::is_constant(a)) {
425  return V_MulScalar( r, a.fastAccessCoeff(0), x );
426  }
427  else {
428  Impl::raise_error("V_MulScalar not implemented for non-constant a");
429  }
430  return r;
431 }
432 
433 // Rank-2 vector scale with Sacado::UQ::PCE scalar type, constant a, b
434 template <typename RS, typename RL, typename RD, typename RM,
435  typename XS, typename XL, typename XD, typename XM>
436 Kokkos::View< Sacado::UQ::PCE<RS>**, RL, RD, RM>
437 MV_MulScalar( const Kokkos::View< Sacado::UQ::PCE<RS>**, RL, RD, RM >& r,
438  const typename Sacado::UQ::PCE<XS>::value_type& a,
439  const Kokkos::View< Sacado::UQ::PCE<XS>**, XL, XD, XM >& x )
440 {
441  typedef Kokkos::View< Sacado::UQ::PCE<RS>**, RL, RD, RM > RVector;
442  typedef Kokkos::View< Sacado::UQ::PCE<XS>**, XL, XD, XM > XVector;
443 
444  typename RVector::flat_array_type r_flat = r;
445  typename XVector::flat_array_type x_flat = x;
446 
447  MV_MulScalar( r_flat, a, x_flat );
448 
449  return r;
450 }
451 
452 // Rank-2 vector scale with Sacado::UQ::PCE scalar type, non-constant a, b
453 template <typename RS, typename RL, typename RD, typename RM,
454  typename XS, typename XL, typename XD, typename XM>
455 Kokkos::View< Sacado::UQ::PCE<RS>**, RL, RD, RM>
456 MV_MulScalar( const Kokkos::View< Sacado::UQ::PCE<RS>**, RL, RD, RM >& r,
457  const Sacado::UQ::PCE<XS>& a,
458  const Kokkos::View< Sacado::UQ::PCE<XS>**, XL, XD, XM >& x )
459 {
460  if (Sacado::is_constant(a)) {
461  return MV_MulScalar( r, a.fastAccessCoeff(0), x );
462  }
463  else {
464  Impl::raise_error("MV_MulScalar not implemented for non-constant a or b");
465  }
466  return r;
467 }
468 
469 template <typename T>
470 struct V_ReciprocalThresholdSelfFunctor;
471 
472 template <typename T, typename D, typename M>
473 struct V_ReciprocalThresholdSelfFunctor<
474  View< T,LayoutLeft,D,M,Impl::ViewPCEContiguous > >
475 {
476  typedef View< T,LayoutLeft,D,M,Impl::ViewPCEContiguous > XVector;
477  typedef typename XVector::array_type array_type;
478 
479  typedef typename array_type::execution_space execution_space;
480  typedef typename array_type::size_type size_type;
481  typedef typename array_type::non_const_value_type value_type;
482  typedef Kokkos::Details::ArithTraits<value_type> KAT;
483  typedef typename KAT::mag_type mag_type;
484 
485  const array_type m_x;
486  const value_type m_min_val;
487  const value_type m_min_val_mag;
488  const size_type m_n_pce;
489 
490  V_ReciprocalThresholdSelfFunctor(
491  const XVector& x,
492  const typename XVector::non_const_value_type& min_val) :
493  m_x(x),
494  m_min_val(min_val.fastAccessCoeff(0)),
495  m_min_val_mag(KAT::abs(m_min_val)),
496  m_n_pce(x.sacado_size()) {}
497  //--------------------------------------------------------------------------
498 
499  KOKKOS_INLINE_FUNCTION
500  void operator()( const size_type i) const
501  {
502  if (KAT::abs(m_x(0,i)) < m_min_val_mag)
503  m_x(0,i) = m_min_val;
504  else
505  m_x(0,i) = KAT::one() / m_x(0,i);
506  for (size_type l=1; l<m_n_pce; ++l)
507  m_x(l,i) = KAT::zero();
508  }
509 };
510 
511 template <typename T>
512 struct MV_ReciprocalThresholdSelfFunctor;
513 
514 template <typename T, typename D, typename M>
515 struct MV_ReciprocalThresholdSelfFunctor<
516  View< T,LayoutLeft,D,M,Impl::ViewPCEContiguous > >
517 {
518  typedef View< T,LayoutLeft,D,M,Impl::ViewPCEContiguous > XVector;
519  typedef typename XVector::array_type array_type;
520 
521  typedef typename array_type::execution_space execution_space;
522  typedef typename array_type::size_type size_type;
523  typedef typename array_type::non_const_value_type value_type;
524  typedef Kokkos::Details::ArithTraits<value_type> KAT;
525  typedef typename KAT::mag_type mag_type;
526 
527  const array_type m_x;
528  const value_type m_min_val;
529  const value_type m_min_val_mag;
530  const size_type m_n;
531  const size_type m_n_pce;
532 
533  MV_ReciprocalThresholdSelfFunctor(
534  const XVector& x,
535  const typename XVector::non_const_value_type& min_val,
536  const size_type& n) :
537  m_x(x),
538  m_min_val(min_val.fastAccessCoeff(0)),
539  m_min_val_mag(KAT::abs(m_min_val)),
540  m_n(n),
541  m_n_pce(x.sacado_size()) {}
542  //--------------------------------------------------------------------------
543 
544  KOKKOS_INLINE_FUNCTION
545  void operator()( const size_type i) const
546  {
547  for (size_type k=0; k<m_n; ++k) {
548  if (KAT::abs(m_x(0,i,k)) < m_min_val_mag)
549  m_x(0,i,k) = m_min_val;
550  else
551  m_x(0,i,k) = KAT::one() / m_x(0,i,k);
552  for (size_type l=1; l<m_n_pce; ++l)
553  m_x(l,i,k) = KAT::zero();
554  }
555  }
556 };
557 
558 } // namespace Kokkos
559 */
560 #endif /* #ifndef KOKKOS_MV_UQ_PCE_HPP */
Kokkos_InnerProductSpaceTraits_UQ_PCE.hpp
Kokkos_Blas1_UQ_PCE.hpp
Sacado_UQ_PCE.hpp
Kokkos_View_UQ_PCE.hpp
Generated by   doxygen 1.8.5