45template <
typename ordinal_type,
typename value_type>
51 const Teuchos::ParameterList& params) :
53 name(
"Monomial Proj Gram Schmidt PCE Basis")
55 this->
setup(max_p, pce, quad);
58template <
typename ordinal_type,
typename value_type>
64template <
typename ordinal_type,
typename value_type>
72template <
typename ordinal_type,
typename value_type>
78 const Teuchos::SerialDenseMatrix<ordinal_type,value_type>& A,
79 const Teuchos::SerialDenseMatrix<ordinal_type,value_type>& F,
80 const Teuchos::Array<value_type>& weights,
82 Teuchos::Array<ordinal_type>& num_terms_,
83 Teuchos::SerialDenseMatrix<ordinal_type,value_type>& Qp_,
84 Teuchos::SerialDenseMatrix<ordinal_type,value_type>& Q_)
88 CPBUtils::compute_terms(max_p, this->d, max_sz, terms_, num_terms_);
96 ordinal_type nqp = weights.size();
98 for (ordinal_type i=0; i<nqp; i++) {
99 for (ordinal_type
j=0;
j<max_sz;
j++) {
101 for (ordinal_type k=0; k<this->d; k++)
102 B(i,
j) *= std::pow(F(i,k), terms_[
j][k]);
107 SDM Bp(this->pce_sz, max_sz);
108 const Teuchos::Array<value_type>& basis_norms =
109 this->pce_basis->norm_squared();
110 for (ordinal_type i=0; i<this->pce_sz; i++) {
111 for (ordinal_type
j=0;
j<max_sz;
j++) {
113 for (ordinal_type k=0; k<nqp; k++)
114 Bp(i,
j) += weights[k]*B(k,
j)*A(k,i);
115 Bp(i,
j) /= basis_norms[i];
120 for (ordinal_type
j=0;
j<max_sz;
j++) {
121 value_type nrm = 0.0;
122 for (ordinal_type i=0; i<this->pce_sz; i++)
123 nrm += Bp(i,
j)*Bp(i,
j)*basis_norms[i];
124 nrm = std::sqrt(nrm);
125 for (ordinal_type i=0; i<this->pce_sz; i++)
132 Teuchos::Array<value_type> w(this->pce_sz, 1.0);
134 Teuchos::Array<ordinal_type> piv(max_sz);
135 for (
int i=0; i<this->d+1; i++)
138 ordinal_type sz_ = SOF::createOrthogonalBasis(
139 this->orthogonalization_method, threshold, this->verbose, Bp, w,
143 Q_.reshape(nqp, sz_);
145 Q_.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1.0, A, Qp_, 0.0);
146 TEUCHOS_ASSERT(ret == 0);
UnitTestSetup< Kokkos::Cuda > setup
Generate a basis from a given set of PCE expansions that is orthogonal with respect to the product me...
virtual const std::string & getName() const
Return string name of basis.
virtual ordinal_type buildReducedBasis(ordinal_type max_p, value_type threshold, const Teuchos::SerialDenseMatrix< ordinal_type, value_type > &A, const Teuchos::SerialDenseMatrix< ordinal_type, value_type > &F, const Teuchos::Array< value_type > &weights, Teuchos::Array< Stokhos::MultiIndex< ordinal_type > > &terms_, Teuchos::Array< ordinal_type > &num_terms_, Teuchos::SerialDenseMatrix< ordinal_type, value_type > &Qp_, Teuchos::SerialDenseMatrix< ordinal_type, value_type > &Q_)
Build the reduced basis, parameterized by total order max_p.
MonomialProjGramSchmidtPCEBasis(ordinal_type p, const Teuchos::Array< Stokhos::OrthogPolyApprox< ordinal_type, value_type > > &pce, const Teuchos::RCP< const Stokhos::Quadrature< ordinal_type, value_type > > &quad, const Teuchos::ParameterList ¶ms=Teuchos::ParameterList())
Constructor.
Teuchos::SerialDenseMatrix< ordinal_type, value_type > SDM
virtual ~MonomialProjGramSchmidtPCEBasis()
Destructor.
A multidimensional index.
Class to store coefficients of a projection onto an orthogonal polynomial basis.
Encapsulate various orthogonalization (ie QR) methods.
Abstract base class for quadrature methods.