queso-0.57.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M > Class Template Reference

A class for handling vector Gaussian random fields (GRF). More...

#include <VectorGaussianRandomField.h>

Public Member Functions

Constructor/Destructor methods
 VectorGaussianRandomField (const char *prefix, const VectorSet< P_V, P_M > &indexSet, const VectorSet< Q_V, Q_M > &imageSetPerIndex, const BaseVectorFunction< P_V, P_M, Q_V, Q_M > &meanFunction, const BaseMatrixCovarianceFunction< P_V, P_M, Q_V, Q_M > &covarianceFunction)
 Constructor. More...
 
 VectorGaussianRandomField (const VectorGaussianRandomField &obj)
 TODO: Copy constructor. More...
 
 ~VectorGaussianRandomField ()
 Destructor. More...
 
Set methods
VectorGaussianRandomFieldoperator= (const VectorGaussianRandomField &rhs)
 TODO: Assignment operator; it copies rhs to this. More...
 
Math methods
const VectorSet< P_V, P_M > & indexSet () const
 Index set; access to protected attribute m_indexSet. More...
 
const BaseVectorFunction< P_V,
P_M, Q_V, Q_M > & 
meanFunction () const
 Mean function; access to protected attribute m_meanFunction. More...
 
const
BaseMatrixCovarianceFunction
< P_V, P_M, Q_V, Q_M > & 
covarianceFunction () const
 Covariance function; access to protected attribute m_covarianceFunction. More...
 
void sampleFunction (const std::vector< P_V * > &fieldPositions, Q_V &sampleValues)
 Function that samples from a Gaussian PDF. More...
 

Protected Member Functions

void copy (const VectorGaussianRandomField &src)
 Copy method. More...
 

Protected Attributes

const BaseEnvironmentm_env
 Environment. More...
 
std::string m_prefix
 Prefix. More...
 
const VectorSet< P_V, P_M > & m_indexSet
 Index set. More...
 
const VectorSet< Q_V, Q_M > & m_imageSetPerIndex
 Image set of the RV, per index. More...
 
const BaseVectorFunction< P_V,
P_M, Q_V, Q_M > & 
m_meanFunction
 Mean function. More...
 
const
BaseMatrixCovarianceFunction
< P_V, P_M, Q_V, Q_M > & 
m_covarianceFunction
 Covariance function. More...
 
std::vector< P_V * > m_savedPositions
 Saved positions. More...
 
VectorSpace< Q_V, Q_M > * m_savedRvImageSpace
 Image set of the RV. More...
 
Q_V * m_savedRvLawExpVector
 Vector of the mean value of the RV. More...
 
Q_M * m_savedRvLawCovMatrix
 Covariance matrix of the RV. More...
 
GaussianVectorRV< Q_V, Q_M > * m_savedRv
 My RV. More...
 

Detailed Description

template<class P_V = GslVector, class P_M = GslMatrix, class Q_V = GslVector, class Q_M = GslMatrix>
class QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >

A class for handling vector Gaussian random fields (GRF).

This class implements a vector Gaussian random field (GRF); i.e. a random field involving vector Gaussian probability density functions (PDFs) of the variables.

Definition at line 48 of file VectorGaussianRandomField.h.

Constructor & Destructor Documentation

template<class P_V , class P_M , class Q_V , class Q_M >
QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::VectorGaussianRandomField ( const char *  prefix,
const VectorSet< P_V, P_M > &  indexSet,
const VectorSet< Q_V, Q_M > &  imageSetPerIndex,
const BaseVectorFunction< P_V, P_M, Q_V, Q_M > &  meanFunction,
const BaseMatrixCovarianceFunction< P_V, P_M, Q_V, Q_M > &  covarianceFunction 
)

Constructor.

Constructs a new object, given a prefix, an index set, and both a mean and a covariance function. This method deletes the previous saved positions.

Definition at line 31 of file VectorGaussianRandomField.C.

References QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::m_savedPositions.

37  :
38  m_env (indexSet.env()),
39  m_prefix ((std::string)(prefix)+"grf_"),
40  m_indexSet (indexSet),
41  m_imageSetPerIndex (imageSetPerIndex),
44  m_savedRvImageSpace (NULL),
47  m_savedRv (NULL)
48 {
49  m_savedPositions.clear();
50 }
Q_V * m_savedRvLawExpVector
Vector of the mean value of the RV.
Q_M * m_savedRvLawCovMatrix
Covariance matrix of the RV.
const VectorSet< Q_V, Q_M > & m_imageSetPerIndex
Image set of the RV, per index.
const BaseVectorFunction< P_V, P_M, Q_V, Q_M > & m_meanFunction
Mean function.
const BaseEnvironment & m_env
Environment.
const BaseVectorFunction< P_V, P_M, Q_V, Q_M > & meanFunction() const
Mean function; access to protected attribute m_meanFunction.
GaussianVectorRV< Q_V, Q_M > * m_savedRv
My RV.
const BaseMatrixCovarianceFunction< P_V, P_M, Q_V, Q_M > & covarianceFunction() const
Covariance function; access to protected attribute m_covarianceFunction.
const VectorSet< P_V, P_M > & m_indexSet
Index set.
VectorSpace< Q_V, Q_M > * m_savedRvImageSpace
Image set of the RV.
std::vector< P_V * > m_savedPositions
Saved positions.
const BaseMatrixCovarianceFunction< P_V, P_M, Q_V, Q_M > & m_covarianceFunction
Covariance function.
const BaseEnvironment & env() const
Environment. Access to private attribute m_env.
Definition: VectorSet.C:78
template<class P_V = GslVector, class P_M = GslMatrix, class Q_V = GslVector, class Q_M = GslMatrix>
QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::VectorGaussianRandomField ( const VectorGaussianRandomField< P_V, P_M, Q_V, Q_M > &  obj)

TODO: Copy constructor.

Todo:
: implement me!
template<class P_V , class P_M , class Q_V , class Q_M >
QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::~VectorGaussianRandomField ( )

Destructor.

Definition at line 53 of file VectorGaussianRandomField.C.

54 {
55 }

Member Function Documentation

template<class P_V = GslVector, class P_M = GslMatrix, class Q_V = GslVector, class Q_M = GslMatrix>
void QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::copy ( const VectorGaussianRandomField< P_V, P_M, Q_V, Q_M > &  src)
protected

Copy method.

template<class P_V , class P_M , class Q_V , class Q_M >
const BaseMatrixCovarianceFunction< P_V, P_M, Q_V, Q_M > & QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::covarianceFunction ( ) const

Covariance function; access to protected attribute m_covarianceFunction.

Definition at line 73 of file VectorGaussianRandomField.C.

74 {
75  return m_covarianceFunction;
76 }
const BaseMatrixCovarianceFunction< P_V, P_M, Q_V, Q_M > & m_covarianceFunction
Covariance function.
template<class P_V , class P_M , class Q_V , class Q_M >
const VectorSet< P_V, P_M > & QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::indexSet ( ) const

Index set; access to protected attribute m_indexSet.

Definition at line 59 of file VectorGaussianRandomField.C.

60 {
61  return m_indexSet;
62 }
const VectorSet< P_V, P_M > & m_indexSet
Index set.
template<class P_V , class P_M , class Q_V , class Q_M >
const BaseVectorFunction< P_V, P_M, Q_V, Q_M > & QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::meanFunction ( ) const

Mean function; access to protected attribute m_meanFunction.

Definition at line 66 of file VectorGaussianRandomField.C.

67 {
68  return m_meanFunction;
69 }
const BaseVectorFunction< P_V, P_M, Q_V, Q_M > & m_meanFunction
Mean function.
template<class P_V = GslVector, class P_M = GslMatrix, class Q_V = GslVector, class Q_M = GslMatrix>
VectorGaussianRandomField& QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::operator= ( const VectorGaussianRandomField< P_V, P_M, Q_V, Q_M > &  rhs)

TODO: Assignment operator; it copies rhs to this.

Todo:
: implement me!
template<class P_V , class P_M , class Q_V , class Q_M >
void QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::sampleFunction ( const std::vector< P_V * > &  fieldPositions,
Q_V &  sampleValues 
)

Function that samples from a Gaussian PDF.

Given the field positions, this method performs a number of tests, calculates the mean vector, the covariance matrix and then it samples from a Gaussian random vector as many positions as required.

Definition at line 80 of file VectorGaussianRandomField.C.

References QUESO::queso_require_equal_to_msg.

81 {
82  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
83  *m_env.subDisplayFile() << "Entering VectorGaussianRandomField<P_V,P_M,Q_V,Q_M>::sampleFunction()"
84  << std::endl;
85  }
86 
87  queso_require_equal_to_msg(( sampleValues.sizeLocal() % fieldPositions.size() ), 0, "input data is not multiple of each other");
88 
89  unsigned int numberOfImageValuesPerIndex = sampleValues.sizeLocal()/fieldPositions.size();
90 
91  queso_require_equal_to_msg(numberOfImageValuesPerIndex, m_imageSetPerIndex.vectorSpace().dimLocal(), "invalid input data dimension");
92 
93  if ((m_savedPositions.size() == 0 ) &&
94  (m_savedRvImageSpace == NULL) &&
95  (m_savedRvLawExpVector == NULL) &&
96  (m_savedRvLawCovMatrix == NULL) &&
97  (m_savedRv == NULL)) {
98  // Ok
99  }
100  else if ((m_savedPositions.size() != 0 ) &&
101  (m_savedRvImageSpace != NULL) &&
102  (m_savedRvLawExpVector != NULL) &&
103  (m_savedRvLawCovMatrix != NULL) &&
104  (m_savedRv != NULL)) {
105  // Ok
106  }
107  else {
108  queso_error_msg("invalid combination of pointer values");
109  }
110 
111  unsigned int numberOfPositions = fieldPositions.size();
112  bool instantiate = true;
113  if (m_savedPositions.size() == numberOfPositions) {
114  bool allPositionsAreEqual = true;
115  for (unsigned int i = 0; i < numberOfPositions; ++i) {
116  queso_require_msg(m_savedPositions[i], "m_savedPositions[i] should not be NULL");
117  if ((m_savedPositions[i]->sizeLocal() == fieldPositions[i]->sizeLocal()) &&
118  (*(m_savedPositions[i]) == *(fieldPositions[i]) )) {
119  // Ok
120  }
121  else {
122  allPositionsAreEqual = false;
123  break;
124  }
125  } // for i
126  instantiate = !allPositionsAreEqual;
127  }
128 
129  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
130  *m_env.subDisplayFile() << "In VectorGaussianRandomField<P_V,P_M,Q_V,Q_M>::sampleFunction()"
131  << ": numberOfPositions = " << numberOfPositions
132  << ", instantiate = " << instantiate
133  << std::endl;
134  }
135 
136  if (instantiate) {
137  delete m_savedRv;
138  delete m_savedRvLawCovMatrix;
139  delete m_savedRvLawExpVector;
140  delete m_savedRvImageSpace;
141  for (unsigned int i = 0; i < m_savedPositions.size(); ++i) {
142  delete m_savedPositions[i];
143  }
144  m_savedPositions.clear();
145 
146  // Set m_savedPositions
147  m_savedPositions.resize(numberOfPositions,NULL);
148  for (unsigned int i = 0; i < m_savedPositions.size(); ++i) {
149  m_savedPositions[i] = new P_V(*(fieldPositions[i]));
150  }
151 
152  // Set m_savedRvImageSpace
153  m_savedRvImageSpace = new VectorSpace<Q_V,Q_M>(m_env, "grf_", numberOfPositions*numberOfImageValuesPerIndex, NULL);
154 
155  // Set m_savedRvLawExpVector
156  Q_V tmpVec(m_imageSetPerIndex.vectorSpace().zeroVector());
158  for (unsigned int i = 0; i < numberOfPositions; ++i) {
159  m_meanFunction.compute(*(fieldPositions[i]),NULL,tmpVec,NULL,NULL,NULL);
160  for (unsigned int j = 0; j < numberOfImageValuesPerIndex; ++j) {
161  (*m_savedRvLawExpVector)[i*numberOfImageValuesPerIndex + j] = tmpVec[j];
162  }
163  }
164 
165  // Set m_savedRvLawCovMatrix
166  Q_M tmpMat(m_imageSetPerIndex.vectorSpace().zeroVector());
168  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
169  *m_env.subDisplayFile() << "In VectorGaussianRandomField<P_V,P_M,Q_V,Q_M>::sampleFunction()"
170  << ": m_savedRvLawCovMatrix order = " << m_savedRvLawCovMatrix->numCols()
171  << ", numberOfPositions = " << numberOfPositions
172  << ", tmpMat order = " << tmpMat.numCols()
173  << ", numberOfImageValuesPerIndex = " << numberOfImageValuesPerIndex
174  << std::endl;
175  }
176  for (unsigned int i = 0; i < numberOfPositions; ++i) {
177  for (unsigned int j = 0; j < numberOfPositions; ++j) {
178  m_covarianceFunction.covMatrix(*(fieldPositions[i]),*(fieldPositions[j]),tmpMat);
179 #if 1
180  Q_M testMat(tmpMat);
181  if (testMat.chol() != 0) {
182  *m_env.subDisplayFile() << "In VectorGaussianRandomField<P_V,P_M,Q_V,Q_M>::sampleFunction()"
183  << ": i = " << i
184  << ", j = " << j
185  << ", *(fieldPositions[i]) = " << *(fieldPositions[i])
186  << ", *(fieldPositions[j]) = " << *(fieldPositions[j])
187  << ", tmpMat = " << tmpMat
188  << ", testMat = " << testMat
189  << ", tmpMat is not positive definite"
190  << std::endl;
191  queso_error_msg("tmpMat is not positive definite");
192  }
193 #endif
194  for (unsigned int k1 = 0; k1 < numberOfImageValuesPerIndex; ++k1) {
195  for (unsigned int k2 = 0; k2 < numberOfImageValuesPerIndex; ++k2) {
196  unsigned int tmpI = i*numberOfImageValuesPerIndex + k1;
197  unsigned int tmpJ = j*numberOfImageValuesPerIndex + k2;
198  (*m_savedRvLawCovMatrix)(tmpI,tmpJ) = tmpMat(k1,k2);
199  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
200  *m_env.subDisplayFile() << "In VectorGaussianRandomField<P_V,P_M,Q_V,Q_M>::sampleFunction()"
201  << ": i = " << i
202  << ", j = " << j
203  << ", k1 = " << k1
204  << ", k2 = " << k2
205  << ", tmpI = " << tmpI
206  << ", tmpJ = " << tmpJ
207  << ", *(fieldPositions[i]) = " << *(fieldPositions[i])
208  << ", *(fieldPositions[j]) = " << *(fieldPositions[j])
209  << ", (*m_savedRvLawCovMatrix)(tmpI,tmpJ) = " << (*m_savedRvLawCovMatrix)(tmpI,tmpJ)
210  << std::endl;
211  }
212  }
213  }
214  }
215  }
216 
217  // Set m_savedRv
218  m_savedRv = new GaussianVectorRV<Q_V,Q_M>("grf_",
222 
223  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
224  *m_env.subDisplayFile() << "In VectorGaussianRandomField<P_V,P_M,Q_V,Q_M>::sampleFunction()"
225  << ": just instantiated Gaussian RV"
226  << "\n *m_savedRvLawExpVector = " << *m_savedRvLawExpVector
227  << "\n *m_savedRvLawCovMatrix = " << *m_savedRvLawCovMatrix
228  << std::endl;
229  for (unsigned int i = 0; i < numberOfPositions; ++i) {
230  *m_env.subDisplayFile() << " *(m_savedPositions[" << i
231  << "]) = " << *(m_savedPositions[i])
232  << std::endl;
233  }
234  }
235  } // if (instantiate)
236 
237  // Generate sample function
238  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
239  *m_env.subDisplayFile() << "In VectorGaussianRandomField<P_V,P_M,Q_V,Q_M>::sampleFunction()"
240  << ": about to realize sample values"
241  << std::endl;
242  }
243  m_savedRv->realizer().realization(sampleValues);
244  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
245  *m_env.subDisplayFile() << "In VectorGaussianRandomField<P_V,P_M,Q_V,Q_M>::sampleFunction()"
246  << ": just realized sample values"
247  << std::endl;
248  }
249 
250  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
251  *m_env.subDisplayFile() << "Leaving VectorGaussianRandomField<P_V,P_M,Q_V,Q_M>::sampleFunction()"
252  << std::endl;
253  }
254 
255  return;
256 }
Q_V * m_savedRvLawExpVector
Vector of the mean value of the RV.
Q_M * m_savedRvLawCovMatrix
Covariance matrix of the RV.
const V & zeroVector() const
Returns a vector filled with zeros.
Definition: VectorSpace.C:174
const VectorSet< Q_V, Q_M > & m_imageSetPerIndex
Image set of the RV, per index.
const BaseVectorFunction< P_V, P_M, Q_V, Q_M > & m_meanFunction
Mean function.
const BaseEnvironment & m_env
Environment.
GaussianVectorRV< Q_V, Q_M > * m_savedRv
My RV.
virtual const VectorSpace< V, M > & vectorSpace() const =0
Vector space to which this set belongs to. See template specialization.
VectorSpace< Q_V, Q_M > * m_savedRvImageSpace
Image set of the RV.
std::vector< P_V * > m_savedPositions
Saved positions.
unsigned int dimLocal() const
Definition: VectorSpace.C:155
MonteCarloSGOptions::MonteCarloSGOptions(const BaseEnvironment &env, const char *prefix, const McOptionsValues &alternativeOptionsValues queso_require_equal_to_msg)(m_env.optionsInputFileName(), std::string(""), std::string("this constructor is incompatible with the existence of an options input file"))
const BaseMatrixCovarianceFunction< P_V, P_M, Q_V, Q_M > & m_covarianceFunction
Covariance function.
unsigned int displayVerbosity() const
Definition: Environment.C:450
const BaseVectorRealizer< V, M > & realizer() const
Finds a realization (sample) of the PDF of this vector RV; access to private attribute m_realizer...
Definition: VectorRV.C:106
std::ofstream * subDisplayFile() const
Access function for m_subDisplayFile (displays file on stream).
Definition: Environment.C:320

Member Data Documentation

template<class P_V = GslVector, class P_M = GslMatrix, class Q_V = GslVector, class Q_M = GslMatrix>
const BaseMatrixCovarianceFunction<P_V,P_M,Q_V,Q_M>& QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::m_covarianceFunction
protected

Covariance function.

Definition at line 114 of file VectorGaussianRandomField.h.

template<class P_V = GslVector, class P_M = GslMatrix, class Q_V = GslVector, class Q_M = GslMatrix>
const BaseEnvironment& QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::m_env
protected

Environment.

Definition at line 99 of file VectorGaussianRandomField.h.

template<class P_V = GslVector, class P_M = GslMatrix, class Q_V = GslVector, class Q_M = GslMatrix>
const VectorSet<Q_V,Q_M>& QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::m_imageSetPerIndex
protected

Image set of the RV, per index.

Definition at line 108 of file VectorGaussianRandomField.h.

template<class P_V = GslVector, class P_M = GslMatrix, class Q_V = GslVector, class Q_M = GslMatrix>
const VectorSet<P_V,P_M>& QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::m_indexSet
protected

Index set.

Definition at line 105 of file VectorGaussianRandomField.h.

template<class P_V = GslVector, class P_M = GslMatrix, class Q_V = GslVector, class Q_M = GslMatrix>
const BaseVectorFunction<P_V,P_M,Q_V,Q_M>& QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::m_meanFunction
protected

Mean function.

Definition at line 111 of file VectorGaussianRandomField.h.

template<class P_V = GslVector, class P_M = GslMatrix, class Q_V = GslVector, class Q_M = GslMatrix>
std::string QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::m_prefix
protected

Prefix.

Definition at line 102 of file VectorGaussianRandomField.h.

template<class P_V = GslVector, class P_M = GslMatrix, class Q_V = GslVector, class Q_M = GslMatrix>
std::vector<P_V*> QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::m_savedPositions
protected
template<class P_V = GslVector, class P_M = GslMatrix, class Q_V = GslVector, class Q_M = GslMatrix>
GaussianVectorRV<Q_V,Q_M>* QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::m_savedRv
protected

My RV.

Definition at line 129 of file VectorGaussianRandomField.h.

template<class P_V = GslVector, class P_M = GslMatrix, class Q_V = GslVector, class Q_M = GslMatrix>
VectorSpace<Q_V,Q_M>* QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::m_savedRvImageSpace
protected

Image set of the RV.

Definition at line 120 of file VectorGaussianRandomField.h.

template<class P_V = GslVector, class P_M = GslMatrix, class Q_V = GslVector, class Q_M = GslMatrix>
Q_M* QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::m_savedRvLawCovMatrix
protected

Covariance matrix of the RV.

Definition at line 126 of file VectorGaussianRandomField.h.

template<class P_V = GslVector, class P_M = GslMatrix, class Q_V = GslVector, class Q_M = GslMatrix>
Q_V* QUESO::VectorGaussianRandomField< P_V, P_M, Q_V, Q_M >::m_savedRvLawExpVector
protected

Vector of the mean value of the RV.

Definition at line 123 of file VectorGaussianRandomField.h.


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

Generated on Tue Jun 5 2018 19:49:38 for queso-0.57.1 by  doxygen 1.8.5