42 #ifndef STOKHOS_UNIT_TEST_HELPERS_HPP 43 #define STOKHOS_UNIT_TEST_HELPERS_HPP 48 #include "Teuchos_SerialDenseMatrix.hpp" 49 #include "Teuchos_FancyOStream.hpp" 53 template<
class OrdinalType,
class ValueType>
55 const std::string& a1_name,
57 const std::string& a2_name,
58 const ValueType& rel_tol,
const ValueType& abs_tol,
59 Teuchos::FancyOStream& out)
63 out <<
"Comparing " << a1_name <<
" == " << a2_name <<
" ... ";
65 const OrdinalType n = a1.
size();
69 out <<
"\nError, "<<a1_name<<
".size() = "<<a1.
size()<<
" == " 70 << a2_name<<
".size() = "<<a2.
size()<<
" : failed!\n";
75 for( OrdinalType i = 0; i < n; ++i ) {
77 ValueType err =
std::abs(a1[i] - a2[i]) / nrm;
82 <<
"\nError, relErr("<<a1_name<<
"["<<i<<
"]," 83 <<a2_name<<
"["<<i<<
"]) = relErr("<<a1[i]<<
","<<a2[i]<<
") = " 84 <<err<<
" <= tol = "<<tol<<
": failed!\n";
93 << a1_name <<
" = " << a1 << std::endl
94 << a2_name <<
" = " << a2 << std::endl;
100 template<
class ValueType>
102 const std::string& a1_name,
104 const std::string& a2_name,
105 const ValueType& rel_tol,
const ValueType& abs_tol,
106 Teuchos::FancyOStream& out)
113 out <<
"\nError, relErr(" << a1_name <<
"," 114 << a2_name <<
") = relErr(" << a1 <<
"," << a2 <<
") = " 115 << err <<
" <= tol = " << tol <<
": failed!\n";
122 template<
class ValueType,
class VectorType1,
class VectorType2>
124 const std::string& a1_name,
125 const VectorType2&a2,
126 const std::string& a2_name,
127 const ValueType& rel_tol,
const ValueType& abs_tol,
128 Teuchos::FancyOStream& out)
132 out <<
"Comparing " << a1_name <<
" == " << a2_name <<
" ... ";
134 const int n = a1.size();
137 if (a2.size() != n) {
138 out <<
"\nError, "<<a1_name<<
".size() = "<<a1.size()<<
" == " 139 << a2_name<<
".size() = "<<a2.size()<<
" : failed!\n";
144 for(
int i = 0; i < n; ++i ) {
145 ValueType err =
std::abs(a1.coeff(i) - a2.coeff(i));
151 <<
"\nError, relErr("<<a1_name<<
"["<<i<<
"]," 152 <<a2_name<<
"["<<i<<
"]) = relErr("<<a1.coeff(i)<<
","<<a2.coeff(i)
153 <<
") = "<<err<<
" <= tol = "<<tol<<
": failed!\n";
162 << a1_name <<
" = " << a1 << std::endl
163 << a2_name <<
" = " << a2 << std::endl;
169 template<
class Array1,
class Array2,
class ValueType>
171 const Array2& a2,
const std::string& a2_name,
172 const ValueType& rel_tol,
173 const ValueType& abs_tol,
174 Teuchos::FancyOStream& out)
179 out <<
"Comparing " << a1_name <<
" == " << a2_name <<
" ... ";
181 const int n = a1.size();
184 if (as<int>(a2.size()) != n) {
185 out <<
"\nError, "<<a1_name<<
".size() = "<<a1.size()<<
" == " 186 << a2_name<<
".size() = "<<a2.size()<<
" : failed!\n";
191 for(
int i = 0; i < n; ++i ) {
192 ValueType err =
std::abs(a1[i] - a2[i]);
196 out <<
"\nError, relErr(" << a1_name <<
"[" << i <<
"]," << a2_name
197 <<
"[" << i <<
"]) = relErr(" << a1[i] <<
"," <<a2[i] <<
") = " 198 << err <<
" <= tol = " << tol <<
": failed!\n";
209 template<
class ordinal_type,
class scalar_type>
211 const Teuchos::SerialDenseMatrix<ordinal_type,scalar_type>& a1,
212 const std::string& a1_name,
213 const Teuchos::SerialDenseMatrix<ordinal_type,scalar_type>& a2,
214 const std::string& a2_name,
217 Teuchos::FancyOStream& out)
222 out <<
"Comparing " << a1_name <<
" == " << a2_name <<
" ... ";
228 if (a2.numRows() != m) {
229 out <<
"\nError, "<<a1_name<<
".numRows() = "<<a1.numRows()<<
" == " 230 << a2_name<<
".numRows() = "<<a2.numRows()<<
" : failed!\n";
235 if (a2.numCols() != n) {
236 out <<
"\nError, "<<a1_name<<
".numCols() = "<<a1.numCols()<<
" == " 237 << a2_name<<
".numCols() = "<<a2.numCols()<<
" : failed!\n";
248 out <<
"\nError, relErr(" << a1_name <<
"(" << i <<
"," <<
j <<
"), " 249 << a2_name <<
"(" << i <<
"," <<
j <<
")) = relErr(" 250 << a1(i,
j) <<
", " <<a2(i,
j) <<
") = " 251 << err <<
" <= tol = " << tol <<
": failed!\n";
263 template<
class ordinal_type,
class scalar_type>
266 const std::string& cijk1_name,
268 const std::string& cijk2_name,
271 Teuchos::FancyOStream& out)
276 out <<
"Comparing " << cijk1_name <<
" == " << cijk2_name <<
" ... ";
283 for (
typename Cijk_type::k_iterator k_it=Cijk2.
k_begin();
284 k_it!=Cijk2.
k_end(); ++k_it) {
286 for (
typename Cijk_type::kj_iterator j_it = Cijk2.
j_begin(k_it);
287 j_it != Cijk2.
j_end(k_it); ++j_it) {
289 for (
typename Cijk_type::kji_iterator i_it = Cijk2.
i_begin(j_it);
290 i_it != Cijk2.
i_end(j_it); ++i_it) {
295 double tol = abs_tol + c2*rel_tol;
300 <<
"Check: rel_err( C(" << i <<
"," <<
j <<
"," << k <<
") )" 301 <<
" = " <<
"rel_err( " << c1 <<
", " << c2 <<
" ) = " << err
302 <<
" <= " << tol <<
" : ";
303 if (s) out <<
"Passed.";
308 success = success && s;
316 template<
class ordinal_type,
class scalar_type>
323 Teuchos::FancyOStream& out,
328 Teuchos::Array< Teuchos::RCP<const Stokhos::OneDOrthogPolyBasis<ordinal_type,scalar_type> > > bases = basis.
getCoordinateBases();
330 Teuchos::Array< Teuchos::RCP<Stokhos::Dense3Tensor<ordinal_type,scalar_type> > > Cijk_1d(d);
332 Cijk_1d[i] = bases[i]->computeTripleProductTensor();
342 c *= (*Cijk_1d[l])(terms_i[l],terms_j[l],terms_k[l]);
355 rel_tol, abs_tol, out);
360 template <
typename operator_type1,
typename operator_type2,
363 const operator_type2& op2,
366 Teuchos::FancyOStream& out) {
369 typedef typename operator_type1::const_set_iterator point_iterator_type;
372 TEUCHOS_TEST_EQUALITY(op1.point_size(), op2.point_size(), out, success);
373 if (!success)
return false;
376 point_iterator_type pi1 = op1.set_begin();
377 point_iterator_type pi2 = op2.set_begin();
378 while (pi1 != op1.set_end() || pi2 != op2.set_end()) {
379 std::stringstream ss1, ss2;
384 pi2->first, ss2.str(),
385 rel_tol, abs_tol, out);
387 std::stringstream ss3, ss4;
388 ss3 << pi1->second.first;
389 ss4 << pi2->second.first;
392 pi2->second.first, ss4.str(),
393 rel_tol, abs_tol, out);
401 template <
typename operator_type1,
typename operator_type2,
402 typename func_type1,
typename func_type2,
405 const operator_type2& op2,
406 const func_type1& func1,
407 const func_type2& func2,
410 Teuchos::FancyOStream& out) {
415 typedef typename operator_type1::const_iterator point_iterator_type;
418 TEUCHOS_TEST_EQUALITY(op1.coeff_size(), op2.coeff_size(), out, success);
427 Teuchos::SerialDenseMatrix<ordinal_type,value_type>
f(point_sz1,2);
429 for (point_iterator_type pi = op1.begin(); pi != op1.end(); ++pi) {
430 f(idx,0) = func1(*pi);
431 f(idx,1) = func2(*pi);
435 Teuchos::SerialDenseMatrix<ordinal_type,value_type> f2(point_sz2,2);
437 for (point_iterator_type pi = op2.begin(); pi != op2.end(); ++pi) {
438 f2(idx,0) = func1(*pi);
439 f2(idx,1) = func2(*pi);
444 if (point_sz1 == point_sz2)
449 Teuchos::SerialDenseMatrix<ordinal_type,value_type>
x(coeff_sz,2);
450 op1.transformQP2PCE(1.0,
f,
x, 0.0);
452 Teuchos::SerialDenseMatrix<ordinal_type,value_type> x2(coeff_sz,2);
453 op2.transformQP2PCE(1.0, f2, x2, 0.0);
462 template <
typename operator_type1,
typename operator_type2,
463 typename func_type1,
typename func_type2,
466 const operator_type2& op2,
467 const func_type1& func1,
468 const func_type2& func2,
471 Teuchos::FancyOStream& out) {
476 typedef typename operator_type1::const_iterator point_iterator_type;
479 TEUCHOS_TEST_EQUALITY(op1.coeff_size(), op2.coeff_size(), out, success);
488 Teuchos::SerialDenseMatrix<ordinal_type,value_type>
f(2,point_sz1);
490 for (point_iterator_type pi = op1.begin(); pi != op1.end(); ++pi) {
491 f(0,idx) = func1(*pi);
492 f(1,idx) = func2(*pi);
496 Teuchos::SerialDenseMatrix<ordinal_type,value_type> f2(2,point_sz2);
498 for (point_iterator_type pi = op2.begin(); pi != op2.end(); ++pi) {
499 f2(0,idx) = func1(*pi);
500 f2(1,idx) = func2(*pi);
505 if (point_sz1 == point_sz2)
510 Teuchos::SerialDenseMatrix<ordinal_type,value_type>
x(2,coeff_sz);
511 op1.transformQP2PCE(1.0,
f,
x, 0.0,
true);
513 Teuchos::SerialDenseMatrix<ordinal_type,value_type> x2(2,coeff_sz);
514 op2.transformQP2PCE(1.0, f2, x2, 0.0,
true);
523 template <
typename basis_type,
typename operator_type,
typename scalar_type>
525 const operator_type& op,
528 Teuchos::FancyOStream& out) {
534 Teuchos::SerialDenseMatrix<ordinal_type,value_type> eye(coeff_sz,coeff_sz);
540 Teuchos::SerialDenseMatrix<ordinal_type,value_type>
f(point_sz,coeff_sz);
541 op.transformPCE2QP(1.0, eye,
f, 0.0);
544 Teuchos::SerialDenseMatrix<ordinal_type,value_type>
x(coeff_sz,coeff_sz);
545 op.transformQP2PCE(1.0,
f,
x, 0.0);
552 Teuchos::SerialDenseMatrix<ordinal_type,value_type> z(coeff_sz,coeff_sz);
555 out <<
"Discrete orthogonality error = " <<
x.normInf() << std::endl;
610 #endif // STOKHOS_UNIT_TEST_HELPERS_HPP
KOKKOS_INLINE_FUNCTION PCE< Storage > sqrt(const PCE< Storage > &a)
kji_iterator i_end(const kj_iterator &j) const
Iterator pointing to last i entry for given j and k.
bool compareVecs(const VectorType1 &a1, const std::string &a1_name, const VectorType2 &a2, const std::string &a2_name, const ValueType &rel_tol, const ValueType &abs_tol, Teuchos::FancyOStream &out)
Data structure storing a sparse 3-tensor C(i,j,k) in a a compressed format.
bool compareArrays(const Array1 &a1, const std::string &a1_name, const Array2 &a2, const std::string &a2_name, const ValueType &rel_tol, const ValueType &abs_tol, Teuchos::FancyOStream &out)
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)
void fillComplete()
Signal all terms have been added.
bool compareValues(const ValueType &a1, const std::string &a1_name, const ValueType &a2, const std::string &a2_name, const ValueType &rel_tol, const ValueType &abs_tol, Teuchos::FancyOStream &out)
virtual ordinal_type dimension() const =0
Return dimension of basis.
ordinal_type order() const
Compute total order of index.
ordinal_type size() const
Return size.
bool compareSparse3Tensor(const Stokhos::Sparse3Tensor< ordinal_type, scalar_type > &Cijk1, const std::string &cijk1_name, const Stokhos::Sparse3Tensor< ordinal_type, scalar_type > &Cijk2, const std::string &cijk2_name, const scalar_type &rel_tol, const scalar_type &abs_tol, Teuchos::FancyOStream &out)
bool testSparse3Tensor(const Stokhos::Sparse3Tensor< ordinal_type, scalar_type > &Cijk, const Stokhos::ProductBasis< ordinal_type, scalar_type > &basis, const scalar_type &sparse_tol, const scalar_type &rel_tol, const scalar_type &abs_tol, Teuchos::FancyOStream &out, bool linear=false)
virtual const MultiIndex< ordinal_type > & term(ordinal_type i) const =0
Get orders of each coordinate polynomial given an index i.
const IndexType const IndexType const IndexType const IndexType const ValueType const ValueType * x
bool compareSDM(const Teuchos::SerialDenseMatrix< ordinal_type, scalar_type > &a1, const std::string &a1_name, const Teuchos::SerialDenseMatrix< ordinal_type, scalar_type > &a2, const std::string &a2_name, const scalar_type &rel_tol, const scalar_type &abs_tol, Teuchos::FancyOStream &out)
void add_term(ordinal_type i, ordinal_type j, ordinal_type k, const value_type &c)
Add new term for given (i,j,k)
KOKKOS_INLINE_FUNCTION PCE< Storage > max(const typename PCE< Storage >::value_type &a, const PCE< Storage > &b)
Top-level namespace for Stokhos classes and functions.
value_type getValue(ordinal_type i, ordinal_type j, ordinal_type k) const
Get Cijk value for a given i, j, k indices.
ordinal_type num_entries() const
Return number of non-zero entries.
KOKKOS_INLINE_FUNCTION PCE< Storage > abs(const PCE< Storage > &a)
Abstract base class for multivariate orthogonal polynomials generated from tensor products of univari...
virtual const Teuchos::Array< value_type > & norm_squared() const =0
Return array storing norm-squared of each basis polynomial.
virtual Teuchos::Array< Teuchos::RCP< const OneDOrthogPolyBasis< ordinal_type, value_type > > > getCoordinateBases() const =0
Return array of coordinate bases.
Legendre polynomial basis.
Stokhos::Sparse3Tensor< int, double > Cijk_type
bool testPseudoSpectralDiscreteOrthogonality(const basis_type &basis, const operator_type &op, const scalar_type &rel_tol, const scalar_type &abs_tol, Teuchos::FancyOStream &out)
ScalarType f(const Teuchos::Array< ScalarType > &x, double a, double b)
k_iterator k_end() const
Iterator pointing to last k entry.
kj_iterator j_begin(const k_iterator &k) const
Iterator pointing to first j entry for given k.
kj_iterator j_end(const k_iterator &k) const
Iterator pointing to last j entry for given k.
bool testPseudoSpectralPoints(const operator_type1 &op1, const operator_type2 &op2, const scalar_type &rel_tol, const scalar_type &abs_tol, Teuchos::FancyOStream &out)
bool testPseudoSpectralApply(const operator_type1 &op1, const operator_type2 &op2, const func_type1 &func1, const func_type2 &func2, const scalar_type &rel_tol, const scalar_type &abs_tol, Teuchos::FancyOStream &out)
Teuchos::RCP< const Stokhos::OrthogPolyBasis< ordinal_type, value_type > > basis() const
Return basis.
kji_iterator i_begin(const kj_iterator &j) const
Iterator pointing to first i entry for given j and k.
virtual ordinal_type size() const =0
Return total size of basis.
bool testPseudoSpectralApplyTrans(const operator_type1 &op1, const operator_type2 &op2, const func_type1 &func1, const func_type2 &func2, const scalar_type &rel_tol, const scalar_type &abs_tol, Teuchos::FancyOStream &out)
k_iterator k_begin() const
Iterator pointing to first k entry.