65 const unsigned int true_quad_order = 200;
67 Teuchos::Array<double> true_quad_points, true_quad_weights;
68 Teuchos::Array< Teuchos::Array<double> > true_quad_values;
70 true_quad_weights, true_quad_values);
73 for (
unsigned int qp=0; qp<true_quad_points.size(); qp++) {
74 double t = std::exp(true_quad_points[qp]);
75 mean_1d += t*true_quad_weights[qp];
76 sd_1d += t*t*true_quad_weights[qp];
79 const unsigned int dmin = 1;
80 const unsigned int dmax = 4;
81 const unsigned int pmin = 1;
82 const unsigned int pmax = 5;
85 for (
unsigned int d=dmin; d<=dmax; d++) {
88 double true_mean = std::pow(mean_1d,
static_cast<double>(d));
89 double true_sd = std::pow(sd_1d,
static_cast<double>(d)) -
91 true_sd = std::sqrt(true_sd);
92 std::cout.precision(12);
93 std::cout <<
"true mean = " << true_mean <<
"\t true std. dev. = "
94 << true_sd << std::endl;
96 Teuchos::Array< Teuchos::RCP<const Stokhos::OneDOrthogPolyBasis<int,double> > > bases(d);
99 for (
unsigned int p=pmin; p<=pmax; p++) {
102 for (
unsigned int i=0; i<d; i++)
104 Teuchos::RCP<const Stokhos::CompletePolynomialBasis<int,double> > basis =
108 int sz = basis->size();
110 for (
unsigned int i=0; i<d; i++) {
115 Teuchos::RCP<const Stokhos::Quadrature<int,double> > quad =
119 Teuchos::RCP<Stokhos::Sparse3Tensor<int,double> > Cijk =
120 basis->computeTripleProductTensor();
128 double mean = u.
mean();
131 std::cout.precision(4);
132 std::cout.setf(std::ios::scientific);
133 std::cout <<
"d = " << d <<
" p = " << p
136 << std::fabs(true_mean-mean) <<
"\tstd. dev. err = "
137 << std::fabs(true_sd-sd) << std::endl;
142 catch (std::exception& e) {
143 std::cout << e.what() << std::endl;
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.
value_type mean() const
Compute mean of expansion.
value_type standard_deviation() const
Compute standard deviation of expansion.
Orthogonal polynomial expansions based on numerical quadrature.
void exp(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
virtual void getQuadPoints(ordinal_type quad_order, Teuchos::Array< value_type > &points, Teuchos::Array< value_type > &weights, Teuchos::Array< Teuchos::Array< value_type > > &values) const
Compute quadrature points, weights, and values of basis polynomials at given set of points points.
Defines quadrature for a tensor product basis by tensor products of 1-D quadrature rules.
int main(int argc, char **argv)
Stokhos::LegendreBasis< int, double > basis_type