44template <
typename ordinal_type,
typename value_type>
45Stokhos::PecosOneDOrthogPolyBasis<ordinal_type, value_type>::
46PecosOneDOrthogPolyBasis(
47 const Teuchos::RCP<Pecos::OrthogonalPolynomial>& pecosPoly_,
48 const std::string& name_, ordinal_type p_) :
49 pecosPoly(pecosPoly_),
52 sparse_grid_growth_rule(webbur::level_to_order_linear_wn),
55 for (ordinal_type i=0; i<=p; i++)
56 norms[i] = pecosPoly->norm_squared(i);
59template <
typename ordinal_type,
typename value_type>
60Stokhos::PecosOneDOrthogPolyBasis<ordinal_type, value_type>::
61PecosOneDOrthogPolyBasis(
62 ordinal_type p_,
const PecosOneDOrthogPolyBasis& basis) :
63 pecosPoly(basis.pecosPoly),
66 sparse_grid_growth_rule(basis.sparse_grid_growth_rule),
69 for (ordinal_type i=0; i<=p; i++)
70 norms[i] = pecosPoly->norm_squared(i);
73template <
typename ordinal_type,
typename value_type>
74Stokhos::PecosOneDOrthogPolyBasis<ordinal_type, value_type>::
75~PecosOneDOrthogPolyBasis()
79template <
typename ordinal_type,
typename value_type>
81Stokhos::PecosOneDOrthogPolyBasis<ordinal_type, value_type>::
87template <
typename ordinal_type,
typename value_type>
89Stokhos::PecosOneDOrthogPolyBasis<ordinal_type, value_type>::
95template <
typename ordinal_type,
typename value_type>
96const Teuchos::Array<value_type>&
97Stokhos::PecosOneDOrthogPolyBasis<ordinal_type, value_type>::
103template <
typename ordinal_type,
typename value_type>
105Stokhos::PecosOneDOrthogPolyBasis<ordinal_type, value_type>::
106norm_squared(ordinal_type i)
const
111template <
typename ordinal_type,
typename value_type>
112Teuchos::RCP< Stokhos::Dense3Tensor<ordinal_type, value_type> >
113Stokhos::PecosOneDOrthogPolyBasis<ordinal_type, value_type>::
114computeTripleProductTensor()
const
118 Teuchos::RCP< Stokhos::Dense3Tensor<ordinal_type, value_type> > Cijk =
119 Teuchos::rcp(
new Dense3Tensor<ordinal_type, value_type>(sz));
120 Teuchos::Array<value_type> points, weights;
121 Teuchos::Array< Teuchos::Array<value_type> > values;
122 getQuadPoints(3*p, points, weights, values);
124 for (ordinal_type i=0; i<sz; i++) {
125 for (ordinal_type
j=0;
j<sz;
j++) {
126 for (ordinal_type k=0; k<sz; k++) {
128 for (ordinal_type l=0; l<static_cast<ordinal_type>(points.size());
131 weights[l]*(values[l][i])*(values[l][
j])*(values[l][k]);
133 (*Cijk)(i,
j,k) = triple_product;
141template <
typename ordinal_type,
typename value_type>
142Teuchos::RCP< Stokhos::Sparse3Tensor<ordinal_type, value_type> >
143Stokhos::PecosOneDOrthogPolyBasis<ordinal_type, value_type>::
144computeSparseTripleProductTensor(ordinal_type order)
const
149 Teuchos::RCP< Stokhos::Sparse3Tensor<ordinal_type, value_type> > Cijk =
150 Teuchos::rcp(
new Sparse3Tensor<ordinal_type, value_type>());
151 Teuchos::Array<value_type> points, weights;
152 Teuchos::Array< Teuchos::Array<value_type> > values;
153 getQuadPoints(3*p, points, weights, values);
155 for (ordinal_type i=0; i<sz; i++) {
156 for (ordinal_type
j=0;
j<sz;
j++) {
157 for (ordinal_type k=0; k<order; k++) {
159 for (ordinal_type l=0; l<static_cast<ordinal_type>(points.size());
162 weights[l]*(values[l][i])*(values[l][
j])*(values[l][k]);
164 if (std::abs(triple_product/norms[i]) > sparse_tol)
165 Cijk->add_term(i,
j,k,triple_product);
169 Cijk->fillComplete();
174template <
typename ordinal_type,
typename value_type>
175Teuchos::RCP< Teuchos::SerialDenseMatrix<ordinal_type, value_type> >
176Stokhos::PecosOneDOrthogPolyBasis<ordinal_type, value_type>::
177computeDerivDoubleProductTensor()
const
180 Teuchos::Array<value_type> points, weights;
181 Teuchos::Array< Teuchos::Array<value_type> > values, derivs;
182 getQuadPoints(2*p, points, weights, values);
186 for (ordinal_type i=0; i<nqp; i++) {
187 derivs[i].resize(sz);
188 evaluateBasesAndDerivatives(points[i], values[i], derivs[i]);
190 Teuchos::RCP< Teuchos::SerialDenseMatrix<ordinal_type, value_type> > Bij =
191 Teuchos::rcp(
new Teuchos::SerialDenseMatrix<ordinal_type, value_type>(sz,sz));
192 for (ordinal_type i=0; i<sz; i++) {
193 for (ordinal_type
j=0;
j<sz;
j++) {
195 for (
int qp=0; qp<nqp; qp++)
196 b += weights[qp]*derivs[qp][i]*values[qp][
j];
204template <
typename ordinal_type,
typename value_type>
206Stokhos::PecosOneDOrthogPolyBasis<ordinal_type,value_type>::
207evaluateBases(
const value_type& x, Teuchos::Array<value_type>& basis_pts)
const
209 for (ordinal_type i=0; i<=p; i++)
210 basis_pts[i] = pecosPoly->type1_value(x, i);
213template <
typename ordinal_type,
typename value_type>
215Stokhos::PecosOneDOrthogPolyBasis<ordinal_type,value_type>::
216evaluateBasesAndDerivatives(
const value_type& x,
217 Teuchos::Array<value_type>& vals,
218 Teuchos::Array<value_type>& derivs)
const
220 for (ordinal_type i=0; i<=p; i++) {
221 vals[i] = pecosPoly->type1_value(x, i);
222 derivs[i] = pecosPoly->type1_gradient(x, i);
226template <
typename ordinal_type,
typename value_type>
228Stokhos::PecosOneDOrthogPolyBasis<ordinal_type,value_type>::
229evaluate(
const value_type& x, ordinal_type k)
const
231 return pecosPoly->type1_value(x, k);
234template <
typename ordinal_type,
typename value_type>
236Stokhos::PecosOneDOrthogPolyBasis<ordinal_type, value_type>::
237print(std::ostream& os)
const
239 os <<
"Pecos " << name <<
" basis of order " << p <<
"." << std::endl;
240 os <<
"Basis polynomial norms (squared):\n\t";
241 for (ordinal_type i=0; i<=p; i++)
242 os << norms[i] <<
" ";
246template <
typename ordinal_type,
typename value_type>
248Stokhos::PecosOneDOrthogPolyBasis<ordinal_type, value_type>::
254template <
typename ordinal_type,
typename value_type>
256Stokhos::PecosOneDOrthogPolyBasis<ordinal_type,value_type>::
257getQuadPoints(ordinal_type quad_order,
258 Teuchos::Array<value_type>& quad_points,
259 Teuchos::Array<value_type>& quad_weights,
260 Teuchos::Array< Teuchos::Array<value_type> >& quad_values)
const
264 static_cast<ordinal_type>(std::ceil((quad_order+1)/2.0));
265 const Pecos::RealArray& gp = pecosPoly->collocation_points(num_points);
266 const Pecos::RealArray& gw = pecosPoly->type1_collocation_weights(num_points);
267 quad_points.resize(num_points);
268 quad_weights.resize(num_points);
269 for (ordinal_type i=0; i<num_points; i++) {
270 quad_points[i] = gp[i];
271 quad_weights[i] = gw[i];
275 quad_values.resize(num_points);
276 for (ordinal_type i=0; i<num_points; i++) {
277 quad_values[i].resize(p+1);
278 evaluateBases(quad_points[i], quad_values[i]);
282template <
typename ordinal_type,
typename value_type>
284Stokhos::PecosOneDOrthogPolyBasis<ordinal_type,value_type>::
285quadDegreeOfExactness(ordinal_type n)
const
290template <
typename ordinal_type,
typename value_type>
291Teuchos::RCP<Stokhos::OneDOrthogPolyBasis<ordinal_type,value_type> >
292Stokhos::PecosOneDOrthogPolyBasis<ordinal_type,value_type>::
293cloneWithOrder(ordinal_type p)
const
296 Teuchos::rcp(
new Stokhos::PecosOneDOrthogPolyBasis<ordinal_type,value_type>(p,*
this));
299template <
typename ordinal_type,
typename value_type>
301Stokhos::PecosOneDOrthogPolyBasis<ordinal_type,value_type>::
302coefficientGrowth(ordinal_type n)
const
307template <
typename ordinal_type,
typename value_type>
309Stokhos::PecosOneDOrthogPolyBasis<ordinal_type,value_type>::
310pointGrowth(ordinal_type n)
const