Stokhos Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
Stokhos_SacadoUQPCEUnitTest.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 "Teuchos_UnitTestHarness.hpp"
43#include "Teuchos_TestingHelpers.hpp"
44#include "Teuchos_UnitTestRepository.hpp"
45#include "Teuchos_GlobalMPISession.hpp"
46
47#include "Stokhos.hpp"
50
51// Tests only run on the host, so use host execution space
52typedef Kokkos::DefaultHostExecutionSpace execution_space;
55
56namespace Stokhos {
57
58 template<class PCEType, class OrdinalType, class ValueType>
59 bool comparePCEs(const PCEType& a1,
60 const std::string& a1_name,
62 const std::string& a2_name,
63 const ValueType& rel_tol, const ValueType& abs_tol,
64 Teuchos::FancyOStream& out)
65 {
66 bool success = true;
67
68 out << "Comparing " << a1_name << " == " << a2_name << " ... ";
69
70 const OrdinalType n = a1.size();
71
72 // Compare sizes
73 if (a2.size() != n) {
74 out << "\nError, "<<a1_name<<".size() = "<<a1.size()<<" == "
75 << a2_name<<".size() = "<<a2.size()<<" : failed!\n";
76 return false;
77 }
78
79 // Compare elements
80 for( OrdinalType i = 0; i < n; ++i ) {
81 ValueType nrm = std::sqrt(a2.basis()->norm_squared(i));
82 ValueType err = std::abs(a1.coeff(i) - a2[i]) / nrm;
83 ValueType tol =
84 abs_tol + rel_tol*std::max(std::abs(a1.coeff(i)),std::abs(a2[i]))/nrm;
85 if (err > tol) {
86 out
87 <<"\nError, relErr("<<a1_name<<"["<<i<<"],"
88 <<a2_name<<"["<<i<<"]) = relErr("<<a1.coeff(i)<<","<<a2[i]<<") = "
89 <<err<<" <= tol = "<<tol<<": failed!\n";
90 success = false;
91 }
92 }
93 if (success) {
94 out << "passed\n";
95 }
96 else {
97 out << std::endl
98 << a1_name << " = " << a1 << std::endl
99 << a2_name << " = " << a2 << std::endl;
100 }
101
102 return success;
103 }
104
105}
106
107namespace SacadoPCEUnitTest {
108
109 // Common setup for unit tests
110 template <typename PCEType>
111 struct UnitTestSetup {
112
114 typedef typename pce_type::ordinal_type ordinal_type;
115 typedef typename pce_type::value_type value_type;
116 typedef typename pce_type::execution_space execution_space;
117 typedef typename pce_type::cijk_type kokkos_cijk_type;
122 Teuchos::RCP<const Stokhos::CompletePolynomialBasis<ordinal_type,value_type> > basis;
123 Teuchos::RCP<Stokhos::Sparse3Tensor<ordinal_type,value_type> > Cijk;
125 Teuchos::RCP< Stokhos::AlgebraicOrthogPolyExpansion<ordinal_type,value_type> > exp;
129
131 rtol = 1e-4;
132 atol = 1e-5;
133 crtol = 1e-12;
134 catol = 1e-12;
135 a = 3.1;
136 const ordinal_type d = 2;
137 const ordinal_type p = 7;
138
139 // Create product basis
140 Teuchos::Array< Teuchos::RCP<const Stokhos::OneDOrthogPolyBasis<ordinal_type,value_type> > > bases(d);
141 for (ordinal_type i=0; i<d; i++)
142 bases[i] =
143 Teuchos::rcp(new Stokhos::LegendreBasis<ordinal_type,value_type>(p, true));
144 basis =
146
147 // Triple product tensor
148 Cijk = basis->computeTripleProductTensor();
149
150 // Kokkos triple product tensor
151 cijk = Stokhos::create_product_tensor<execution_space>(*basis, *Cijk);
152
153 // Algebraic expansion
155
156 // Quad expansion for initialization
157 Teuchos::RCP<const Stokhos::Quadrature<int,double> > quad =
159 Teuchos::RCP< Stokhos::QuadOrthogPolyExpansion<int,double> > quad_exp =
161
162 // Create approximation
167 cx_opa.reset(basis,1);
168 x_opa.term(0, 0) = 1.0;
169 y_opa.term(0, 0) = 2.0;
170 cx_opa.term(0, 0) = a;
171 for (int i=0; i<d; i++) {
172 x_opa.term(i, 1) = 0.1;
173 y_opa.term(i, 1) = 0.25;
174 }
175 quad_exp->sin(sin_x_opa, x_opa);
176 quad_exp->cos(cos_y_opa, y_opa);
177
178 // Create PCEs
179 x.reset(cijk);
180 y.reset(cijk);
181 sin_x.reset(cijk);
182 cos_y.reset(cijk);
183 cx.reset(cijk, 1);
184 x.load(x_opa.coeff());
185 y.load(y_opa.coeff());
186 sin_x.load(sin_x_opa.coeff());
187 cos_y.load(cos_y_opa.coeff());
188 cx.load(cx_opa.coeff());
189
190 u.reset(cijk);
191 u2.reset(cijk);
192 cu.reset(cijk);
193 cu2.reset(cijk, 1);
194 sx.reset(cijk, d+1);
195 su.reset(cijk, d+1);
196 su2.reset(cijk, d+1);
197 for (ordinal_type i=0; i<d; i++) {
198 sx.fastAccessCoeff(i+1) = 0.0;
199 }
200 }
201 };
202
203 typedef UnitTestSetup<pce_type> UTS;
204
205 TEUCHOS_UNIT_TEST( Stokhos_PCE, UMinus) {
206 UTS setup;
207 UTS::pce_type u = -setup.sin_x;
208 UTS::opa_type u_opa(setup.basis);
209 setup.exp->unaryMinus(u_opa, setup.sin_x_opa);
210 success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa",
211 setup.rtol, setup.atol, out);
212 }
213
214
215#define UNARY_UNIT_TEST(OP) \
216 TEUCHOS_UNIT_TEST( Stokhos_PCE, OP##_const) { \
217 UTS setup; \
218 UTS::pce_type u = OP(setup.cx); \
219 UTS::opa_type u_opa(setup.basis); \
220 setup.exp->OP(u_opa, setup.cx_opa); \
221 success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
222 setup.rtol, setup.atol, out); \
223 } \
224 TEUCHOS_UNIT_TEST( Stokhos_PCE, OP##_resize) { \
225 UTS setup; \
226 UTS::pce_type u; \
227 u = OP(setup.cx); \
228 UTS::opa_type u_opa(setup.basis); \
229 setup.exp->OP(u_opa, setup.cx_opa); \
230 success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
231 setup.rtol, setup.atol, out); \
232 }
233
234 UNARY_UNIT_TEST(exp)
235 UNARY_UNIT_TEST(log)
236 UNARY_UNIT_TEST(log10)
237 UNARY_UNIT_TEST(sqrt)
238 UNARY_UNIT_TEST(cbrt)
239 UNARY_UNIT_TEST(sin)
240 UNARY_UNIT_TEST(cos)
241 UNARY_UNIT_TEST(tan)
242 UNARY_UNIT_TEST(sinh)
243 UNARY_UNIT_TEST(cosh)
244 UNARY_UNIT_TEST(tanh)
245 UNARY_UNIT_TEST(asin)
246 UNARY_UNIT_TEST(acos)
247 UNARY_UNIT_TEST(atan)
248 // UNARY_UNIT_TEST(asinh)
249 // UNARY_UNIT_TEST(acosh)
250 // UNARY_UNIT_TEST(atanh)
251
252#define BINARY_UNIT_TEST(OP, EXPOP) \
253 TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP) { \
254 UTS setup; \
255 UTS::pce_type v = setup.sin_x; \
256 UTS::pce_type w = setup.cos_y; \
257 UTS::pce_type u = OP(v,w); \
258 UTS::opa_type u_opa(setup.basis); \
259 setup.exp->EXPOP(u_opa, setup.sin_x_opa, setup.cos_y_opa); \
260 success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
261 setup.rtol, setup.atol, out); \
262 } \
263 TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_left_const) { \
264 UTS setup; \
265 UTS::pce_type w = setup.sin_x; \
266 UTS::pce_type u = OP(setup.a, w); \
267 UTS::opa_type u_opa(setup.basis); \
268 setup.exp->EXPOP(u_opa, setup.a, setup.sin_x_opa); \
269 success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
270 setup.rtol, setup.atol, out); \
271 } \
272 TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_right_const) { \
273 UTS setup; \
274 UTS::pce_type v = setup.sin_x; \
275 UTS::pce_type u = OP(v, setup.a); \
276 UTS::opa_type u_opa(setup.basis); \
277 setup.exp->EXPOP(u_opa, setup.sin_x_opa, setup.a); \
278 success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
279 setup.rtol, setup.atol, out); \
280 } \
281 TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_both_const) { \
282 UTS setup; \
283 UTS::pce_type u = OP(setup.cx, setup.cx); \
284 UTS::opa_type u_opa(setup.basis); \
285 setup.exp->EXPOP(u_opa, setup.cx_opa, setup.cx_opa); \
286 success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
287 setup.rtol, setup.atol, out); \
288 } \
289 TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_left_const2) { \
290 UTS setup; \
291 UTS::pce_type w = setup.sin_x; \
292 UTS::pce_type u = OP(setup.cx, w); \
293 UTS::opa_type u_opa(setup.basis); \
294 setup.exp->EXPOP(u_opa, setup.cx_opa, setup.sin_x_opa); \
295 success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
296 setup.rtol, setup.atol, out); \
297 } \
298 TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_right_const2) { \
299 UTS setup; \
300 UTS::pce_type v = setup.sin_x; \
301 UTS::pce_type u = OP(v, setup.cx); \
302 UTS::opa_type u_opa(setup.basis); \
303 setup.exp->EXPOP(u_opa, setup.sin_x_opa, setup.cx_opa); \
304 success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
305 setup.rtol, setup.atol, out); \
306 } \
307 TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_resize) { \
308 UTS setup; \
309 UTS::pce_type v = setup.sin_x; \
310 UTS::pce_type w = setup.cos_y; \
311 UTS::pce_type u; \
312 u = OP(v, w); \
313 UTS::opa_type u_opa(setup.basis); \
314 setup.exp->EXPOP(u_opa, setup.sin_x_opa, setup.cos_y_opa); \
315 success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
316 setup.rtol, setup.atol, out); \
317 } \
318 TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_left_const_resize) { \
319 UTS setup; \
320 UTS::pce_type w = setup.sin_x; \
321 UTS::pce_type u; \
322 u = OP(setup.a, w); \
323 UTS::opa_type u_opa(setup.basis); \
324 setup.exp->EXPOP(u_opa, setup.a, setup.sin_x_opa); \
325 success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
326 setup.rtol, setup.atol, out); \
327 } \
328 TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_right_const_resize) { \
329 UTS setup; \
330 UTS::pce_type v = setup.sin_x; \
331 UTS::pce_type u; \
332 u = OP(v, setup.a); \
333 UTS::opa_type u_opa(setup.basis); \
334 setup.exp->EXPOP(u_opa, setup.sin_x_opa, setup.a); \
335 success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
336 setup.rtol, setup.atol, out); \
337 }
338
339 BINARY_UNIT_TEST(operator+, plus)
340 BINARY_UNIT_TEST(operator-, minus)
341 BINARY_UNIT_TEST(operator*, times)
342 BINARY_UNIT_TEST(operator/, divide)
343
344#define OPASSIGN_UNIT_TEST(OP, EXPOP) \
345 TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP) { \
346 UTS setup; \
347 UTS::pce_type v = setup.sin_x; \
348 UTS::pce_type u = setup.cos_y; \
349 u OP v; \
350 UTS::opa_type u_opa = setup.cos_y_opa; \
351 setup.exp->EXPOP(u_opa, setup.sin_x_opa); \
352 success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
353 setup.rtol, setup.atol, out); \
354 } \
355 TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_const) { \
356 UTS setup; \
357 UTS::pce_type u = setup.sin_x; \
358 u OP setup.a; \
359 UTS::opa_type u_opa = setup.sin_x_opa; \
360 setup.exp->EXPOP(u_opa, setup.a); \
361 success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
362 setup.rtol, setup.atol, out); \
363 } \
364 TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_const2) { \
365 UTS setup; \
366 UTS::pce_type u = setup.sin_x; \
367 u OP setup.cx; \
368 UTS::opa_type u_opa = setup.sin_x_opa; \
369 setup.exp->EXPOP(u_opa, setup.cx_opa); \
370 success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
371 setup.rtol, setup.atol, out); \
372 } \
373 TEUCHOS_UNIT_TEST( Stokhos_PCE, EXPOP##_resize) { \
374 UTS setup; \
375 UTS::pce_type v = setup.sin_x; \
376 UTS::pce_type u = setup.a; \
377 u OP v; \
378 UTS::opa_type u_opa = setup.cx_opa; \
379 setup.exp->EXPOP(u_opa, setup.sin_x_opa); \
380 success = Stokhos::comparePCEs(u, "u", u_opa, "u_opa", \
381 setup.rtol, setup.atol, out); \
382 }
383
384 OPASSIGN_UNIT_TEST(+=, plusEqual)
385 OPASSIGN_UNIT_TEST(-=, minusEqual)
386 OPASSIGN_UNIT_TEST(*=, timesEqual)
387 OPASSIGN_UNIT_TEST(/=, divideEqual)
388
389 TEUCHOS_UNIT_TEST( Stokhos_PCE, initializer_list_copy ) {
390 UTS setup;
391 UTS::pce_type u = { 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0 };
392 UTS::opa_type v(setup.basis, 8);
393 for (int i=0; i<8; i++)
394 v[i] = i+1;
395 success = comparePCEs(u, "u", v, "v",
396 setup.rtol, setup.atol, out);
397 }
398 TEUCHOS_UNIT_TEST( Stokhos_PCE, initializer_list_assign ) {
399 UTS setup;
400 UTS::pce_type u;
401 u = { 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0 };
402 UTS::opa_type v(setup.basis, 8);
403 for (int i=0; i<8; i++)
404 v[i] = i+1;
405 success = comparePCEs(u, "u", v, "v",
406 setup.rtol, setup.atol, out);
407 }
408 TEUCHOS_UNIT_TEST( Stokhos_PCE, range_based_for ) {
409 UTS setup;
410 UTS::pce_type u;
411 u.reset(UTS::pce_type::cijk_type(), 8);
412 for (auto& z : u) { z = 3.0; }
413 UTS::opa_type v(setup.basis, 8);
414 for (int i=0; i<8; i++)
415 v[i] = 3.0;
416 success = comparePCEs(u, "u", v, "v",
417 setup.rtol, setup.atol, out);
418 }
419}
#define UNARY_UNIT_TEST(VEC, SCALAR_T, OP, OPNAME, USING_OP)
#define OPASSIGN_UNIT_TEST(VEC, SCALAR_T, OP, OPNAME)
#define BINARY_UNIT_TEST(VEC, SCALAR_T, OP, OPNAME)
Sacado::UQ::PCE< storage_type > pce_type
Stokhos::DynamicStorage< int, double, execution_space > storage_type
Kokkos::DefaultHostExecutionSpace execution_space
Orthogonal polynomial expansions limited to algebraic operations.
Multivariate orthogonal polynomial basis generated from a total-order complete-polynomial tensor prod...
Legendre polynomial basis.
Class to store coefficients of a projection onto an orthogonal polynomial basis.
void reset(const Teuchos::RCP< const Stokhos::OrthogPolyBasis< ordinal_type, value_type > > &new_basis, ordinal_type sz=0)
Reset to a new basis.
pointer coeff()
Return coefficient array.
reference term(ordinal_type dimension, ordinal_type order)
Get coefficient term for given dimension and order.
Orthogonal polynomial expansions based on numerical quadrature.
Defines quadrature for a tensor product basis by tensor products of 1-D quadrature rules.
TEUCHOS_UNIT_TEST(Stokhos_PCE, UMinus)
UnitTestSetup< pce_type > UTS
Top-level namespace for Stokhos classes and functions.
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< const Stokhos::Quadrature< int, double > > quad
Teuchos::RCP< Stokhos::Sparse3Tensor< ordinal_type, value_type > > Cijk
Teuchos::RCP< Stokhos::Sparse3Tensor< int, double > > Cijk
Teuchos::RCP< const Stokhos::CompletePolynomialBasis< int, double > > basis
Teuchos::RCP< const Stokhos::CompletePolynomialBasis< ordinal_type, value_type > > basis
Stokhos::OrthogPolyApprox< ordinal_type, value_type > opa_type
Teuchos::RCP< Stokhos::AlgebraicOrthogPolyExpansion< ordinal_type, value_type > > exp