A class for representing and using expansions of orthogonal multivariate polynomials. More...
#include <PolynomialChaosExpansion.h>
A class for representing and using expansions of orthogonal multivariate polynomials.
A particular polynomial chaos expansion for a function from R^n->R^m. This class uses some the MultiIndexSet class to define PCE terms that are in the expansion. For each PCE term and output, there is a coefficient. The PCE is built from 1D polynomials specified so each dimension may have independent polynomial choices.
Definition at line 28 of file PolynomialChaosExpansion.h.
Public Member Functions | |
PolynomialChaosExpansion (std::shared_ptr< OrthogonalPolynomial > const &basisCompsIn, std::shared_ptr< muq::Utilities::MultiIndexSet > const &multisIn, Eigen::MatrixXd const &coeffsIn) | |
PolynomialChaosExpansion (std::shared_ptr< OrthogonalPolynomial > const &basisCompsIn, std::shared_ptr< muq::Utilities::MultiIndexSet > const &multisIn, unsigned int outputDim) | |
PolynomialChaosExpansion (std::vector< std::shared_ptr< IndexedScalarBasis >> const &basisCompsIn, std::shared_ptr< muq::Utilities::MultiIndexSet > const &multisIn, Eigen::MatrixXd const &coeffsIn) | |
PolynomialChaosExpansion (std::vector< std::shared_ptr< IndexedScalarBasis >> const &basisCompsIn, std::shared_ptr< muq::Utilities::MultiIndexSet > const &multisIn, unsigned int outputDim) | |
virtual | ~PolynomialChaosExpansion ()=default |
Eigen::VectorXd | Variance () const |
compute the variance of the current expansion More... | |
Eigen::MatrixXd | Covariance () const |
Eigen::VectorXd | Mean () const |
Eigen::VectorXd | Magnitude () const |
Compute the L2 norm of each output. More... | |
Eigen::VectorXd | TotalSensitivity (unsigned int targetDim) const |
Compute the Sobol total sensitivity index for the input dimension, for each output dimension. More... | |
Eigen::MatrixXd | TotalSensitivity () const |
Compute all Sobol total sensitivities. Rows are outputs, each column is an input. More... | |
Eigen::VectorXd | SobolSensitivity (unsigned int targetDim) const |
Compute the main sensitivity index for the input dimension, for each output dimension. More... | |
Eigen::VectorXd | SobolSensitivity (std::vector< unsigned int > const &targetDims) const |
Eigen::MatrixXd | MainSensitivity () const |
virtual void | ToHDF5 (muq::Utilities::H5Object &group) const override |
Write PCE to HDF5 file. More... | |
Public Member Functions inherited from muq::Approximation::BasisExpansion | |
BasisExpansion (std::vector< std::shared_ptr< IndexedScalarBasis >> const &basisCompsIn, bool coeffInput=false) | |
BasisExpansion (std::vector< std::shared_ptr< IndexedScalarBasis >> const &basisCompsIn, std::shared_ptr< muq::Utilities::MultiIndexSet > multisIn, bool coeffInput=false) | |
BasisExpansion (std::vector< std::shared_ptr< IndexedScalarBasis >> const &basisCompsIn, std::shared_ptr< muq::Utilities::MultiIndexSet > multisIn, Eigen::MatrixXd const &coeffsIn, bool coeffInput=false) | |
virtual | ~BasisExpansion ()=default |
virtual unsigned | NumTerms () const |
Eigen::MatrixXd | BuildVandermonde (Eigen::MatrixXd const &evalPts) const |
Eigen::MatrixXd | BuildDerivMatrix (Eigen::MatrixXd const &evalPts, int wrtDim) const |
Eigen::MatrixXd | SecondDerivative (unsigned outputDim, unsigned wrtDim1, unsigned wrtDim2, Eigen::VectorXd const &evalPt, Eigen::MatrixXd const &coeffs) |
Eigen::MatrixXd | SecondDerivative (unsigned outputDim, unsigned wrtDim1, unsigned wrtDim2, Eigen::VectorXd const &evalPt) |
Eigen::MatrixXd | GetCoeffs () const |
void | SetCoeffs (Eigen::MatrixXd const &allCoeffs) |
const std::shared_ptr< muq::Utilities::MultiIndexSet > | Multis () const |
virtual void | ToHDF5 (std::string filename, std::string groupName="/") const |
Saves the expansion to group in an HDF5 file. More... | |
Public Member Functions inherited from muq::Modeling::ModPiece | |
ModPiece (std::vector< int > const &inputSizes, std::vector< int > const &outputSizes) | |
ModPiece (Eigen::VectorXi const &inputSizes, Eigen::VectorXi const &outputSizes) | |
virtual | ~ModPiece ()=default |
virtual double | GetRunTime (const std::string &method="Evaluate") const override |
Get the average run time for one of the implemented methods. More... | |
virtual void | ResetCallTime () override |
Resets the number of call and times. More... | |
virtual unsigned long int | GetNumCalls (const std::string &method="Evaluate") const override |
get the number of times one of the implemented methods has been called. More... | |
virtual std::vector< Eigen::VectorXd > const & | Evaluate (std::vector< Eigen::VectorXd > const &input) |
Evaluate the ModPiece. More... | |
virtual std::vector< Eigen::VectorXd > const & | Evaluate (ref_vector< Eigen::VectorXd > const &input) |
VARIADIC_TO_REFVECTOR (Evaluate, Eigen::VectorXd, std::vector< Eigen::VectorXd > const &) | |
virtual Eigen::VectorXd const & | Gradient (unsigned int const outputDimWrt, unsigned int const inputDimWrt, std::vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sensitivity) |
Compute the Gradient \(J^Tv\). More... | |
virtual Eigen::VectorXd const & | Gradient (unsigned int const outputDimWrt, unsigned int const inputDimWrt, ref_vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sensitivity) |
Eigen::VectorXd const & | Gradient (unsigned int outWrt, unsigned int inWrt, Eigen::VectorXd const &last, Eigen::VectorXd const &sens) |
template<typename... Args> | |
Eigen::VectorXd const & | Gradient (unsigned int wrtOut, unsigned int wrtIn, Args const &... args) |
Eigen::VectorXd const & | ApplyHessian (unsigned int outWrt, unsigned int inWrt1, unsigned int inWrt2, Eigen::VectorXd const &last, Eigen::VectorXd const &sens, Eigen::VectorXd const &vec) |
template<typename... Args> | |
Eigen::VectorXd const & | ApplyHessian (unsigned int wrtOut, unsigned int wrtIn1, unsigned int wrtIn2, Args const &... args) |
virtual Eigen::MatrixXd const & | Jacobian (unsigned int const outputDimWrt, unsigned int const inputDimWrt, std::vector< Eigen::VectorXd > const &input) |
Compute the Jacobian of this ModPiece. More... | |
virtual Eigen::MatrixXd const & | Jacobian (unsigned int const outputDimWrt, unsigned int const inputDimWrt, ref_vector< Eigen::VectorXd > const &input) |
template<typename... Args> | |
Eigen::MatrixXd const & | Jacobian (unsigned int outWrt, unsigned int inWrt, Args const &... args) |
template<typename... Args> | |
Eigen::MatrixXd | JacobianByFD (unsigned int outWrt, unsigned int inWrt, Args const &... args) |
template<typename... Args> | |
Eigen::MatrixXd | ApplyJacobianByFD (unsigned int outWrt, unsigned int inWrt, Args const &... args) |
virtual Eigen::VectorXd const & | ApplyJacobian (unsigned int const outputDimWrt, unsigned int const inputDimWrt, std::vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &vec) |
Apply the Jacobian of this ModPiece to a vector. More... | |
virtual Eigen::VectorXd const & | ApplyJacobian (unsigned int const outputDimWrt, unsigned int const inputDimWrt, ref_vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &vec) |
Eigen::VectorXd const & | ApplyJacobian (unsigned int outWrt, unsigned int inWrt, Eigen::VectorXd const &last, Eigen::VectorXd const &sens) |
template<typename... Args> | |
Eigen::VectorXd const & | ApplyJacobian (unsigned int wrtOut, unsigned int wrtIn, Args const &... args) |
virtual Eigen::VectorXd | GradientByFD (unsigned int const outputDimWrt, unsigned int const inputDimWrt, std::vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sensitivity) |
virtual Eigen::VectorXd | GradientByFD (unsigned int const outputDimWrt, unsigned int const inputDimWrt, ref_vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sensitivity) |
virtual Eigen::MatrixXd | JacobianByFD (unsigned int const outputDimWrt, unsigned int const inputDimWrt, std::vector< Eigen::VectorXd > const &input) |
virtual Eigen::MatrixXd | JacobianByFD (unsigned int const outputDimWrt, unsigned int const inputDimWrt, ref_vector< Eigen::VectorXd > const &input) |
virtual Eigen::VectorXd | ApplyJacobianByFD (unsigned int const outputDimWrt, unsigned int const inputDimWrt, std::vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &vec) |
virtual Eigen::VectorXd | ApplyJacobianByFD (unsigned int const outputDimWrt, unsigned int const inputDimWrt, ref_vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &vec) |
virtual Eigen::VectorXd const & | ApplyHessian (unsigned int const outWrt, unsigned int const inWrt1, unsigned int const inWrt2, std::vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sens, Eigen::VectorXd const &vec) |
virtual Eigen::VectorXd const & | ApplyHessian (unsigned int const outWrt, unsigned int const inWrt1, unsigned int const inWrt2, ref_vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sens, Eigen::VectorXd const &vec) |
virtual Eigen::VectorXd | ApplyHessianByFD (unsigned int const outWrt, unsigned int const inWrt1, unsigned int const inWrt2, std::vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sens, Eigen::VectorXd const &vec) |
virtual Eigen::VectorXd | ApplyHessianByFD (unsigned int const outWrt, unsigned int const inWrt1, unsigned int const inWrt2, ref_vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sens, Eigen::VectorXd const &vec) |
void | EnableCache () |
void | DisableCache () |
bool | CacheStatus () const |
virtual void | SetWarnLevel (unsigned int newLevel) |
Public Member Functions inherited from muq::Modeling::WorkPiece | |
WorkPiece () | |
Create a muq::Modeling::WorkPiece with no fixed number of inputs and outputs and variable input/output types. More... | |
WorkPiece (int const num, WorkPiece::Fix const fix=WorkPiece::Fix::Inputs) | |
Create a muq::Modeling::WorkPiece with either a fixed number of inputs or outputs and variable input/output types. More... | |
WorkPiece (int const numIns, int const numOuts) | |
Create a muq::Modeling::WorkPiece with a fixed number of inputs and outputs but variable input/output types. More... | |
WorkPiece (std::vector< std::string > const &types, WorkPiece::Fix const fix=WorkPiece::Fix::Inputs) | |
Create a muq::Modeling::WorkPiece with either a fixed number of inputs with specified types or a fixed number of outputs with specified types. More... | |
WorkPiece (std::map< unsigned int, std::string > const &types, WorkPiece::Fix const fix=WorkPiece::Fix::Inputs) | |
Create a muq::Modeling::WorkPiece where either some of the inputs have specified types or some of the outputs have specified types. More... | |
WorkPiece (std::map< unsigned int, std::string > const &types, int const num, WorkPiece::Fix const fixTypes=WorkPiece::Fix::Inputs, WorkPiece::Fix const fixNum=WorkPiece::Fix::Inputs) | |
Create a muq::Modeling::WorkPiece where either some of the inputs have specified types or some of the outputs have specified types and either the number of inputs or the number of outputs is fixed. More... | |
WorkPiece (std::vector< std::string > const &types, int const num) | |
Create a muq::Modeling::WorkPiece with a fixed number of inputs with specified types and a fixed number of outputs (of uknown type) More... | |
WorkPiece (int const num, std::vector< std::string > const &types) | |
Create a muq::Modeling::WorkPiece with a fixed number of outputs with specified types and a fixed number of inputs (of uknown type) More... | |
WorkPiece (std::map< unsigned int, std::string > const &inTypes, int const numIns, int const numOuts) | |
Create a muq::Modeling::WorkPiece where some of the inputs are known and we know the input and output numbers. More... | |
WorkPiece (int const numIns, std::map< unsigned int, std::string > const &outTypes, int const numOuts) | |
Create a muq::Modeling::WorkPiece where some of the outputs are known and we know the input and output numbers. More... | |
WorkPiece (std::vector< std::string > const &inTypes, std::vector< std::string > const &outTypes) | |
Create a muq::Modeling::WorkPiece with a fixed number of inputs and outputs with specified types. More... | |
WorkPiece (std::map< unsigned int, std::string > const &inTypes, std::vector< std::string > const &outTypes) | |
Create a muq::Modeling::WorkPiece where some of the inputs are known and all of the outputs have specified types. More... | |
WorkPiece (std::map< unsigned int, std::string > const &inTypes, int const num, std::vector< std::string > const &outTypes) | |
Create a muq::Modeling::WorkPiece where some of the inputs are known with a known number of inputs and all of the outputs have specified types. More... | |
WorkPiece (std::vector< std::string > const &inTypes, std::map< unsigned int, std::string > const &outTypes) | |
Create a muq::Modeling::WorkPiece where some of the outputs and all of the inputs have specified types. More... | |
WorkPiece (std::vector< std::string > const &inTypes, std::map< unsigned int, std::string > const &outTypes, int const num) | |
Create a muq::Modeling::WorkPiece where some of the outputs with a known number of outputs and all of the inputs have specified types. More... | |
WorkPiece (std::map< unsigned int, std::string > const &inTypes, std::map< unsigned int, std::string > const &outTypes) | |
Create a muq::Mdoeling::WorkPiece where some of the inputs and some of the outputs have specified types. More... | |
WorkPiece (std::map< unsigned int, std::string > const &inTypes, int const numIn, std::map< unsigned int, std::string > const &outTypes) | |
Create a muq::Mdoeling::WorkPiece where some of the inputs and some of the outputs have specified types with a fixed number of inputs. More... | |
WorkPiece (std::map< unsigned int, std::string > const &inTypes, std::map< unsigned int, std::string > const &outTypes, int const numOut) | |
Create a muq::Mdoeling::WorkPiece where some of the inputs and some of the outputs have specified types with a fixed number of outputs. More... | |
WorkPiece (std::map< unsigned int, std::string > const &inTypes, int const numIn, std::map< unsigned int, std::string > const &outTypes, int const numOut) | |
Create a muq::Mdoeling::WorkPiece where some of the inputs and some of the outputs have specified types with a fixed number of inputs and outputs. More... | |
virtual | ~WorkPiece () |
Default destructor. More... | |
std::vector< boost::any > const & | Evaluate (std::vector< boost::any > const &ins) |
Evaluate this muq::Modeling::WorkPiece. More... | |
std::vector< boost::any > const & | Evaluate (ref_vector< boost::any > const &ins) |
Evaluate this muq::Modeling::WorkPiece using references to the inputs. More... | |
std::vector< boost::any > const & | Evaluate () |
Evaluate this muq::Modeling::WorkPiece in the case that there are no inputs. More... | |
template<typename... Args> | |
std::vector< boost::any > const & | Evaluate (Args... args) |
Evalaute this muq::Modeling::WorkPiece using multiple arguments. More... | |
std::string const & | Name () |
Get the (unique) name of this work piece. More... | |
void | SetName (std::string const &newName) |
Set the name of this work piece. More... | |
std::string | InputType (unsigned int inputNum, bool const demangle=true) const |
Get the input type (if we know it) for a specific input. More... | |
int | InputSize (unsigned int inputNum) const |
Get the length of a vector valued input with fixed size. More... | |
void | SetInputSize (unsigned int inputNum, int newSize) |
std::string | OutputType (unsigned int outputNum, bool const demangle=true) const |
Get the output type (if we know it) for a specific output. More... | |
std::map< unsigned int, std::string > | OutputTypes () const |
Get the output types. More... | |
std::map< unsigned int, std::string > | InputTypes () const |
Get the input types. More... | |
unsigned int | ID () const |
Get the unique ID number. More... | |
Static Public Member Functions | |
static std::shared_ptr< PolynomialChaosExpansion > | ComputeWeightedSum (std::vector< std::shared_ptr< PolynomialChaosExpansion >> expansions, Eigen::VectorXd const &weights) |
static std::shared_ptr< PolynomialChaosExpansion > | FromHDF5 (std::string filename, std::string groupName="/") |
Loads an expansion from an HDF5 file. More... | |
static std::shared_ptr< PolynomialChaosExpansion > | FromHDF5 (muq::Utilities::H5Object &group) |
Loads the expansion from an existing HDF5 group. More... | |
Static Public Member Functions inherited from muq::Approximation::BasisExpansion | |
static std::shared_ptr< BasisExpansion > | FromHDF5 (std::string filename, std::string groupName="/") |
Loads an expansion from an HDF5 file. More... | |
static std::shared_ptr< BasisExpansion > | FromHDF5 (muq::Utilities::H5Object &group) |
Loads the expansion from an existing HDF5 group. More... | |
Static Public Member Functions inherited from muq::Modeling::WorkPiece | |
static ref_vector< const boost::any > | ToRefVector (std::vector< boost::any > const &anyVec) |
Create vector of references from a vector of boost::any's. More... | |
static ref_vector< const Eigen::VectorXd > | ToRefVector (std::vector< Eigen::VectorXd > const &anyVec) |
Friends | |
class | PCEFactory |
Additional Inherited Members | |
Public Attributes inherited from muq::Modeling::ModPiece | |
const Eigen::VectorXi | inputSizes |
const Eigen::VectorXi | outputSizes |
Public Attributes inherited from muq::Modeling::WorkPiece | |
int | numInputs |
The number of inputs. More... | |
int | numOutputs |
The number of outputs. More... | |
PolynomialChaosExpansion::PolynomialChaosExpansion | ( | std::shared_ptr< OrthogonalPolynomial > const & | basisCompsIn, |
std::shared_ptr< muq::Utilities::MultiIndexSet > const & | multisIn, | ||
Eigen::MatrixXd const & | coeffsIn | ||
) |
Definition at line 13 of file PolynomialChaosExpansion.cpp.
PolynomialChaosExpansion::PolynomialChaosExpansion | ( | std::shared_ptr< OrthogonalPolynomial > const & | basisCompsIn, |
std::shared_ptr< muq::Utilities::MultiIndexSet > const & | multisIn, | ||
unsigned int | outputDim | ||
) |
Definition at line 20 of file PolynomialChaosExpansion.cpp.
PolynomialChaosExpansion::PolynomialChaosExpansion | ( | std::vector< std::shared_ptr< IndexedScalarBasis >> const & | basisCompsIn, |
std::shared_ptr< muq::Utilities::MultiIndexSet > const & | multisIn, | ||
Eigen::MatrixXd const & | coeffsIn | ||
) |
Definition at line 27 of file PolynomialChaosExpansion.cpp.
PolynomialChaosExpansion::PolynomialChaosExpansion | ( | std::vector< std::shared_ptr< IndexedScalarBasis >> const & | basisCompsIn, |
std::shared_ptr< muq::Utilities::MultiIndexSet > const & | multisIn, | ||
unsigned int | outputDim | ||
) |
Definition at line 33 of file PolynomialChaosExpansion.cpp.
|
virtualdefault |
|
static |
Compute the weighted sum of polynomial expansions. Slow, because it assumes you haven't been tracking which polynomials to keep, so it does that, then calls the private method. If you know already, your class should be a friend and use the private method directly, as that will be faster.
Definition at line 162 of file PolynomialChaosExpansion.cpp.
Referenced by muq::Approximation::AdaptiveSmolyakPCE::AddEstimates().
|
staticprivate |
An internal function to compute the weighted sum of polynomial expansions if you already know which polynomials are included, where polynomials is a clean copy not used by other expansions. Actually does the addition.
Definition at line 134 of file PolynomialChaosExpansion.cpp.
References muq::Approximation::BasisExpansion::coeffs, and nlohmann::detail::dtoa_impl::k.
Eigen::MatrixXd PolynomialChaosExpansion::Covariance | ( | ) | const |
Definition at line 84 of file PolynomialChaosExpansion.cpp.
References muq::Approximation::BasisExpansion::basisComps, muq::Approximation::BasisExpansion::coeffs, GetNormalizationVec(), muq::Modeling::ModPiece::inputSizes, and muq::Modeling::ModPiece::outputSizes.
|
static |
Loads the expansion from an existing HDF5 group.
This function will read the multiindices in an an HDF5 dataset and construct an instance of the MultiIndexSet class.
[in] | group | An HDF5 group containing the "multiindices", "coefficients", and "basis_types" datasets. |
Definition at line 298 of file PolynomialChaosExpansion.cpp.
References muq::Approximation::BasisExpansion::FromHDF5(), and muq::Utilities::H5Object::Path().
|
static |
Loads an expansion from an HDF5 file.
This function works in tandem with the BasisExpansion::ToHDF5 function. It will read the multiindices, coefficients, and scalar basis type from the HDF5 file and construct a BasisExpansion. See the BasisExpansion::ToHDF5 function for the details of these datasets.
[in] | filename | A string to an HDF5 file. If the file doesn't exist or the correct datasets don't exist, an exception will be thrown. |
[in] | dsetName | The path to the HDF5 group containing expansion datasets. |
Definition at line 290 of file PolynomialChaosExpansion.cpp.
References muq::Approximation::BasisExpansion::FromHDF5().
|
private |
This function returns the sqrt of the normalization, sqrt(<p*p>), for each PCE basis function p. NB: This function does not change depending on the coefficients.
Definition at line 40 of file PolynomialChaosExpansion.cpp.
References muq::Approximation::BasisExpansion::basisComps, nlohmann::detail::dtoa_impl::k, and muq::Approximation::BasisExpansion::multis.
Referenced by Covariance(), Magnitude(), SobolSensitivity(), TotalSensitivity(), and Variance().
Eigen::VectorXd PolynomialChaosExpansion::Magnitude | ( | ) | const |
Compute the L2 norm of each output.
Definition at line 124 of file PolynomialChaosExpansion.cpp.
References muq::Approximation::BasisExpansion::coeffs, and GetNormalizationVec().
Eigen::MatrixXd PolynomialChaosExpansion::MainSensitivity | ( | ) | const |
Compute all the main sensitivities (e.g., one term Sobol sensitivity indices). Rows are outputs, each column is an input.
Definition at line 273 of file PolynomialChaosExpansion.cpp.
References muq::Approximation::BasisExpansion::coeffs, muq::Modeling::ModPiece::inputSizes, and SobolSensitivity().
Eigen::VectorXd PolynomialChaosExpansion::Mean | ( | ) | const |
Definition at line 118 of file PolynomialChaosExpansion.cpp.
References muq::Approximation::BasisExpansion::coeffs.
Eigen::VectorXd PolynomialChaosExpansion::SobolSensitivity | ( | std::vector< unsigned int > const & | targetDims | ) | const |
Computes the Sobol sensitivity for a group of input parameters.
Definition at line 240 of file PolynomialChaosExpansion.cpp.
References muq::Approximation::BasisExpansion::basisComps, muq::Approximation::BasisExpansion::coeffs, GetNormalizationVec(), muq::Modeling::ModPiece::inputSizes, nlohmann::detail::dtoa_impl::k, muq::Approximation::BasisExpansion::multis, muq::Approximation::BasisExpansion::NumTerms(), and Variance().
Eigen::VectorXd PolynomialChaosExpansion::SobolSensitivity | ( | unsigned int | targetDim | ) | const |
Compute the main sensitivity index for the input dimension, for each output dimension.
Definition at line 235 of file PolynomialChaosExpansion.cpp.
Referenced by MainSensitivity().
|
overridevirtual |
Write PCE to HDF5 file.
This is the same as BasisExpansion::ToHDF5, but adds additional attribute specifiying that this is a PCE.
Reimplemented from muq::Approximation::BasisExpansion.
Definition at line 284 of file PolynomialChaosExpansion.cpp.
References muq::Utilities::H5Object::attrs, and muq::Approximation::BasisExpansion::ToHDF5().
Eigen::MatrixXd PolynomialChaosExpansion::TotalSensitivity | ( | ) | const |
Compute all Sobol total sensitivities. Rows are outputs, each column is an input.
Definition at line 225 of file PolynomialChaosExpansion.cpp.
References muq::Approximation::BasisExpansion::coeffs, and muq::Modeling::ModPiece::inputSizes.
Eigen::VectorXd PolynomialChaosExpansion::TotalSensitivity | ( | unsigned int | targetDim | ) | const |
Compute the Sobol total sensitivity index for the input dimension, for each output dimension.
Definition at line 204 of file PolynomialChaosExpansion.cpp.
References muq::Approximation::BasisExpansion::basisComps, muq::Approximation::BasisExpansion::coeffs, GetNormalizationVec(), muq::Modeling::ModPiece::inputSizes, muq::Approximation::BasisExpansion::multis, and Variance().
Eigen::VectorXd PolynomialChaosExpansion::Variance | ( | ) | const |
compute the variance of the current expansion
Definition at line 60 of file PolynomialChaosExpansion.cpp.
References muq::Approximation::BasisExpansion::basisComps, muq::Approximation::BasisExpansion::coeffs, GetNormalizationVec(), and muq::Modeling::ModPiece::inputSizes.
Referenced by SobolSensitivity(), and TotalSensitivity().
|
friend |
Definition at line 29 of file PolynomialChaosExpansion.h.