queso-0.57.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
QUESO::ScalarGaussianRandomField< V, M > Class Template Reference

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

#include <ScalarGaussianRandomField.h>

Public Member Functions

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

Protected Member Functions

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

Protected Attributes

const BaseEnvironmentm_env
 Environment. More...
 
std::string m_prefix
 Prefix. More...
 
const VectorSet< V, M > & m_indexSet
 Index set. More...
 
const BaseScalarFunction< V, M > & m_meanFunction
 Mean function. More...
 
const
BaseScalarCovarianceFunction
< V, M > & 
m_covarianceFunction
 Covariance function. More...
 
std::vector< V * > m_savedPositions
 Saved positions. More...
 
VectorSpace< V, M > * m_savedRvImageSpace
 Image set of the RV. More...
 
V * m_savedRvLawExpVector
 Vector of the mean value of the RV. More...
 
M * m_savedRvLawCovMatrix
 Covariance matrix of the RV. More...
 
GaussianVectorRV< V, M > * m_savedRv
 My RV. More...
 

Detailed Description

template<class V = GslVector, class M = GslMatrix>
class QUESO::ScalarGaussianRandomField< V, M >

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

This class implements a scalar Gaussian random field (GRF); i.e. a random field involving Gaussian probability density functions (PDFs) of the variables. A one-dimensional GRF is also called a Gaussian process.

Definition at line 49 of file ScalarGaussianRandomField.h.

Constructor & Destructor Documentation

template<class V , class M >
QUESO::ScalarGaussianRandomField< V, M >::ScalarGaussianRandomField ( const char *  prefix,
const VectorSet< V, M > &  indexSet,
const BaseScalarFunction< V, M > &  meanFunction,
const BaseScalarCovarianceFunction< V, 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 32 of file ScalarGaussianRandomField.C.

References QUESO::ScalarGaussianRandomField< V, M >::m_savedPositions.

37  :
38  m_env (indexSet.env()),
39  m_prefix ((std::string)(prefix)+"grf_"),
43  m_savedRvImageSpace (NULL),
46  m_savedRv (NULL)
47 {
48  m_savedPositions.clear();
49 }
M * m_savedRvLawCovMatrix
Covariance matrix of the RV.
V * m_savedRvLawExpVector
Vector of the mean value of the RV.
GaussianVectorRV< V, M > * m_savedRv
My RV.
VectorSpace< V, M > * m_savedRvImageSpace
Image set of the RV.
const BaseEnvironment & m_env
Environment.
const BaseScalarCovarianceFunction< V, M > & covarianceFunction() const
Covariance function; access to protected attribute m_covarianceFunction.
const BaseScalarCovarianceFunction< V, M > & m_covarianceFunction
Covariance function.
const BaseScalarFunction< V, M > & m_meanFunction
Mean function.
const BaseScalarFunction< V, M > & meanFunction() const
Mean function; access to protected attribute m_meanFunction.
std::vector< V * > m_savedPositions
Saved positions.
const VectorSet< V, M > & indexSet() const
Index set; access to protected attribute m_indexSet.
const VectorSet< V, M > & m_indexSet
Index set.
template<class V = GslVector, class M = GslMatrix>
QUESO::ScalarGaussianRandomField< V, M >::ScalarGaussianRandomField ( const ScalarGaussianRandomField< V, M > &  obj)

TODO: Copy constructor.

Todo:
: implement me!
template<class V , class M >
QUESO::ScalarGaussianRandomField< V, M >::~ScalarGaussianRandomField ( )

Destructor.

Definition at line 52 of file ScalarGaussianRandomField.C.

53 {
54 }

Member Function Documentation

template<class V = GslVector, class M = GslMatrix>
void QUESO::ScalarGaussianRandomField< V, M >::copy ( const ScalarGaussianRandomField< V, M > &  src)
protected

Copy method.

template<class V , class M >
const BaseScalarCovarianceFunction< V, M > & QUESO::ScalarGaussianRandomField< V, M >::covarianceFunction ( ) const

Covariance function; access to protected attribute m_covarianceFunction.

Definition at line 72 of file ScalarGaussianRandomField.C.

73 {
74  return m_covarianceFunction;
75 }
const BaseScalarCovarianceFunction< V, M > & m_covarianceFunction
Covariance function.
template<class V , class M >
const VectorSet< V, M > & QUESO::ScalarGaussianRandomField< V, M >::indexSet ( ) const

Index set; access to protected attribute m_indexSet.

Definition at line 58 of file ScalarGaussianRandomField.C.

59 {
60  return m_indexSet;
61 }
const VectorSet< V, M > & m_indexSet
Index set.
template<class V , class M >
const BaseScalarFunction< V, M > & QUESO::ScalarGaussianRandomField< V, M >::meanFunction ( ) const

Mean function; access to protected attribute m_meanFunction.

Definition at line 65 of file ScalarGaussianRandomField.C.

66 {
67  return m_meanFunction;
68 }
const BaseScalarFunction< V, M > & m_meanFunction
Mean function.
template<class V = GslVector, class M = GslMatrix>
ScalarGaussianRandomField& QUESO::ScalarGaussianRandomField< V, M >::operator= ( const ScalarGaussianRandomField< V, M > &  rhs)

TODO: Assignment operator; it copies rhs to this.

Todo:
: implement me!
template<class V , class M >
void QUESO::ScalarGaussianRandomField< V, M >::sampleFunction ( const std::vector< V * > &  fieldPositions,
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 79 of file ScalarGaussianRandomField.C.

References QUESO::queso_require_equal_to_msg.

80 {
81  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 99)) {
82  *m_env.subDisplayFile() << "Entering ScalarGaussianRandomField<V,M>::sampleFunction()"
83  << std::endl;
84  }
85 
86  queso_require_equal_to_msg(fieldPositions.size(), sampleValues.sizeLocal(), "invalid input data");
87 
88  if ((m_savedPositions.size() == 0 ) &&
89  (m_savedRvImageSpace == NULL) &&
90  (m_savedRvLawExpVector == NULL) &&
91  (m_savedRvLawCovMatrix == NULL) &&
92  (m_savedRv == NULL)) {
93  // Ok
94  }
95  else if ((m_savedPositions.size() != 0 ) &&
96  (m_savedRvImageSpace != NULL) &&
97  (m_savedRvLawExpVector != NULL) &&
98  (m_savedRvLawCovMatrix != NULL) &&
99  (m_savedRv != NULL)) {
100  // Ok
101  }
102  else {
103  queso_error_msg("invalid combination of pointer values");
104  }
105 
106  unsigned int numberOfPositions = fieldPositions.size();
107  bool instantiate = true;
108  if (m_savedPositions.size() == numberOfPositions) {
109  bool allPositionsAreEqual = true;
110  for (unsigned int i = 0; i < numberOfPositions; ++i) {
111  queso_require_msg(m_savedPositions[i], "m_savedPositions[i] should not be NULL");
112  if ((m_savedPositions[i]->sizeLocal() == fieldPositions[i]->sizeLocal()) &&
113  (*(m_savedPositions[i]) == *(fieldPositions[i]) )) {
114  // Ok
115  }
116  else {
117  allPositionsAreEqual = false;
118  break;
119  }
120  } // for i
121  instantiate = !allPositionsAreEqual;
122  }
123 
124  if (instantiate) {
125  delete m_savedRv;
126  delete m_savedRvLawCovMatrix;
127  delete m_savedRvLawExpVector;
128  delete m_savedRvImageSpace;
129  for (unsigned int i = 0; i < m_savedPositions.size(); ++i) {
130  delete m_savedPositions[i];
131  }
132  m_savedPositions.clear();
133 
134  // Set m_savedPositions
135  m_savedPositions.resize(numberOfPositions,NULL);
136  for (unsigned int i = 0; i < m_savedPositions.size(); ++i) {
137  m_savedPositions[i] = new V(*(fieldPositions[i]));
138  }
139 
140  // Set m_savedRvImageSpace
141  m_savedRvImageSpace = new VectorSpace<V,M>(m_env, "grf_", numberOfPositions, NULL);
142 
143  // Set m_savedRvLawExpVector
144  m_savedRvLawExpVector = new V(m_savedRvImageSpace->zeroVector());
145  for (unsigned int i = 0; i < numberOfPositions; ++i) {
146  (*m_savedRvLawExpVector)[i] = m_meanFunction.actualValue(*(fieldPositions[i]),NULL,NULL,NULL,NULL);
147  }
148 
149  // Set m_savedRvLawCovMatrix
150  m_savedRvLawCovMatrix = new M(m_savedRvImageSpace->zeroVector());
151  for (unsigned int i = 0; i < numberOfPositions; ++i) {
152  for (unsigned int j = 0; j < numberOfPositions; ++j) {
153  (*m_savedRvLawCovMatrix)(i,j) = m_covarianceFunction.value(*(fieldPositions[i]),*(fieldPositions[j]));
154  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
155  *m_env.subDisplayFile() << "In ScalarGaussianRandomField<V,M>::sampleFunction()"
156  << ": i = " << i
157  << ", j = " << j
158  << ", *(fieldPositions[i]) = " << *(fieldPositions[i])
159  << ", *(fieldPositions[j]) = " << *(fieldPositions[j])
160  << ", (*m_savedRvLawCovMatrix)(i,j) = " << (*m_savedRvLawCovMatrix)(i,j)
161  << std::endl;
162  }
163  }
164  }
165 
166  // Set m_savedRv
167  m_savedRv = new GaussianVectorRV<V,M>("grf_",
171 
172  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
173  *m_env.subDisplayFile() << "In ScalarGaussianRandomField<V,M>::sampleFunction()"
174  << ": just instantiated Gaussian RV"
175  << "\n *m_savedRvLawExpVector = " << *m_savedRvLawExpVector
176  << "\n *m_savedRvLawCovMatrix = " << *m_savedRvLawCovMatrix
177  << std::endl;
178  for (unsigned int i = 0; i < numberOfPositions; ++i) {
179  *m_env.subDisplayFile() << " *(m_savedPositions[" << i
180  << "]) = " << *(m_savedPositions[i])
181  << std::endl;
182  }
183  }
184  } // if (instantiate)
185 
186  // Generate sample function
187  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 99)) {
188  *m_env.subDisplayFile() << "In ScalarGaussianRandomField<V,M>::sampleFunction()"
189  << ": about to realize sample values"
190  << std::endl;
191  }
192  m_savedRv->realizer().realization(sampleValues);
193  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 99)) {
194  *m_env.subDisplayFile() << "In ScalarGaussianRandomField<V,M>::sampleFunction()"
195  << ": just realized sample values"
196  << std::endl;
197  }
198 
199  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 99)) {
200  *m_env.subDisplayFile() << "Leaving ScalarGaussianRandomField<V,M>::sampleFunction()"
201  << std::endl;
202  }
203 
204  return;
205 }
M * m_savedRvLawCovMatrix
Covariance matrix of the RV.
V * m_savedRvLawExpVector
Vector of the mean value of the RV.
GaussianVectorRV< V, M > * m_savedRv
My RV.
VectorSpace< V, M > * m_savedRvImageSpace
Image set of the RV.
const BaseEnvironment & m_env
Environment.
const BaseScalarCovarianceFunction< V, M > & m_covarianceFunction
Covariance function.
const BaseScalarFunction< V, M > & m_meanFunction
Mean function.
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"))
std::vector< V * > m_savedPositions
Saved positions.
unsigned int displayVerbosity() const
Definition: Environment.C:450
std::ofstream * subDisplayFile() const
Access function for m_subDisplayFile (displays file on stream).
Definition: Environment.C:320

Member Data Documentation

template<class V = GslVector, class M = GslMatrix>
const BaseScalarCovarianceFunction<V,M>& QUESO::ScalarGaussianRandomField< V, M >::m_covarianceFunction
protected

Covariance function.

Definition at line 112 of file ScalarGaussianRandomField.h.

template<class V = GslVector, class M = GslMatrix>
const BaseEnvironment& QUESO::ScalarGaussianRandomField< V, M >::m_env
protected

Environment.

Definition at line 100 of file ScalarGaussianRandomField.h.

template<class V = GslVector, class M = GslMatrix>
const VectorSet<V,M>& QUESO::ScalarGaussianRandomField< V, M >::m_indexSet
protected

Index set.

Definition at line 106 of file ScalarGaussianRandomField.h.

template<class V = GslVector, class M = GslMatrix>
const BaseScalarFunction<V,M>& QUESO::ScalarGaussianRandomField< V, M >::m_meanFunction
protected

Mean function.

Definition at line 109 of file ScalarGaussianRandomField.h.

template<class V = GslVector, class M = GslMatrix>
std::string QUESO::ScalarGaussianRandomField< V, M >::m_prefix
protected

Prefix.

Definition at line 103 of file ScalarGaussianRandomField.h.

template<class V = GslVector, class M = GslMatrix>
std::vector<V*> QUESO::ScalarGaussianRandomField< V, M >::m_savedPositions
protected
template<class V = GslVector, class M = GslMatrix>
GaussianVectorRV<V,M>* QUESO::ScalarGaussianRandomField< V, M >::m_savedRv
protected

My RV.

Definition at line 127 of file ScalarGaussianRandomField.h.

template<class V = GslVector, class M = GslMatrix>
VectorSpace<V,M>* QUESO::ScalarGaussianRandomField< V, M >::m_savedRvImageSpace
protected

Image set of the RV.

Definition at line 118 of file ScalarGaussianRandomField.h.

template<class V = GslVector, class M = GslMatrix>
M* QUESO::ScalarGaussianRandomField< V, M >::m_savedRvLawCovMatrix
protected

Covariance matrix of the RV.

Definition at line 124 of file ScalarGaussianRandomField.h.

template<class V = GslVector, class M = GslMatrix>
V* QUESO::ScalarGaussianRandomField< V, M >::m_savedRvLawExpVector
protected

Vector of the mean value of the RV.

Definition at line 121 of file ScalarGaussianRandomField.h.


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

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