https://mooseframework.inl.gov
Public Types | Public Member Functions | Static Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | Private Attributes | List of all members
LMC Class Reference

Covariance function for multi-output Gaussian Processes based on the linear model of coregionalization (LMC) More...

#include <LMC.h>

Inheritance diagram for LMC:
[legend]

Public Types

typedef DataFileName DataFileParameterType
 

Public Member Functions

 LMC (const InputParameters &parameters)
 
void computeCovarianceMatrix (RealEigenMatrix &K, const RealEigenMatrix &x, const RealEigenMatrix &xp, const bool is_self_covariance) const override
 Generates the Covariance Matrix given two sets of points in the parameter space. More...
 
bool computedKdhyper (RealEigenMatrix &dKdhp, const RealEigenMatrix &x, const std::string &hyper_param_name, unsigned int ind) const override
 Redirect dK/dhp for hyperparameter "hp". More...
 
void loadHyperParamMap (const std::unordered_map< std::string, Real > &map, const std::unordered_map< std::string, std::vector< Real >> &vec_map)
 Load some hyperparameters into the local maps contained in this object. More...
 
void buildHyperParamMap (std::unordered_map< std::string, Real > &map, std::unordered_map< std::string, std::vector< Real >> &vec_map) const
 Populates the input maps with the owned hyperparameters. More...
 
virtual bool getTuningData (const std::string &name, unsigned int &size, Real &min, Real &max) const
 Get the default minimum and maximum and size of a hyperparameter. More...
 
void dependentCovarianceTypes (std::map< UserObjectName, std::string > &name_type_map) const
 Populate a map with the names and types of the dependent covariance functions. More...
 
const std::vector< UserObjectName > & dependentCovarianceNames () const
 Get the names of the dependent covariances. More...
 
virtual bool isTunable (const std::string &name) const
 Check if a given parameter is tunable. More...
 
unsigned int numOutputs () const
 Return the number of outputs assumed for this covariance function. More...
 
std::unordered_map< std::string, Real > & hyperParamMapReal ()
 Get the map of scalar parameters. More...
 
std::unordered_map< std::string, std::vector< Real > > & hyperParamMapVectorReal ()
 Get the map of vector parameters. More...
 
virtual bool enabled () const
 
std::shared_ptr< MooseObjectgetSharedPtr ()
 
std::shared_ptr< const MooseObjectgetSharedPtr () const
 
MooseAppgetMooseApp () const
 
const std::string & type () const
 
virtual const std::string & name () const
 
std::string typeAndName () const
 
std::string errorPrefix (const std::string &error_type) const
 
void callMooseError (std::string msg, const bool with_prefix) const
 
MooseObjectParameterName uniqueParameterName (const std::string &parameter_name) const
 
const InputParametersparameters () const
 
MooseObjectName uniqueName () const
 
const T & getParam (const std::string &name) const
 
std::vector< std::pair< T1, T2 > > getParam (const std::string &param1, const std::string &param2) const
 
const T * queryParam (const std::string &name) const
 
const T & getRenamedParam (const std::string &old_name, const std::string &new_name) const
 
getCheckedPointerParam (const std::string &name, const std::string &error_string="") const
 
bool isParamValid (const std::string &name) const
 
bool isParamSetByUser (const std::string &nm) const
 
void paramError (const std::string &param, Args... args) const
 
void paramWarning (const std::string &param, Args... args) const
 
void paramInfo (const std::string &param, Args... args) const
 
void connectControllableParams (const std::string &parameter, const std::string &object_type, const std::string &object_name, const std::string &object_parameter) const
 
void mooseError (Args &&... args) const
 
void mooseErrorNonPrefixed (Args &&... args) const
 
void mooseDocumentedError (const std::string &repo_name, const unsigned int issue_num, Args &&... args) const
 
void mooseWarning (Args &&... args) const
 
void mooseWarningNonPrefixed (Args &&... args) const
 
void mooseDeprecated (Args &&... args) const
 
void mooseInfo (Args &&... args) const
 
std::string getDataFileName (const std::string &param) const
 
std::string getDataFileNameByName (const std::string &relative_path) const
 
std::string getDataFilePath (const std::string &relative_path) const
 
const Parallel::Communicator & comm () const
 
processor_id_type n_processors () const
 
processor_id_type processor_id () const
 

Static Public Member Functions

static InputParameters validParams ()
 

Public Attributes

const ConsoleStream _console
 

Protected Member Functions

void computeBMatrix (RealEigenMatrix &Bmat, const unsigned int exp_i) const
 Computes the covariance matrix for the outputs (using the latent coefficients) We use a $B = a_i a_i^T + diag(lambda_i)$ expansion here where $a_i$ and $lambda_i$ are vectors. More...
 
void computeAGradient (RealEigenMatrix &grad, const unsigned int exp_i, const unsigned int index) const
 Computes the gradient of $B$ with respect to the entries in $a_i$ in the following expression: $B = a_i a_i^T + diag(lambda_i)$. More...
 
void computeLambdaGradient (RealEigenMatrix &grad, const unsigned int exp_i, const unsigned int index) const
 Computes the gradient of $B$ with respect to the entries in $lambda_i$ in the following expression: $B = a_i a_i^T + diag(lambda_i)$. More...
 
const RealaddRealHyperParameter (const std::string &name, const Real value, const bool is_tunable)
 Register a scalar hyperparameter to this covariance function. More...
 
const std::vector< Real > & addVectorRealHyperParameter (const std::string &name, const std::vector< Real > value, const bool is_tunable)
 Register a vector hyperparameter to this covariance function. More...
 
CovarianceFunctionBasegetCovarianceFunctionByName (const UserObjectName &name) const
 Lookup a CovarianceFunction object by name and return pointer. More...
 

Protected Attributes

const unsigned int _num_expansion_terms
 The number of expansion terms in the output ovariance matrix. More...
 
std::unordered_map< std::string, Real_hp_map_real
 Map of real-valued hyperparameters. More...
 
std::unordered_map< std::string, std::vector< Real > > _hp_map_vector_real
 Map of vector-valued hyperparameters. More...
 
std::unordered_set< std::string > _tunable_hp
 list of tunable hyper-parameters More...
 
const unsigned int _num_outputs
 The number of outputs this covariance function is used to describe. More...
 
const std::vector< UserObjectName > _dependent_covariance_names
 The names of the dependent covariance functions. More...
 
std::vector< std::string > _dependent_covariance_types
 The types of the dependent covariance functions. More...
 
std::vector< CovarianceFunctionBase * > _covariance_functions
 Vector of pointers to the dependent covariance functions. More...
 
const bool & _enabled
 
MooseApp_app
 
const std::string _type
 
const std::string _name
 
const InputParameters_pars
 
Factory_factory
 
ActionFactory_action_factory
 
const Parallel::Communicator & _communicator
 

Private Attributes

std::vector< const std::vector< Real > * > _a_coeffs
 The vectors in the $B = a_i a_i^T + diag(lambda_i)$ expansion. More...
 
std::vector< const std::vector< Real > * > _lambdas
 

Detailed Description

Covariance function for multi-output Gaussian Processes based on the linear model of coregionalization (LMC)

Definition at line 18 of file LMC.h.

Constructor & Destructor Documentation

◆ LMC()

LMC::LMC ( const InputParameters parameters)

Definition at line 29 of file LMC.C.

31  _num_expansion_terms(getParam<unsigned int>("num_latent_funcs"))
32 {
33  // We use a random number generator to obtain the initial guess for the
34  // hyperparams
35  MooseRandom generator_latent;
36  generator_latent.seed(0, 1980);
37 
38  // First add and initialize the a A coefficients in the (aa^T+lambda*I) matrix
39  for (const auto exp_i : make_range(_num_expansion_terms))
40  {
41  const std::string a_coeff_name = "acoeff_" + std::to_string(exp_i);
42  const std::string a_coeff_name_with_prefix = _name + ":" + a_coeff_name;
43  _a_coeffs.push_back(
44  &addVectorRealHyperParameter(a_coeff_name, std::vector(_num_outputs, 1.0), true));
45 
46  auto & acoeff_vector = _hp_map_vector_real[a_coeff_name_with_prefix];
47  for (const auto out_i : make_range(_num_outputs))
48  acoeff_vector[out_i] = 3.0 * generator_latent.rand(0) + 1.0;
49  }
50 
51  // Then add and initialize the lambda coefficients in the (aa^T+lambda*I) matrix
52  for (const auto exp_i : make_range(_num_expansion_terms))
53  {
54  const std::string lambda_name = "lambda_" + std::to_string(exp_i);
55  const std::string lambda_name_with_prefix = _name + ":" + lambda_name;
56  _lambdas.push_back(
57  &addVectorRealHyperParameter(lambda_name, std::vector(_num_outputs, 1.0), true));
58 
59  auto & lambda_vector = _hp_map_vector_real[lambda_name_with_prefix];
60  for (const auto out_i : make_range(_num_outputs))
61  lambda_vector[out_i] = 3.0 * generator_latent.rand(0) + 1.0;
62  }
63 }
const unsigned int _num_expansion_terms
The number of expansion terms in the output ovariance matrix.
Definition: LMC.h:67
void seed(std::size_t i, unsigned int seed)
CovarianceFunctionBase(const InputParameters &parameters)
const std::vector< Real > & addVectorRealHyperParameter(const std::string &name, const std::vector< Real > value, const bool is_tunable)
Register a vector hyperparameter to this covariance function.
std::vector< const std::vector< Real > * > _a_coeffs
The vectors in the $B = a_i a_i^T + diag(lambda_i)$ expansion.
Definition: LMC.h:72
const unsigned int _num_outputs
The number of outputs this covariance function is used to describe.
std::vector< const std::vector< Real > * > _lambdas
Definition: LMC.h:73
const std::string _name
IntRange< T > make_range(T beg, T end)
std::unordered_map< std::string, std::vector< Real > > _hp_map_vector_real
Map of vector-valued hyperparameters.
const InputParameters & parameters() const
Real rand(std::size_t i)

Member Function Documentation

◆ addRealHyperParameter()

const Real & CovarianceFunctionBase::addRealHyperParameter ( const std::string &  name,
const Real  value,
const bool  is_tunable 
)
protectedinherited

Register a scalar hyperparameter to this covariance function.

Parameters
nameThe name of the parameter
valueThe initial value of the parameter
is_tunableIf the parameter is tunable during optimization

Definition at line 52 of file CovarianceFunctionBase.C.

55 {
56  const auto prefixed_name = _name + ":" + name;
57  if (is_tunable)
58  _tunable_hp.insert(prefixed_name);
59  return _hp_map_real.emplace(prefixed_name, value).first->second;
60 }
std::unordered_map< std::string, Real > _hp_map_real
Map of real-valued hyperparameters.
std::unordered_set< std::string > _tunable_hp
list of tunable hyper-parameters
virtual const std::string & name() const
const std::string _name

◆ addVectorRealHyperParameter()

const std::vector< Real > & CovarianceFunctionBase::addVectorRealHyperParameter ( const std::string &  name,
const std::vector< Real value,
const bool  is_tunable 
)
protectedinherited

Register a vector hyperparameter to this covariance function.

Parameters
nameThe name of the parameter
valueThe initial value of the parameter
is_tunableIf the parameter is tunable during optimization

Definition at line 63 of file CovarianceFunctionBase.C.

Referenced by LMC().

66 {
67  const auto prefixed_name = _name + ":" + name;
68  if (is_tunable)
69  _tunable_hp.insert(prefixed_name);
70  return _hp_map_vector_real.emplace(prefixed_name, value).first->second;
71 }
std::unordered_set< std::string > _tunable_hp
list of tunable hyper-parameters
virtual const std::string & name() const
const std::string _name
std::unordered_map< std::string, std::vector< Real > > _hp_map_vector_real
Map of vector-valued hyperparameters.

◆ buildHyperParamMap()

void CovarianceFunctionBase::buildHyperParamMap ( std::unordered_map< std::string, Real > &  map,
std::unordered_map< std::string, std::vector< Real >> &  vec_map 
) const
inherited

Populates the input maps with the owned hyperparameters.

Parameters
mapMap of scalar hyperparameters that should be populated
vec_mapMap of vector hyperparameters that should be populated

Definition at line 115 of file CovarianceFunctionBase.C.

Referenced by StochasticTools::GaussianProcess::setupCovarianceMatrix(), and StochasticTools::GaussianProcess::tuneHyperParamsAdam().

118 {
119  // First, add the hyperparameters of the dependent covariance functions
120  for (const auto dependent_covar : _covariance_functions)
121  dependent_covar->buildHyperParamMap(map, vec_map);
122 
123  // At the end we just append the hyperparameters this object owns
124  for (const auto & iter : _hp_map_real)
125  map[iter.first] = iter.second;
126  for (const auto & iter : _hp_map_vector_real)
127  vec_map[iter.first] = iter.second;
128 }
std::unordered_map< std::string, Real > _hp_map_real
Map of real-valued hyperparameters.
std::vector< CovarianceFunctionBase * > _covariance_functions
Vector of pointers to the dependent covariance functions.
std::unordered_map< std::string, std::vector< Real > > _hp_map_vector_real
Map of vector-valued hyperparameters.

◆ computeAGradient()

void LMC::computeAGradient ( RealEigenMatrix grad,
const unsigned int  exp_i,
const unsigned int  index 
) const
protected

Computes the gradient of $B$ with respect to the entries in $a_i$ in the following expression: $B = a_i a_i^T + diag(lambda_i)$.

Parameters
gradThe gradient matrix that should be populated
exp_iThe index of the expansion of B
indexThe index within the vector $a_i$

Definition at line 176 of file LMC.C.

Referenced by computedKdhyper().

179 {
180  const auto & a_coeffs = *_a_coeffs[exp_i];
181  // Add asserts here
182  grad.setZero();
183  for (const auto col_i : make_range(_num_outputs))
184  grad(index, col_i) = a_coeffs[col_i];
185  const RealEigenMatrix transpose = grad.transpose();
186  grad = grad + transpose;
187 }
std::vector< const std::vector< Real > * > _a_coeffs
The vectors in the $B = a_i a_i^T + diag(lambda_i)$ expansion.
Definition: LMC.h:72
const unsigned int _num_outputs
The number of outputs this covariance function is used to describe.
Eigen::Matrix< Real, Eigen::Dynamic, Eigen::Dynamic > RealEigenMatrix
std::string grad(const std::string &var)
Definition: NS.h:91
IntRange< T > make_range(T beg, T end)

◆ computeBMatrix()

void LMC::computeBMatrix ( RealEigenMatrix Bmat,
const unsigned int  exp_i 
) const
protected

Computes the covariance matrix for the outputs (using the latent coefficients) We use a $B = a_i a_i^T + diag(lambda_i)$ expansion here where $a_i$ and $lambda_i$ are vectors.

Parameters
BmatThe matrix which should be populated by the covariance values
exp_iThe expansion index in the latent space

Definition at line 161 of file LMC.C.

Referenced by computeCovarianceMatrix(), and computedKdhyper().

162 {
163  const auto & a_coeffs = *_a_coeffs[exp_i];
164  const auto & lambda_coeffs = *_lambdas[exp_i];
165 
166  for (const auto row_i : make_range(_num_outputs))
167  for (const auto col_i : make_range(_num_outputs))
168  {
169  Bmat(row_i, col_i) = a_coeffs[row_i] * a_coeffs[col_i];
170  if (row_i == col_i)
171  Bmat(row_i, col_i) += lambda_coeffs[col_i];
172  }
173 }
std::vector< const std::vector< Real > * > _a_coeffs
The vectors in the $B = a_i a_i^T + diag(lambda_i)$ expansion.
Definition: LMC.h:72
const unsigned int _num_outputs
The number of outputs this covariance function is used to describe.
std::vector< const std::vector< Real > * > _lambdas
Definition: LMC.h:73
IntRange< T > make_range(T beg, T end)

◆ computeCovarianceMatrix()

void LMC::computeCovarianceMatrix ( RealEigenMatrix K,
const RealEigenMatrix x,
const RealEigenMatrix xp,
const bool  is_self_covariance 
) const
overridevirtual

Generates the Covariance Matrix given two sets of points in the parameter space.

Parameters
KReference to a matrix which should be populated by the covariance entries
xReference to the first set of points
xpReference to the second set of points
is_self_covarianceSwitch to enable adding the noise variance to the diagonal of the covariance matrix

Implements CovarianceFunctionBase.

Definition at line 66 of file LMC.C.

70 {
71  // Create temporary vectors for constructing the covariance matrix
72  RealEigenMatrix K_params = RealEigenMatrix::Zero(x.rows(), xp.rows());
73  RealEigenMatrix B = RealEigenMatrix::Zero(_num_outputs, _num_outputs);
74  K = RealEigenMatrix::Zero(x.rows() * _num_outputs, xp.rows() * _num_outputs);
75  RealEigenMatrix K_working =
76  RealEigenMatrix::Zero(x.rows() * _num_outputs, xp.rows() * _num_outputs);
77 
78  // For every expansion term we add the contribution to the covariance matrix
79  for (const auto exp_i : make_range(_num_expansion_terms))
80  {
81  _covariance_functions[exp_i]->computeCovarianceMatrix(K_params, x, xp, is_self_covariance);
82  computeBMatrix(B, exp_i);
83  MathUtils::kron(K_working, B, K_params);
84  K += K_working;
85  }
86 }
void computeBMatrix(RealEigenMatrix &Bmat, const unsigned int exp_i) const
Computes the covariance matrix for the outputs (using the latent coefficients) We use a $B = a_i a_i...
Definition: LMC.C:161
const unsigned int _num_expansion_terms
The number of expansion terms in the output ovariance matrix.
Definition: LMC.h:67
void kron(RealEigenMatrix &product, const RealEigenMatrix &mat_A, const RealEigenMatrix &mat_B)
static const std::string K
Definition: NS.h:170
std::vector< CovarianceFunctionBase * > _covariance_functions
Vector of pointers to the dependent covariance functions.
const unsigned int _num_outputs
The number of outputs this covariance function is used to describe.
const std::vector< double > x
Eigen::Matrix< Real, Eigen::Dynamic, Eigen::Dynamic > RealEigenMatrix
IntRange< T > make_range(T beg, T end)

◆ computedKdhyper()

bool LMC::computedKdhyper ( RealEigenMatrix dKdhp,
const RealEigenMatrix x,
const std::string &  hyper_param_name,
unsigned int  ind 
) const
overridevirtual

Redirect dK/dhp for hyperparameter "hp".

Returns false is the parameter has not been found in this covariance object.

Parameters
dKdhpThe matrix which should be populated with the derivatives
xThe input vector for which the derivatives of the covariance matrix is computed
hyper_param_nameThe name of the hyperparameter
indThe index within the hyperparameter. 0 if it is a scalar parameter. If it is a vector parameter, it should be the index within the vector.

Reimplemented from CovarianceFunctionBase.

Definition at line 89 of file LMC.C.

93 {
94  // Early return in the paramter name is longer than the expected [name] prefix.
95  // We prefix the parameter names with the name of the covariance function.
96  if (name().length() + 1 > hyper_param_name.length())
97  return false;
98 
99  // Strip the prefix from the given parameter name
100  const std::string name_without_prefix = hyper_param_name.substr(name().length() + 1);
101 
102  // Check if the parameter is tunable
103  if (_tunable_hp.find(hyper_param_name) != _tunable_hp.end())
104  {
105  const std::string acoeff_prefix = "acoeff_";
106  const std::string lambda_prefix = "lambda_";
107 
108  // Allocate storage for the factors of the total gradient matrix
109  RealEigenMatrix dBdhp = RealEigenMatrix::Zero(_num_outputs, _num_outputs);
110  RealEigenMatrix K_params = RealEigenMatrix::Zero(x.rows(), x.rows());
111 
112  if (name_without_prefix.find(acoeff_prefix) != std::string::npos)
113  {
114  // Automatically grab the expansion index
115  const int number = std::stoi(name_without_prefix.substr(acoeff_prefix.length()));
116  computeAGradient(dBdhp, number, ind);
117  _covariance_functions[number]->computeCovarianceMatrix(K_params, x, x, true);
118  }
119  else if (name_without_prefix.find(lambda_prefix) != std::string::npos)
120  {
121  // Automatically grab the expansion index
122  const int number = std::stoi(name_without_prefix.substr(lambda_prefix.length()));
123  computeLambdaGradient(dBdhp, number, ind);
124  _covariance_functions[number]->computeCovarianceMatrix(K_params, x, x, true);
125  }
126  MathUtils::kron(dKdhp, dBdhp, K_params);
127  return true;
128  }
129  else
130  {
131  // Allocate storage for the matrix factors
132  RealEigenMatrix B_tmp = RealEigenMatrix::Zero(_num_outputs, _num_outputs);
133  RealEigenMatrix B = RealEigenMatrix::Zero(_num_outputs, _num_outputs);
134  RealEigenMatrix dKdhp_sub = RealEigenMatrix::Zero(x.rows(), x.rows());
135 
136  // First, check the dependent covariances
137  bool found = false;
138  for (const auto dependent_covar : _covariance_functions)
139  if (!found)
140  found = dependent_covar->computedKdhyper(dKdhp_sub, x, hyper_param_name, ind);
141 
142  if (!found)
143  mooseError("Hyperparameter ", hyper_param_name, "not found!");
144 
145  // Then we compute the output covariance
146  for (const auto exp_i : make_range(_num_expansion_terms))
147  {
148  computeBMatrix(B_tmp, exp_i);
149  B += B_tmp;
150  }
151 
152  MathUtils::kron(dKdhp, B, dKdhp_sub);
153 
154  return true;
155  }
156 
157  return false;
158 }
void computeBMatrix(RealEigenMatrix &Bmat, const unsigned int exp_i) const
Computes the covariance matrix for the outputs (using the latent coefficients) We use a $B = a_i a_i...
Definition: LMC.C:161
std::unordered_set< std::string > _tunable_hp
list of tunable hyper-parameters
void computeAGradient(RealEigenMatrix &grad, const unsigned int exp_i, const unsigned int index) const
Computes the gradient of $B$ with respect to the entries in $a_i$ in the following expression: $B = ...
Definition: LMC.C:176
const unsigned int _num_expansion_terms
The number of expansion terms in the output ovariance matrix.
Definition: LMC.h:67
void kron(RealEigenMatrix &product, const RealEigenMatrix &mat_A, const RealEigenMatrix &mat_B)
std::vector< CovarianceFunctionBase * > _covariance_functions
Vector of pointers to the dependent covariance functions.
virtual const std::string & name() const
const unsigned int _num_outputs
The number of outputs this covariance function is used to describe.
const std::vector< double > x
Eigen::Matrix< Real, Eigen::Dynamic, Eigen::Dynamic > RealEigenMatrix
IntRange< T > make_range(T beg, T end)
void mooseError(Args &&... args) const
void computeLambdaGradient(RealEigenMatrix &grad, const unsigned int exp_i, const unsigned int index) const
Computes the gradient of $B$ with respect to the entries in $lambda_i$ in the following expression: $...
Definition: LMC.C:190

◆ computeLambdaGradient()

void LMC::computeLambdaGradient ( RealEigenMatrix grad,
const unsigned int  exp_i,
const unsigned int  index 
) const
protected

Computes the gradient of $B$ with respect to the entries in $lambda_i$ in the following expression: $B = a_i a_i^T + diag(lambda_i)$.

Parameters
gradThe gradient matrix that should be populated
exp_iThe index of the expansion of B
indexThe index within the vector $lambda_i$

Definition at line 190 of file LMC.C.

Referenced by computedKdhyper().

193 {
194  grad.setZero();
195  grad(index, index) = 1.0;
196 }
std::string grad(const std::string &var)
Definition: NS.h:91

◆ dependentCovarianceNames()

const std::vector<UserObjectName>& CovarianceFunctionBase::dependentCovarianceNames ( ) const
inlineinherited

Get the names of the dependent covariances.

Definition at line 62 of file CovarianceFunctionBase.h.

Referenced by StochasticTools::GaussianProcess::linkCovarianceFunction().

63  {
65  }
const std::vector< UserObjectName > _dependent_covariance_names
The names of the dependent covariance functions.

◆ dependentCovarianceTypes()

void CovarianceFunctionBase::dependentCovarianceTypes ( std::map< UserObjectName, std::string > &  name_type_map) const
inherited

Populate a map with the names and types of the dependent covariance functions.

Parameters
name_type_mapReference to the map which should be populated

Definition at line 163 of file CovarianceFunctionBase.C.

Referenced by StochasticTools::GaussianProcess::linkCovarianceFunction().

165 {
166  for (const auto dependent_covar : _covariance_functions)
167  {
168  dependent_covar->dependentCovarianceTypes(name_type_map);
169  name_type_map.insert(std::make_pair(dependent_covar->name(), dependent_covar->type()));
170  }
171 }
std::vector< CovarianceFunctionBase * > _covariance_functions
Vector of pointers to the dependent covariance functions.

◆ getCovarianceFunctionByName()

CovarianceFunctionBase * CovarianceInterface::getCovarianceFunctionByName ( const UserObjectName &  name) const
protectedinherited

Lookup a CovarianceFunction object by name and return pointer.

Definition at line 25 of file CovarianceInterface.C.

Referenced by ActiveLearningGaussianProcess::ActiveLearningGaussianProcess(), CovarianceFunctionBase::CovarianceFunctionBase(), GaussianProcessTrainer::GaussianProcessTrainer(), and GaussianProcessSurrogate::setupCovariance().

26 {
27  std::vector<CovarianceFunctionBase *> models;
29  .query()
30  .condition<AttribName>(name)
31  .condition<AttribSystem>("CovarianceFunction")
32  .queryInto(models);
33  if (models.empty())
34  mooseError("Unable to find a CovarianceFunction object with the name '" + name + "'");
35  return models[0];
36 }
void mooseError(Args &&... args)
TheWarehouse & theWarehouse() const
const std::string name
Definition: Setup.h:20
Query query()
FEProblemBase & _covar_feproblem
Reference to FEProblemBase instance.

◆ getTuningData()

bool CovarianceFunctionBase::getTuningData ( const std::string &  name,
unsigned int size,
Real min,
Real max 
) const
virtualinherited

Get the default minimum and maximum and size of a hyperparameter.

Returns false is the parameter has not been found in this covariance object.

Parameters
nameThe name of the hyperparameter
sizeReference to an unsigned int that will contain the size of the hyperparameter (will be populated with 1 if it is scalar)
minReference to a number which will be populated by the maximum allowed value of the hyperparameter
maxReference to a number which will be populated by the minimum allowed value of the hyperparameter

Definition at line 131 of file CovarianceFunctionBase.C.

Referenced by StochasticTools::GaussianProcess::generateTuningMap().

135 {
136  // First, check the dependent covariances
137  for (const auto dependent_covar : _covariance_functions)
138  if (dependent_covar->getTuningData(name, size, min, max))
139  return true;
140 
141  min = 1e-9;
142  max = 1e9;
143 
144  if (_hp_map_real.find(name) != _hp_map_real.end())
145  {
146  size = 1;
147  return true;
148  }
149  else if (_hp_map_vector_real.find(name) != _hp_map_vector_real.end())
150  {
151  const auto & vector_value = _hp_map_vector_real.find(name);
152  size = vector_value->second.size();
153  return true;
154  }
155  else
156  {
157  size = 0;
158  return false;
159  }
160 }
std::unordered_map< std::string, Real > _hp_map_real
Map of real-valued hyperparameters.
std::vector< CovarianceFunctionBase * > _covariance_functions
Vector of pointers to the dependent covariance functions.
virtual const std::string & name() const
auto max(const L &left, const R &right)
std::unordered_map< std::string, std::vector< Real > > _hp_map_vector_real
Map of vector-valued hyperparameters.
auto min(const L &left, const R &right)

◆ hyperParamMapReal()

std::unordered_map<std::string, Real>& CovarianceFunctionBase::hyperParamMapReal ( )
inlineinherited

Get the map of scalar parameters.

Definition at line 88 of file CovarianceFunctionBase.h.

88 { return _hp_map_real; }
std::unordered_map< std::string, Real > _hp_map_real
Map of real-valued hyperparameters.

◆ hyperParamMapVectorReal()

std::unordered_map<std::string, std::vector<Real> >& CovarianceFunctionBase::hyperParamMapVectorReal ( )
inlineinherited

Get the map of vector parameters.

Definition at line 91 of file CovarianceFunctionBase.h.

92  {
93  return _hp_map_vector_real;
94  }
std::unordered_map< std::string, std::vector< Real > > _hp_map_vector_real
Map of vector-valued hyperparameters.

◆ isTunable()

bool CovarianceFunctionBase::isTunable ( const std::string &  name) const
virtualinherited

Check if a given parameter is tunable.

Parameters
Thename of the hyperparameter

Definition at line 74 of file CovarianceFunctionBase.C.

Referenced by StochasticTools::GaussianProcess::generateTuningMap().

75 {
76  // First, we check if the dependent covariances have the parameter
77  for (const auto dependent_covar : _covariance_functions)
78  if (dependent_covar->isTunable(name))
79  return true;
80 
81  if (_tunable_hp.find(name) != _tunable_hp.end())
82  return true;
83  else if (_hp_map_real.find(name) != _hp_map_real.end() ||
85  mooseError("We found hyperparameter ", name, " but it was not declared tunable!");
86 
87  return false;
88 }
std::unordered_map< std::string, Real > _hp_map_real
Map of real-valued hyperparameters.
std::unordered_set< std::string > _tunable_hp
list of tunable hyper-parameters
std::vector< CovarianceFunctionBase * > _covariance_functions
Vector of pointers to the dependent covariance functions.
virtual const std::string & name() const
std::unordered_map< std::string, std::vector< Real > > _hp_map_vector_real
Map of vector-valued hyperparameters.
void mooseError(Args &&... args) const

◆ loadHyperParamMap()

void CovarianceFunctionBase::loadHyperParamMap ( const std::unordered_map< std::string, Real > &  map,
const std::unordered_map< std::string, std::vector< Real >> &  vec_map 
)
inherited

Load some hyperparameters into the local maps contained in this object.

Parameters
mapInput map of scalar hyperparameters
vec_mapInput map of vector hyperparameters

Definition at line 91 of file CovarianceFunctionBase.C.

Referenced by LoadCovarianceDataAction::load(), and StochasticTools::GaussianProcess::tuneHyperParamsAdam().

94 {
95  // First, load the hyperparameters of the dependent covariance functions
96  for (const auto dependent_covar : _covariance_functions)
97  dependent_covar->loadHyperParamMap(map, vec_map);
98 
99  // Then we load the hyperparameters of this object
100  for (auto & iter : _hp_map_real)
101  {
102  const auto & map_iter = map.find(iter.first);
103  if (map_iter != map.end())
104  iter.second = map_iter->second;
105  }
106  for (auto & iter : _hp_map_vector_real)
107  {
108  const auto & map_iter = vec_map.find(iter.first);
109  if (map_iter != vec_map.end())
110  iter.second = map_iter->second;
111  }
112 }
std::unordered_map< std::string, Real > _hp_map_real
Map of real-valued hyperparameters.
std::vector< CovarianceFunctionBase * > _covariance_functions
Vector of pointers to the dependent covariance functions.
std::unordered_map< std::string, std::vector< Real > > _hp_map_vector_real
Map of vector-valued hyperparameters.

◆ numOutputs()

unsigned int CovarianceFunctionBase::numOutputs ( ) const
inlineinherited

Return the number of outputs assumed for this covariance function.

Definition at line 85 of file CovarianceFunctionBase.h.

Referenced by GaussianProcessSurrogate::evaluate(), GaussianProcessTrainer::GaussianProcessTrainer(), and StochasticTools::GaussianProcess::linkCovarianceFunction().

85 { return _num_outputs; }
const unsigned int _num_outputs
The number of outputs this covariance function is used to describe.

◆ validParams()

InputParameters LMC::validParams ( )
static

Definition at line 17 of file LMC.C.

18 {
20  params.addClassDescription("Covariance function for multioutput Gaussian Processes based on the "
21  "Linear Model of Coregionalization (LMC).");
22  params.addParam<unsigned int>(
23  "num_latent_funcs", 1., "The number of latent functions for the expansion of the outputs.");
24  params.makeParamRequired<unsigned int>("num_outputs");
25  params.makeParamRequired<std::vector<UserObjectName>>("covariance_functions");
26  return params;
27 }
void addParam(const std::string &name, const std::initializer_list< typename T::value_type > &value, const std::string &doc_string)
static InputParameters validParams()
void makeParamRequired(const std::string &name)
void addClassDescription(const std::string &doc_string)

Member Data Documentation

◆ _a_coeffs

std::vector<const std::vector<Real> *> LMC::_a_coeffs
private

The vectors in the $B = a_i a_i^T + diag(lambda_i)$ expansion.

Definition at line 72 of file LMC.h.

Referenced by computeAGradient(), computeBMatrix(), and LMC().

◆ _covariance_functions

std::vector<CovarianceFunctionBase *> CovarianceFunctionBase::_covariance_functions
protectedinherited

◆ _dependent_covariance_names

const std::vector<UserObjectName> CovarianceFunctionBase::_dependent_covariance_names
protectedinherited

The names of the dependent covariance functions.

Definition at line 125 of file CovarianceFunctionBase.h.

Referenced by CovarianceFunctionBase::CovarianceFunctionBase(), and CovarianceFunctionBase::dependentCovarianceNames().

◆ _dependent_covariance_types

std::vector<std::string> CovarianceFunctionBase::_dependent_covariance_types
protectedinherited

The types of the dependent covariance functions.

Definition at line 128 of file CovarianceFunctionBase.h.

Referenced by CovarianceFunctionBase::CovarianceFunctionBase().

◆ _hp_map_real

std::unordered_map<std::string, Real> CovarianceFunctionBase::_hp_map_real
protectedinherited

◆ _hp_map_vector_real

std::unordered_map<std::string, std::vector<Real> > CovarianceFunctionBase::_hp_map_vector_real
protectedinherited

◆ _lambdas

std::vector<const std::vector<Real> *> LMC::_lambdas
private

Definition at line 73 of file LMC.h.

Referenced by computeBMatrix(), and LMC().

◆ _num_expansion_terms

const unsigned int LMC::_num_expansion_terms
protected

The number of expansion terms in the output ovariance matrix.

Definition at line 67 of file LMC.h.

Referenced by computeCovarianceMatrix(), computedKdhyper(), and LMC().

◆ _num_outputs

const unsigned int CovarianceFunctionBase::_num_outputs
protectedinherited

The number of outputs this covariance function is used to describe.

Definition at line 122 of file CovarianceFunctionBase.h.

Referenced by computeAGradient(), computeBMatrix(), computeCovarianceMatrix(), computedKdhyper(), LMC(), and CovarianceFunctionBase::numOutputs().

◆ _tunable_hp

std::unordered_set<std::string> CovarianceFunctionBase::_tunable_hp
protectedinherited

The documentation for this class was generated from the following files: