27 #include <queso/GslVector.h>
28 #include <queso/GslMatrix.h>
29 #include <queso/VectorSet.h>
30 #include <queso/GaussianLikelihoodBlockDiagonalCovariance.h>
34 template<
class V,
class M>
39 m_covarianceCoefficients(covariance.numBlocks(), 1.0),
40 m_covariance(covariance)
42 unsigned int totalDim = 0;
48 if (totalDim != observations.sizeLocal()) {
49 queso_error_msg(
"Covariance matrix not same dimension as observation vector");
53 template<
class V,
class M>
58 template<
class V,
class M>
63 return this->m_covarianceCoefficients[i];
66 template<
class V,
class M>
71 return this->m_covarianceCoefficients[i];
74 template<
class V,
class M>
77 const V & domainVector,
const V * domainDirection, V * gradVector,
78 M * hessianMatrix, V * hessianEffect)
const
80 return std::exp(this->lnValue(domainVector, domainDirection, gradVector,
81 hessianMatrix, hessianEffect));
84 template<
class V,
class M>
87 const V & domainVector,
const V * domainDirection, V * gradVector,
88 M * hessianMatrix, V * hessianEffect)
const
90 V modelOutput(this->m_observations, 0, 0);
91 V weightedMisfit(this->m_observations, 0, 0);
93 this->evaluateModel(domainVector, domainDirection, modelOutput, gradVector,
94 hessianMatrix, hessianEffect);
97 modelOutput -= this->m_observations;
100 this->m_covariance.invertMultiply(modelOutput, weightedMisfit);
103 unsigned int offset = 0;
106 for (
unsigned int i = 0; i < this->m_covariance.numBlocks(); i++) {
108 unsigned int blockDim = this->m_covariance.getBlock(i).numRowsLocal();
109 for (
unsigned int j = 0; j < blockDim; j++) {
111 modelOutput[offset+j] /= this->m_covarianceCoefficients[i];
117 modelOutput *= weightedMisfit;
119 double norm2_squared = modelOutput.sumOfComponents();
121 return -0.5 * norm2_squared;
double & blockCoefficient(unsigned int i)
Get (non-const) multiplicative coefficient for block i.
A templated class for handling sets.
virtual double lnValue(const V &domainVector, const V *domainDirection, V *gradVector, M *hessianMatrix, V *hessianEffect) const
Logarithm of the value of the scalar function.
unsigned int numRowsLocal() const
Returns the local row dimension of this matrix.
unsigned int numBlocks() const
Return the number of blocks in the block diagonal matrix.
virtual ~GaussianLikelihoodBlockDiagonalCovariance()
Destructor.
#define queso_error_msg(msg)
GaussianLikelihoodBlockDiagonalCovariance(const char *prefix, const VectorSet< V, M > &domainSet, const V &observations, const GslBlockMatrix &covariance)
Default constructor.
const double & getBlockCoefficient(unsigned int i) const
Get (const) multiplicative coefficient for block i.
GslMatrix & getBlock(unsigned int i) const
Return block i in the block diagonal matrix.
virtual double actualValue(const V &domainVector, const V *domainDirection, V *gradVector, M *hessianMatrix, V *hessianEffect) const
Actual value of the scalar function.
Class for representing block matrices using GSL library.
Base class for canned Gaussian likelihoods.
A class representing a Gaussian likelihood with block-diagonal covariance matrix. ...
const GslBlockMatrix & m_covariance