44 #include "Teuchos_BLAS.hpp" 46 template <
typename ordinal_type,
typename value_type>
49 const Teuchos::RCP<
const OrthogPolyBasis<ordinal_type,value_type> >& basis_,
50 const Teuchos::Array< Teuchos::Array<value_type> >& points,
51 const Teuchos::Array<value_type>& weights_,
53 name(
"Gram Schmidt Basis"),
56 sparse_tol(sparse_tol_),
58 d(basis->dimension()),
66 Teuchos::Array< Teuchos::Array<value_type> > values(nqp);
69 basis->evaluateBases(points[k], values[k]);
73 Teuchos::SerialDenseMatrix<ordinal_type, value_type> inner_product(sz,sz);
74 inner_product.putScalar(0.0);
79 t += weights[k]*values[k][i]*values[k][
j];
80 inner_product(i,
j) = t;
98 t += gs_mat(
j,k)*inner_product(i,k);
103 gs_mat(i,k) -= t*gs_mat(
j,k);
110 nrm += gs_mat(i,
j)*gs_mat(i,k)*inner_product(
j,k);
112 nrm += gs_mat(i,
j)*gs_mat(i,k)*inner_product(k,
j);
118 basis_values.resize(nqp);
120 basis_values[k].resize(sz);
124 t += gs_mat(i,
j)*values[k][
j];
125 basis_values[k][i] = t;
130 template <
typename ordinal_type,
typename value_type>
136 template <
typename ordinal_type,
typename value_type>
144 template <
typename ordinal_type,
typename value_type>
152 template <
typename ordinal_type,
typename value_type>
160 template <
typename ordinal_type,
typename value_type>
161 const Teuchos::Array<value_type>&
168 template <
typename ordinal_type,
typename value_type>
176 template <
typename ordinal_type,
typename value_type>
177 Teuchos::RCP< Stokhos::Sparse3Tensor<ordinal_type, value_type> >
181 Teuchos::RCP< Stokhos::Sparse3Tensor<ordinal_type, value_type> > Cijk =
182 Teuchos::rcp(
new Sparse3Tensor<ordinal_type, value_type>);
190 weights[l]*basis_values[l][i]*basis_values[l][
j]*basis_values[l][k];
192 Cijk->add_term(i,
j,k,t);
197 Cijk->fillComplete();
202 template <
typename ordinal_type,
typename value_type>
203 Teuchos::RCP< Stokhos::Sparse3Tensor<ordinal_type, value_type> >
207 Teuchos::RCP< Stokhos::Sparse3Tensor<ordinal_type, value_type> > Cijk =
208 Teuchos::rcp(
new Sparse3Tensor<ordinal_type, value_type>);
216 weights[l]*basis_values[l][i]*basis_values[l][
j]*basis_values[l][k];
218 Cijk->add_term(i,
j,k,t);
223 Cijk->fillComplete();
228 template <
typename ordinal_type,
typename value_type>
235 z += gs_mat(i,
j)*basis->evaluateZero(
j);
240 template <
typename ordinal_type,
typename value_type>
243 evaluateBases(
const Teuchos::ArrayView<const value_type>& point,
244 Teuchos::Array<value_type>& basis_vals)
const 246 basis->evaluateBases(point, basis_vals_tmp);
250 t += gs_mat(i,
j)*basis_vals_tmp[
j];
255 template <
typename ordinal_type,
typename value_type>
258 print(std::ostream& os)
const 260 os <<
"Gram-Schmidt basis of order " << p <<
", dimension " << d
261 <<
", and size " << sz <<
". Matrix coefficients:\n";
262 os << gs_mat << std::endl;
263 os <<
"Basis vector norms (squared):\n\t";
265 os << norms[i] <<
" ";
267 os <<
"Underlying basis:\n";
271 template <
typename ordinal_type,
typename value_type>
279 template <
typename ordinal_type,
typename value_type>
284 Teuchos::BLAS<ordinal_type, value_type> blas;
287 blas.TRSM(Teuchos::LEFT_SIDE, Teuchos::LOWER_TRI, Teuchos::TRANS,
288 Teuchos::UNIT_DIAG, sz, 1, 1.0, gs_mat.values(), sz, out, sz);
virtual Teuchos::RCP< Stokhos::Sparse3Tensor< ordinal_type, value_type > > computeTripleProductTensor() const =0
Compute triple product tensor.
virtual ordinal_type order() const =0
Return order of basis.
virtual void evaluateBases(const Teuchos::ArrayView< const value_type > &point, Teuchos::Array< value_type > &basis_vals) const =0
Evaluate basis polynomials at given point point.
virtual ordinal_type dimension() const =0
Return dimension of basis.
virtual value_type evaluateZero(ordinal_type i) const =0
Evaluate basis polynomial i at zero.
Transforms a non-orthogonal multivariate basis to an orthogonal one using the Gram-Schmit procedure...
KOKKOS_INLINE_FUNCTION PCE< Storage > abs(const PCE< Storage > &a)
virtual const Teuchos::Array< value_type > & norm_squared() const =0
Return array storing norm-squared of each basis polynomial.
virtual void print(std::ostream &os) const =0
Print basis to stream os.
virtual Teuchos::RCP< Stokhos::Sparse3Tensor< ordinal_type, value_type > > computeLinearTripleProductTensor() const =0
Compute linear triple product tensor where k = 0,1.
virtual ordinal_type size() const =0
Return total size of basis.
virtual const std::string & getName() const =0
Return string name of basis.