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;
virtual ~GaussianLikelihoodBlockDiagonalCovariance()
Destructor.
Class for representing block matrices using GSL library.
virtual double actualValue(const V &domainVector, const V *domainDirection, V *gradVector, M *hessianMatrix, V *hessianEffect) const
Actual value of the scalar function.
GslMatrix & getBlock(unsigned int i) const
Return block i in the block diagonal matrix.
GaussianLikelihoodBlockDiagonalCovariance(const char *prefix, const VectorSet< V, M > &domainSet, const V &observations, const GslBlockMatrix &covariance)
Default constructor.
A class representing a Gaussian likelihood with block-diagonal covariance matrix. ...
const double & getBlockCoefficient(unsigned int i) const
Get (const) multiplicative coefficient for block i.
unsigned int numBlocks() const
Return the number of blocks in the block diagonal matrix.
unsigned int numRowsLocal() const
Returns the local row dimension of this matrix.
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.
const GslBlockMatrix & m_covariance
double & blockCoefficient(unsigned int i)
Get (non-const) multiplicative coefficient for block i.
#define queso_error_msg(msg)
Base class for canned Gaussian likelihoods.