Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Stokhos_SacadoUQPCEUnitTest.cpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Stokhos Package
4 //
5 // Copyright 2009 NTESS and the Stokhos contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
14 
15 #include "Stokhos.hpp"
18 
19 // Tests only run on the host, so use host execution space
20 typedef Kokkos::DefaultHostExecutionSpace execution_space;
23 
24 namespace Stokhos {
25 
26  template<class PCEType, class OrdinalType, class ValueType>
27  bool comparePCEs(const PCEType& a1,
28  const std::string& a1_name,
30  const std::string& a2_name,
31  const ValueType& rel_tol, const ValueType& abs_tol,
33  {
34  bool success = true;
35 
36  out << "Comparing " << a1_name << " == " << a2_name << " ... ";
37 
38  const OrdinalType n = a1.size();
39 
40  // Compare sizes
41  if (a2.size() != n) {
42  out << "\nError, "<<a1_name<<".size() = "<<a1.size()<<" == "
43  << a2_name<<".size() = "<<a2.size()<<" : failed!\n";
44  return false;
45  }
46 
47  // Compare elements
48  for( OrdinalType i = 0; i < n; ++i ) {
49  ValueType nrm = std::sqrt(a2.basis()->norm_squared(i));
50  ValueType err = std::abs(a1.coeff(i) - a2[i]) / nrm;
51  ValueType tol =
52  abs_tol + rel_tol*std::max(std::abs(a1.coeff(i)),std::abs(a2[i]))/nrm;
53  if (err > tol) {
54  out
55  <<"\nError, relErr("<<a1_name<<"["<<i<<"],"
56  <<a2_name<<"["<<i<<"]) = relErr("<<a1.coeff(i)<<","<<a2[i]<<") = "
57  <<err<<" <= tol = "<<tol<<": failed!\n";
58  success = false;
59  }
60  }
61  if (success) {
62  out << "passed\n";
63  }
64  else {
65  out << std::endl
66  << a1_name << " = " << a1 << std::endl
67  << a2_name << " = " << a2 << std::endl;
68  }
69 
70  return success;
71  }
72 
73 }
74 
75 namespace SacadoPCEUnitTest {
76 
77  // Common setup for unit tests
78  template <typename PCEType>
79  struct UnitTestSetup {
80 
81  typedef PCEType pce_type;
83  typedef typename pce_type::value_type value_type;
85  typedef typename pce_type::cijk_type kokkos_cijk_type;
94  pce_type x, y, sin_x, cos_y, cx, u, u2, cu, cu2, sx, su, su2;
97 
99  rtol = 1e-4;
100  atol = 1e-5;
101  crtol = 1e-12;
102  catol = 1e-12;
103  a = 3.1;
104  const ordinal_type d = 2;
105  const ordinal_type p = 7;
106 
107  // Create product basis
109  for (ordinal_type i=0; i<d; i++)
110  bases[i] =
112  basis =
114 
115  // Triple product tensor
117 
118  // Kokkos triple product tensor
119  cijk = Stokhos::create_product_tensor<execution_space>(*basis, *Cijk);
120 
121  // Algebraic expansion
123 
124  // Quad expansion for initialization
129 
130  // Create approximation
131  x_opa.reset(basis);
132  y_opa.reset(basis);
133  sin_x_opa.reset(basis);
134  cos_y_opa.reset(basis);
135  cx_opa.reset(basis,1);
136  x_opa.term(0, 0) = 1.0;
137  y_opa.term(0, 0) = 2.0;
138  cx_opa.term(0, 0) = a;
139  for (int i=0; i<d; i++) {
140  x_opa.term(i, 1) = 0.1;
141  y_opa.term(i, 1) = 0.25;
142  }
143  quad_exp->sin(sin_x_opa, x_opa);
144  quad_exp->cos(cos_y_opa, y_opa);
145 
146  // Create PCEs
147  x.reset(cijk);
148  y.reset(cijk);
149  sin_x.reset(cijk);
150  cos_y.reset(cijk);
151  cx.reset(cijk, 1);
152  x.load(x_opa.coeff());
153  y.load(y_opa.coeff());
154  sin_x.load(sin_x_opa.coeff());
155  cos_y.load(cos_y_opa.coeff());
156  cx.load(cx_opa.coeff());
157 
158  u.reset(cijk);
159  u2.reset(cijk);
160  cu.reset(cijk);
161  cu2.reset(cijk, 1);
162  sx.reset(cijk, d+1);
163  su.reset(cijk, d+1);
164  su2.reset(cijk, d+1);
165  for (ordinal_type i=0; i<d; i++) {
166  sx.fastAccessCoeff(i+1) = 0.0;
167  }
168  }
169  };
170 
172 
173  TEUCHOS_UNIT_TEST( Stokhos_PCE, UMinus) {
174  UTS setup;
175  UTS::pce_type u = -setup.sin_x;
176  UTS::opa_type u_opa(setup.basis);
177  setup.exp->unaryMinus(u_opa, setup.sin_x_opa);
178  success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa",
179  setup.rtol, setup.atol, out);
180  }
181 
182 
183 #define UNARY_UNIT_TEST(OP) \
184  TEUCHOS_UNIT_TEST( Stokhos_PCE, OP##_const) { \
185  UTS setup; \
186  UTS::pce_type u = OP(setup.cx); \
187  UTS::opa_type u_opa(setup.basis); \
188  setup.exp->OP(u_opa, setup.cx_opa); \
189  success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
190  setup.rtol, setup.atol, out); \
191  } \
192  TEUCHOS_UNIT_TEST( Stokhos_PCE, OP##_resize) { \
193  UTS setup; \
194  UTS::pce_type u; \
195  u = OP(setup.cx); \
196  UTS::opa_type u_opa(setup.basis); \
197  setup.exp->OP(u_opa, setup.cx_opa); \
198  success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
199  setup.rtol, setup.atol, out); \
200  }
201 
216  // UNARY_UNIT_TEST(asinh)
217  // UNARY_UNIT_TEST(acosh)
218  // UNARY_UNIT_TEST(atanh)
219 
220 #define BINARY_UNIT_TEST(OP, EXPOP) \
221  TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP) { \
222  UTS setup; \
223  UTS::pce_type v = setup.sin_x; \
224  UTS::pce_type w = setup.cos_y; \
225  UTS::pce_type u = OP(v,w); \
226  UTS::opa_type u_opa(setup.basis); \
227  setup.exp->EXPOP(u_opa, setup.sin_x_opa, setup.cos_y_opa); \
228  success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
229  setup.rtol, setup.atol, out); \
230  } \
231  TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_left_const) { \
232  UTS setup; \
233  UTS::pce_type w = setup.sin_x; \
234  UTS::pce_type u = OP(setup.a, w); \
235  UTS::opa_type u_opa(setup.basis); \
236  setup.exp->EXPOP(u_opa, setup.a, setup.sin_x_opa); \
237  success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
238  setup.rtol, setup.atol, out); \
239  } \
240  TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_right_const) { \
241  UTS setup; \
242  UTS::pce_type v = setup.sin_x; \
243  UTS::pce_type u = OP(v, setup.a); \
244  UTS::opa_type u_opa(setup.basis); \
245  setup.exp->EXPOP(u_opa, setup.sin_x_opa, setup.a); \
246  success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
247  setup.rtol, setup.atol, out); \
248  } \
249  TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_both_const) { \
250  UTS setup; \
251  UTS::pce_type u = OP(setup.cx, setup.cx); \
252  UTS::opa_type u_opa(setup.basis); \
253  setup.exp->EXPOP(u_opa, setup.cx_opa, setup.cx_opa); \
254  success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
255  setup.rtol, setup.atol, out); \
256  } \
257  TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_left_const2) { \
258  UTS setup; \
259  UTS::pce_type w = setup.sin_x; \
260  UTS::pce_type u = OP(setup.cx, w); \
261  UTS::opa_type u_opa(setup.basis); \
262  setup.exp->EXPOP(u_opa, setup.cx_opa, setup.sin_x_opa); \
263  success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
264  setup.rtol, setup.atol, out); \
265  } \
266  TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_right_const2) { \
267  UTS setup; \
268  UTS::pce_type v = setup.sin_x; \
269  UTS::pce_type u = OP(v, setup.cx); \
270  UTS::opa_type u_opa(setup.basis); \
271  setup.exp->EXPOP(u_opa, setup.sin_x_opa, setup.cx_opa); \
272  success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
273  setup.rtol, setup.atol, out); \
274  } \
275  TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_resize) { \
276  UTS setup; \
277  UTS::pce_type v = setup.sin_x; \
278  UTS::pce_type w = setup.cos_y; \
279  UTS::pce_type u; \
280  u = OP(v, w); \
281  UTS::opa_type u_opa(setup.basis); \
282  setup.exp->EXPOP(u_opa, setup.sin_x_opa, setup.cos_y_opa); \
283  success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
284  setup.rtol, setup.atol, out); \
285  } \
286  TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_left_const_resize) { \
287  UTS setup; \
288  UTS::pce_type w = setup.sin_x; \
289  UTS::pce_type u; \
290  u = OP(setup.a, w); \
291  UTS::opa_type u_opa(setup.basis); \
292  setup.exp->EXPOP(u_opa, setup.a, setup.sin_x_opa); \
293  success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
294  setup.rtol, setup.atol, out); \
295  } \
296  TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_right_const_resize) { \
297  UTS setup; \
298  UTS::pce_type v = setup.sin_x; \
299  UTS::pce_type u; \
300  u = OP(v, setup.a); \
301  UTS::opa_type u_opa(setup.basis); \
302  setup.exp->EXPOP(u_opa, setup.sin_x_opa, setup.a); \
303  success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
304  setup.rtol, setup.atol, out); \
305  }
306 
307  BINARY_UNIT_TEST(operator+, plus)
308  BINARY_UNIT_TEST(operator-, minus)
309  BINARY_UNIT_TEST(operator*, times)
310  BINARY_UNIT_TEST(operator/, divide)
311 
312 #define OPASSIGN_UNIT_TEST(OP, EXPOP) \
313  TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP) { \
314  UTS setup; \
315  UTS::pce_type v = setup.sin_x; \
316  UTS::pce_type u = setup.cos_y; \
317  u OP v; \
318  UTS::opa_type u_opa = setup.cos_y_opa; \
319  setup.exp->EXPOP(u_opa, setup.sin_x_opa); \
320  success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
321  setup.rtol, setup.atol, out); \
322  } \
323  TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_const) { \
324  UTS setup; \
325  UTS::pce_type u = setup.sin_x; \
326  u OP setup.a; \
327  UTS::opa_type u_opa = setup.sin_x_opa; \
328  setup.exp->EXPOP(u_opa, setup.a); \
329  success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
330  setup.rtol, setup.atol, out); \
331  } \
332  TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_const2) { \
333  UTS setup; \
334  UTS::pce_type u = setup.sin_x; \
335  u OP setup.cx; \
336  UTS::opa_type u_opa = setup.sin_x_opa; \
337  setup.exp->EXPOP(u_opa, setup.cx_opa); \
338  success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
339  setup.rtol, setup.atol, out); \
340  } \
341  TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_resize) { \
342  UTS setup; \
343  UTS::pce_type v = setup.sin_x; \
344  UTS::pce_type u = setup.a; \
345  u OP v; \
346  UTS::opa_type u_opa = setup.cx_opa; \
347  setup.exp->EXPOP(u_opa, setup.sin_x_opa); \
348  success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
349  setup.rtol, setup.atol, out); \
350  }
351 
352  OPASSIGN_UNIT_TEST(+=, plusEqual)
353  OPASSIGN_UNIT_TEST(-=, minusEqual)
354  OPASSIGN_UNIT_TEST(*=, timesEqual)
355  OPASSIGN_UNIT_TEST(/=, divideEqual)
356 
357  TEUCHOS_UNIT_TEST( Stokhos_PCE, initializer_list_copy ) {
358  UTS setup;
359  UTS::pce_type u = { 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0 };
360  UTS::opa_type v(setup.basis, 8);
361  for (int i=0; i<8; i++)
362  v[i] = i+1;
363  success = comparePCEs(u, "u", v, "v",
364  setup.rtol, setup.atol, out);
365  }
366  TEUCHOS_UNIT_TEST( Stokhos_PCE, initializer_list_assign ) {
367  UTS setup;
368  UTS::pce_type u;
369  u = { 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0 };
370  UTS::opa_type v(setup.basis, 8);
371  for (int i=0; i<8; i++)
372  v[i] = i+1;
373  success = comparePCEs(u, "u", v, "v",
374  setup.rtol, setup.atol, out);
375  }
376  TEUCHOS_UNIT_TEST( Stokhos_PCE, range_based_for ) {
377  UTS setup;
378  UTS::pce_type u;
379  u.reset(UTS::pce_type::cijk_type(), 8);
380  for (auto& z : u) { z = 3.0; }
381  UTS::opa_type v(setup.basis, 8);
382  for (int i=0; i<8; i++)
383  v[i] = 3.0;
384  success = comparePCEs(u, "u", v, "v",
385  setup.rtol, setup.atol, out);
386  }
387 }
KOKKOS_INLINE_FUNCTION PCE< Storage > sqrt(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > tan(const PCE< Storage > &a)
Stokhos::OrthogPolyApprox< ordinal_type, value_type > opa_type
Stokhos::StandardStorage< int, double > storage_type
KOKKOS_INLINE_FUNCTION PCE< Storage > sinh(const PCE< Storage > &a)
Sacado::ETPCE::OrthogPoly< double, Stokhos::StandardStorage< int, double > > pce_type
Teuchos::RCP< const Stokhos::Quadrature< int, double > > quad
Teuchos::RCP< Stokhos::Sparse3Tensor< ordinal_type, value_type > > Cijk
#define UNARY_UNIT_TEST(OP)
Kokkos::DefaultExecutionSpace execution_space
Stokhos::OrthogPolyApprox< int, double > opa_type
bool comparePCEs(const PCEType &a1, const std::string &a1_name, const Stokhos::OrthogPolyApprox< OrdinalType, ValueType > &a2, const std::string &a2_name, const ValueType &rel_tol, const ValueType &abs_tol, Teuchos::FancyOStream &out)
Teuchos::RCP< Stokhos::Sparse3Tensor< int, double > > Cijk
Orthogonal polynomial expansions limited to algebraic operations.
pointer coeff()
Return coefficient array.
KOKKOS_INLINE_FUNCTION PCE< Storage > tanh(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > cbrt(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > acos(const PCE< Storage > &a)
UnitTestSetup< pce_type > UTS
Teuchos::RCP< Stokhos::AlgebraicOrthogPolyExpansion< ordinal_type, value_type > > exp
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
void reset(const Teuchos::RCP< const Stokhos::OrthogPolyBasis< ordinal_type, value_type > > &new_basis, ordinal_type sz=0)
Reset to a new basis.
KOKKOS_INLINE_FUNCTION PCE< Storage > max(const typename PCE< Storage >::value_type &a, const PCE< Storage > &b)
TEUCHOS_UNIT_TEST(Stokhos_PCE, UMinus)
virtual Teuchos::RCP< Stokhos::Sparse3Tensor< ordinal_type, value_type > > computeTripleProductTensor() const
Compute triple product tensor.
KOKKOS_INLINE_FUNCTION PCE< Storage > cosh(const PCE< Storage > &a)
Teuchos::RCP< const Stokhos::OrthogPolyBasis< ordinal_type, value_type > > basis() const
Return basis.
#define BINARY_UNIT_TEST(OP, EXPOP)
KOKKOS_INLINE_FUNCTION PCE< Storage > abs(const PCE< Storage > &a)
Multivariate orthogonal polynomial basis generated from a total-order complete-polynomial tensor prod...
virtual const Teuchos::Array< value_type > & norm_squared() const =0
Return array storing norm-squared of each basis polynomial.
KOKKOS_INLINE_FUNCTION PCE< Storage > atan(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > exp(const PCE< Storage > &a)
Legendre polynomial basis.
Teuchos::RCP< const Stokhos::CompletePolynomialBasis< ordinal_type, value_type > > basis
#define OPASSIGN_UNIT_TEST(OP, EXPOP)
KOKKOS_INLINE_FUNCTION PCE< Storage > sin(const PCE< Storage > &a)
ordinal_type size() const
Return size.
Teuchos::RCP< Stokhos::QuadOrthogPolyExpansion< int, double > > exp
Teuchos::RCP< const Stokhos::CompletePolynomialBasis< int, double > > basis
KOKKOS_INLINE_FUNCTION PCE< Storage > log(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > log10(const PCE< Storage > &a)
int n
KOKKOS_INLINE_FUNCTION PCE< Storage > asin(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > cos(const PCE< Storage > &a)
Defines quadrature for a tensor product basis by tensor products of 1-D quadrature rules...
reference term(ordinal_type dimension, ordinal_type order)
Get coefficient term for given dimension and order.