MUQ  0.4.3
AdaptiveSmolyakPCE.cpp
Go to the documentation of this file.
2 
3 
4 using namespace muq::Approximation;
5 using namespace muq::Utilities;
6 
7 
8 
9 
10 
11 AdaptiveSmolyakPCE::AdaptiveSmolyakPCE(std::shared_ptr<muq::Modeling::ModPiece> const& modelIn,
12  std::vector<std::shared_ptr<Quadrature>> const& quad1dIn,
13  std::vector<std::shared_ptr<IndexedScalarBasis>> const& polys1dIn)
14  : SmolyakEstimator<std::shared_ptr<PolynomialChaosExpansion>>(modelIn),
15  tensFactory(quad1dIn, polys1dIn)
16 {}
17 
18 
19 
20 std::vector<Eigen::VectorXd> AdaptiveSmolyakPCE::OneTermPoints(std::shared_ptr<MultiIndex> const& multi)
21 {
22  return tensFactory.QuadPts(multi);
23 }
24 
25 std::shared_ptr<PolynomialChaosExpansion> AdaptiveSmolyakPCE::ComputeOneTerm(std::shared_ptr<MultiIndex> const& multi,
26  std::vector<std::reference_wrapper<const Eigen::VectorXd>> const& modEvals)
27 {
28  return tensFactory.Compute(modEvals,multi);
29 }
30 
31 std::shared_ptr<PolynomialChaosExpansion> AdaptiveSmolyakPCE::AddEstimates(double w1,
32  std::shared_ptr<PolynomialChaosExpansion> const& part1,
33  double w2,
34  std::shared_ptr<PolynomialChaosExpansion> const& part2) const
35 {
36  Eigen::VectorXd wts(2);
37  wts << w1, w2;
38  return PolynomialChaosExpansion::ComputeWeightedSum({part1,part2},wts);
39 }
40 
41 
42 double AdaptiveSmolyakPCE::ComputeMagnitude(std::shared_ptr<PolynomialChaosExpansion> const& estimate) const
43 {
44  return estimate->Magnitude().norm();
45 }
virtual std::vector< Eigen::VectorXd > OneTermPoints(std::shared_ptr< muq::Utilities::MultiIndex > const &multi) override
virtual double ComputeMagnitude(std::shared_ptr< PolynomialChaosExpansion > const &estimate) const override
AdaptiveSmolyakPCE(std::shared_ptr< muq::Modeling::ModPiece > const &modelIn, std::vector< std::shared_ptr< Quadrature >> const &quad1dIn, std::vector< std::shared_ptr< IndexedScalarBasis >> const &polys1dIn)
virtual std::shared_ptr< PolynomialChaosExpansion > AddEstimates(double w1, std::shared_ptr< PolynomialChaosExpansion > const &part1, double w2, std::shared_ptr< PolynomialChaosExpansion > const &part2) const override
virtual std::shared_ptr< PolynomialChaosExpansion > ComputeOneTerm(std::shared_ptr< muq::Utilities::MultiIndex > const &multi, std::vector< std::reference_wrapper< const Eigen::VectorXd >> const &modEvals) override
std::vector< Eigen::VectorXd > const & QuadPts() const
Definition: PCEFactory.h:107
std::shared_ptr< PolynomialChaosExpansion > Compute(std::shared_ptr< muq::Modeling::ModPiece > const &model)
Definition: PCEFactory.cpp:65
A class for representing and using expansions of orthogonal multivariate polynomials.
static std::shared_ptr< PolynomialChaosExpansion > ComputeWeightedSum(std::vector< std::shared_ptr< PolynomialChaosExpansion >> expansions, Eigen::VectorXd const &weights)