43#include "Teuchos_TimeMonitor.hpp"
45template <
typename ordinal_type,
typename value_type>
49#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
50 TEUCHOS_FUNC_TIME_MONITOR(
"Stokhos::TensorProductQuadrature -- Quad Grid Generation");
52 ordinal_type d = product_basis->dimension();
53 ordinal_type sz = product_basis->size();
55 Teuchos::Array< Teuchos::RCP<const OneDOrthogPolyBasis<ordinal_type,value_type> > > coordinate_bases = product_basis->getCoordinateBases();
58 Teuchos::Array< Teuchos::Array<value_type> > gp(d);
59 Teuchos::Array< Teuchos::Array<value_type> > gw(d);
60 Teuchos::Array< Teuchos::Array< Teuchos::Array<value_type> > > gv(d);
61 Teuchos::Array<ordinal_type> n(d);
62 ordinal_type ntot = 1;
63 for (ordinal_type i=0; i<d; i++) {
64 coordinate_bases[i]->getQuadPoints(2*(coordinate_bases[i]->order()),
69 quad_points.resize(ntot);
70 quad_weights.resize(ntot);
71 quad_values.resize(ntot);
72 Teuchos::Array<ordinal_type> index(d);
73 for (ordinal_type i=0; i<d; i++)
77 quad_points[cnt].resize(d);
78 quad_weights[cnt] = value_type(1.0);
79 quad_values[cnt].resize(sz);
80 for (ordinal_type
j=0;
j<d;
j++) {
81 quad_points[cnt][
j] = gp[
j][index[
j]];
82 quad_weights[cnt] *= gw[
j][index[
j]];
84 for (ordinal_type k=0; k<sz; k++) {
85 quad_values[cnt][k] = value_type(1.0);
87 for (ordinal_type
j=0;
j<d;
j++)
88 quad_values[cnt][k] *= gv[
j][index[
j]][term[
j]];
92 while (i < d-1 && index[i] == n[i]) {
103template <
typename ordinal_type,
typename value_type>
107#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
108 TEUCHOS_FUNC_TIME_MONITOR(
"Stokhos::TensorProductQuadrature -- Quad Grid Generation");
110 ordinal_type d = product_basis->dimension();
111 ordinal_type sz = product_basis->size();
113 Teuchos::Array< Teuchos::RCP<const OneDOrthogPolyBasis<ordinal_type,value_type> > > coordinate_bases = product_basis->getCoordinateBases();
116 Teuchos::Array< Teuchos::Array<value_type> > gp(d);
117 Teuchos::Array< Teuchos::Array<value_type> > gw(d);
118 Teuchos::Array< Teuchos::Array< Teuchos::Array<value_type> > > gv(d);
119 Teuchos::Array<ordinal_type> n(d);
120 ordinal_type ntot = 1;
121 for (ordinal_type i=0; i<d; i++) {
122 coordinate_bases[i]->getQuadPoints(quad_order,
123 gp[i], gw[i], gv[i]);
127 quad_points.resize(ntot);
128 quad_weights.resize(ntot);
129 quad_values.resize(ntot);
130 Teuchos::Array<ordinal_type> index(d);
131 for (ordinal_type i=0; i<d; i++)
133 ordinal_type cnt = 0;
135 quad_points[cnt].resize(d);
136 quad_weights[cnt] = value_type(1.0);
137 quad_values[cnt].resize(sz);
138 for (ordinal_type
j=0;
j<d;
j++) {
139 quad_points[cnt][
j] = gp[
j][index[
j]];
140 quad_weights[cnt] *= gw[
j][index[
j]];
142 for (ordinal_type k=0; k<sz; k++) {
143 quad_values[cnt][k] = value_type(1.0);
145 for (ordinal_type
j=0;
j<d;
j++)
146 quad_values[cnt][k] *= gv[
j][index[
j]][term[
j]];
150 while (i < d-1 && index[i] == n[i]) {
161template <
typename ordinal_type,
typename value_type>
162const Teuchos::Array< Teuchos::Array<value_type> >&
169template <
typename ordinal_type,
typename value_type>
170const Teuchos::Array<value_type>&
177template <
typename ordinal_type,
typename value_type>
178const Teuchos::Array< Teuchos::Array<value_type> >&
185template <
typename ordinal_type,
typename value_type>
188print(std::ostream& os)
const
190 ordinal_type nqp = quad_weights.size();
191 os <<
"Tensor Product Quadrature with " << nqp <<
" points:"
192 << std::endl <<
"Weight : Points" << std::endl;
193 for (ordinal_type i=0; i<nqp; i++) {
194 os << i <<
": " << quad_weights[i] <<
" : ";
195 for (ordinal_type
j=0; j<static_cast<ordinal_type>(quad_points[i].size());
197 os << quad_points[i][
j] <<
" ";
200 os <<
"Basis values at quadrature points:" << std::endl;
201 for (ordinal_type i=0; i<nqp; i++) {
202 os << i <<
" " <<
": ";
203 for (ordinal_type
j=0; j<static_cast<ordinal_type>(quad_values[i].size());
205 os << quad_values[i][
j] <<
" ";
A multidimensional index.
Abstract base class for multivariate orthogonal polynomials generated from tensor products of univari...
virtual const Teuchos::Array< Teuchos::Array< value_type > > & getQuadPoints() const
Get quadrature points.
virtual const Teuchos::Array< value_type > & getQuadWeights() const
Get quadrature weights.
virtual std::ostream & print(std::ostream &os) const
Print quadrature data.
TensorProductQuadrature(const Teuchos::RCP< const ProductBasis< ordinal_type, value_type > > &product_basis)
Constructor.
virtual const Teuchos::Array< Teuchos::Array< value_type > > & getBasisAtQuadPoints() const
Get values of basis at quadrature points.