44#ifndef STOKHOS_SPD_DENSE_DIRECT_DIVISION_EXPANSION_STRATEGY_HPP
45#define STOKHOS_SPD_DENSE_DIRECT_DIVISION_EXPANSION_STRATEGY_HPP
51#include "Teuchos_TimeMonitor.hpp"
52#include "Teuchos_RCP.hpp"
53#include "Teuchos_SerialSymDenseMatrix.hpp"
54#include "Teuchos_SerialSpdDenseSolver.hpp"
62 template <
typename ordinal_type,
typename value_type,
typename node_type>
78 const value_type& alpha,
81 const value_type& beta);
96 Teuchos::RCP<const Stokhos::OrthogPolyBasis<ordinal_type, value_type> >
basis;
102 Teuchos::RCP<const Cijk_type>
Cijk;
105 Teuchos::RCP< Teuchos::SerialSymDenseMatrix<ordinal_type,value_type> >
A;
107 Teuchos::RCP< Teuchos::SerialDenseMatrix<ordinal_type,value_type> >
X,
B;
110 Teuchos::SerialSpdDenseSolver<ordinal_type,value_type>
solver;
116template <
typename ordinal_type,
typename value_type,
typename node_type>
125 ordinal_type sz =
basis->size();
126 A = Teuchos::rcp(
new Teuchos::SerialSymDenseMatrix<ordinal_type,value_type>(
128 B = Teuchos::rcp(
new Teuchos::SerialDenseMatrix<ordinal_type,value_type>(
130 X = Teuchos::rcp(
new Teuchos::SerialDenseMatrix<ordinal_type,value_type>(
134template <
typename ordinal_type,
typename value_type,
typename node_type>
138 const value_type& alpha,
141 const value_type& beta)
143#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
144 TEUCHOS_FUNC_TIME_MONITOR(
"Stokhos::SPDDenseDirectDivisionStrategy::divide()");
147 ordinal_type sz = basis->size();
148 ordinal_type pa = a.
size();
149 ordinal_type pb = b.
size();
158 const value_type* ca = a.
coeff();
159 const value_type* cb = b.
coeff();
160 value_type* cc = c.
coeff();
167 if (pb < Cijk->num_k())
168 k_end = Cijk->find_k(pb);
174 j_it != Cijk->j_end(k_it); ++j_it) {
177 i_it != Cijk->i_end(j_it); ++i_it) {
180 (*A)(i,
j) += cijk*cb[k];
188 for (ordinal_type i=0; i<pa; i++)
189 (*B)(i,0) = ca[i]*basis->norm_squared(i);
193 solver.setVectors(X, B);
194 if (solver.shouldEquilibrate()) {
195 solver.factorWithEquilibration(
true);
196 solver.equilibrateMatrix();
207 for (ordinal_type i=0; i<pc; i++)
208 cc[i] = alpha*(*X)(i,0) + beta*cc[i];
211 for (ordinal_type i=0; i<pc; i++)
212 cc[i] = alpha*ca[i]/cb[0] + beta*cc[i];
Strategy interface for computing PCE of a/b.
Class to store coefficients of a projection onto an orthogonal polynomial basis.
void resize(ordinal_type sz)
Resize coefficient array (coefficients are preserved)
pointer coeff()
Return coefficient array.
ordinal_type size() const
Return size.
Abstract base class for multivariate orthogonal polynomials.
Strategy interface for computing PCE of a/b using only b[0].
Teuchos::RCP< const Stokhos::OrthogPolyBasis< ordinal_type, value_type > > basis
Basis.
Teuchos::SerialSpdDenseSolver< ordinal_type, value_type > solver
Serial dense solver.
Teuchos::RCP< const Cijk_type > Cijk
Triple product.
virtual void divide(Stokhos::OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const value_type &alpha, const Stokhos::OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const Stokhos::OrthogPolyApprox< ordinal_type, value_type, node_type > &b, const value_type &beta)
Teuchos::RCP< Teuchos::SerialSymDenseMatrix< ordinal_type, value_type > > A
Dense matrices for linear system.
Stokhos::Sparse3Tensor< ordinal_type, value_type > Cijk_type
Short-hand for Cijk.
Teuchos::RCP< Teuchos::SerialDenseMatrix< ordinal_type, value_type > > X
SPDDenseDirectDivisionExpansionStrategy(const SPDDenseDirectDivisionExpansionStrategy &)
Teuchos::RCP< Teuchos::SerialDenseMatrix< ordinal_type, value_type > > B
SPDDenseDirectDivisionExpansionStrategy(const Teuchos::RCP< const Stokhos::OrthogPolyBasis< ordinal_type, value_type > > &basis_, const Teuchos::RCP< const Stokhos::Sparse3Tensor< ordinal_type, value_type > > &Cijk_)
Constructor.
virtual ~SPDDenseDirectDivisionExpansionStrategy()
Destructor.
SPDDenseDirectDivisionExpansionStrategy & operator=(const SPDDenseDirectDivisionExpansionStrategy &b)
Data structure storing a sparse 3-tensor C(i,j,k) in a a compressed format.
kji_sparse_array::const_iterator k_iterator
Iterator for looping over k entries.
j_sparse_array::const_iterator kji_iterator
Iterator for looping over i entries given k and j.
Top-level namespace for Stokhos classes and functions.
Bi-directional iterator for traversing a sparse array.