queso-0.53.0
SimulationModel.C
Go to the documentation of this file.
1 //-----------------------------------------------------------------------bl-
2 //--------------------------------------------------------------------------
3 //
4 // QUESO - a library to support the Quantification of Uncertainty
5 // for Estimation, Simulation and Optimization
6 //
7 // Copyright (C) 2008-2015 The PECOS Development Team
8 //
9 // This library is free software; you can redistribute it and/or
10 // modify it under the terms of the Version 2.1 GNU Lesser General
11 // Public License as published by the Free Software Foundation.
12 //
13 // This library is distributed in the hope that it will be useful,
14 // but WITHOUT ANY WARRANTY; without even the implied warranty of
15 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16 // Lesser General Public License for more details.
17 //
18 // You should have received a copy of the GNU Lesser General Public
19 // License along with this library; if not, write to the Free Software
20 // Foundation, Inc. 51 Franklin Street, Fifth Floor,
21 // Boston, MA 02110-1301 USA
22 //
23 //-----------------------------------------------------------------------el-
24 
25 #include <queso/SimulationModel.h>
26 #include <queso/GslVector.h>
27 #include <queso/GslMatrix.h>
28 
29 namespace QUESO {
30 
31 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
32 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::SimulationModel(
33  const char* prefix,
34  const SmOptionsValues* alternativeOptionsValues, // dakota
35  const SimulationStorage<S_V,S_M,P_V,P_M,Q_V,Q_M>& simulStorage)
36  :
37  m_env (simulStorage.env()),
38  m_optionsObj (alternativeOptionsValues),
39  m_paper_p_x (simulStorage.scenarioSpace().dimLocal()),
40  m_paper_p_t (simulStorage.parameterSpace().dimLocal()),
41  m_paper_m (simulStorage.numSimulations()),
42  m_paper_n_eta (simulStorage.outputSpace().dimLocal()),
43  m_p_x_space (m_env, "p_x_", m_paper_p_x, NULL),
44  m_xSeq_original (m_p_x_space, m_paper_m, "simul_xSeq_original"),
45  m_xSeq_original_mins (m_p_x_space.zeroVector()),
46  m_xSeq_original_maxs (m_p_x_space.zeroVector()),
47  m_xSeq_original_ranges (m_p_x_space.zeroVector()),
48  m_xSeq_standard (m_p_x_space, m_paper_m, "simul_xSeq_standard"),
49  m_xSeq_standard_mins (m_p_x_space.zeroVector()),
50  m_xSeq_standard_maxs (m_p_x_space.zeroVector()),
51  m_xSeq_standard_ranges (m_p_x_space.zeroVector()),
52  m_xs_asterisks_standard (simulStorage.xs_asterisks_original().size(),(const S_V*) NULL), // to be deleted on destructor
53  m_p_t_space (m_env, "p_t_", m_paper_p_t, NULL),
54  m_tSeq_original (m_p_t_space, m_paper_m, "simul_tSeq_original"),
55  m_tSeq_mins (m_p_t_space.zeroVector()),
56  m_tSeq_maxs (m_p_t_space.zeroVector()),
57  m_tSeq_ranges (m_p_t_space.zeroVector()),
58  m_tSeq_standard (m_p_t_space, m_paper_m, "simul_tSeq_standard"),
59  m_ts_asterisks_standard (simulStorage.ts_asterisks_original().size(),(const P_V*) NULL), // to be deleted on destructor
60  m_n_eta_space (m_env, "n_eta_", m_paper_n_eta, NULL),
61  m_etaSeq_original (m_n_eta_space, m_paper_m, "simul_etaSeq_original"),
62  m_etaSeq_original_mean (m_n_eta_space.zeroVector()),
63  m_etaSeq_original_std (m_n_eta_space.zeroVector()),
64 #ifdef _GPMSA_CODE_TREATS_SIMULATION_VECTORS_IN_CHUNKS
65  m_etaSeq_chunkMeans (simulStorage.chunkSizes().size(),0.),
66  m_etaSeq_chunkStds (simulStorage.chunkSizes().size(),0.),
67 #else
68  m_etaSeq_allMean (0.),
69  m_etaSeq_allStd (0.),
70 #endif
71  m_etaSeq_transformed (m_n_eta_space, m_paper_m, "simul_etaSeq_transformed"),
72  m_etaSeq_transformed_mean (m_n_eta_space.zeroVector()),
73  m_etaSeq_transformed_std (m_n_eta_space.zeroVector()),
74  m_eta_space (m_env, "m_eta_simul_model", m_paper_m*m_paper_n_eta, NULL),
75  m_etaVec_transformed (m_eta_space.zeroVector()),
76  m_etaMat_transformed (m_env, m_n_eta_space.map(), m_paper_m),
77  m_m_space (m_env, "m_", m_paper_m, NULL),
78  m_m_unitVec (m_env, m_m_space.map(), 1.),
79  m_m_Imat (m_m_unitVec),
80  m_paper_p_eta (0),
81  m_p_eta_space (NULL), // to be deleted on destructor
82  m_Kmat_eta (NULL), // to be deleted on destructor
83  m_kvec_is (m_paper_p_eta, (Q_V*) NULL), // to be deleted on destructor
84  m_Kmat_is (m_paper_p_eta, (Q_M*) NULL), // to be deleted on destructor
85  m_Kmat (NULL) // to be deleted on destructor
86 {
87  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 2)) {
88  *m_env.subDisplayFile() << "Entering SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
89  << ": prefix = " << prefix
90  << ", alternativeOptionsValues = " << alternativeOptionsValues
91  << ", m_env.optionsInputFileName() = " << m_env.optionsInputFileName()
92  << "\n m_paper_p_x = " << m_paper_p_x
93  << "\n m_paper_p_t = " << m_paper_p_t
94  << "\n m_paper_m = " << m_paper_m
95  << "\n m_paper_n_eta =" << m_paper_n_eta
96  << "\n m_paper_p_eta =" << m_paper_p_eta
97  << "\n m_eta_space.dimLocal() = " << m_eta_space.dimLocal()
98  << std::endl;
99  }
100 
101  // If NULL, we create one
102  if (m_optionsObj == NULL) {
103  SmOptionsValues * tempOptions = new SmOptionsValues(&m_env, prefix);
104 
105  // We did this dance because scanOptionsValues is not a const method, but
106  // m_optionsObj is a pointer to const
107  m_optionsObj = tempOptions;
108  }
109 
110  // We'll need to remove this later
111  m_simulationModelOptions = new SimulationModelOptions(m_env, prefix);
112 
113 
114  queso_require_less_equal_msg(m_paper_p_eta, m_paper_m, "'m_paper_p_eta' cannot be bigger than 'm_paper_m'");
115 
116  queso_require_less_equal_msg(m_paper_p_eta, m_paper_n_eta, "'m_paper_p_eta' cannot be bigger than 'm_paper_n_eta'");
117 
118  queso_require_equal_to_msg(m_paper_m, simulStorage.xs_asterisks_original().size(), "invalid m_paper_m (1)");
119 
120  queso_require_equal_to_msg(m_paper_m, simulStorage.ts_asterisks_original().size(), "invalid m_paper_m (2)");
121 
122  //***********************************************************************
123  // Standardize 'xs_asterisks_original'
124  //***********************************************************************
125  // Form 'm_xSeq_original', compute 'm_xSeq_original_mins' and 'm_xSeq_original_ranges', form 'm_xSeq_standard'
126  for (unsigned int i = 0; i < m_paper_m; ++i) {
127  m_xSeq_original.setPositionValues(i,*(simulStorage.xs_asterisks_original()[i]));
128  }
129  m_xSeq_original_mins = m_xSeq_original.subMinPlain();
130  m_xSeq_original_maxs = m_xSeq_original.subMaxPlain();
131  m_xSeq_original_ranges = m_xSeq_original_maxs - m_xSeq_original_mins;
132  S_V tmpXVec(m_p_x_space.zeroVector());
133  for (unsigned int i = 0; i < m_paper_m; ++i) {
134  m_xSeq_original.getPositionValues(i,tmpXVec);
135  tmpXVec -= m_xSeq_original_mins;
136  for (unsigned int j = 0; j < m_paper_p_x; ++j) {
137  if (m_xSeq_original_ranges[j] == 0.) { // CSRI, 2013-Aug-06, with Laura
138  tmpXVec[j] = 0.5;
139  }
140  else {
141  tmpXVec[j] /= m_xSeq_original_ranges[j];
142  }
143  }
144  m_xSeq_standard.setPositionValues(i,tmpXVec);
145  }
146 
147  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
148  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
149  << ": finished setting 'm_xSeq_original'"
150  << ", computing mins of 'm_xSeq_original'"
151  << ", computing ranges of 'm_xSeq_original'"
152  << ", computing 'm_xSeq_standard'"
153  << std::endl;
154  }
155 
156  for (unsigned int i = 0; i < m_paper_m; ++i) {
157  m_xSeq_standard.getPositionValues(i,tmpXVec);
158  m_xs_asterisks_standard[i] = new S_V(tmpXVec); // to be deleted on destructor
159  }
160 
161  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 99)) {
162  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()";
163  for (unsigned int i = 0; i < m_paper_m; ++i) {
164  m_xSeq_original.getPositionValues(i,tmpXVec);
165  *m_env.subDisplayFile() << "\n m_xSeq_original[" << i << "] = " << tmpXVec;
166  }
167  *m_env.subDisplayFile() << "\n m_xSeq_original_mins = " << m_xSeq_original_mins
168  << "\n m_xSeq_original_maxs = " << m_xSeq_original_maxs
169  << "\n m_xSeq_original_ranges = " << m_xSeq_original_ranges
170  << std::endl;
171  for (unsigned int i = 0; i < m_paper_m; ++i) {
172  m_xSeq_standard.getPositionValues(i,tmpXVec);
173  *m_env.subDisplayFile() << "\n m_xSeq_standard[" << i << "] = " << tmpXVec;
174  }
175  m_xSeq_standard_mins = m_xSeq_standard.subMinPlain();
176  m_xSeq_standard_maxs = m_xSeq_standard.subMaxPlain();
177  m_xSeq_standard_ranges = m_xSeq_standard_maxs - m_xSeq_standard_mins;
178  *m_env.subDisplayFile() << "\n m_xSeq_standard_mins = " << m_xSeq_standard_mins
179  << "\n m_xSeq_standard_maxs = " << m_xSeq_standard_maxs
180  << "\n m_xSeq_standard_ranges = " << m_xSeq_standard_ranges
181  << std::endl;
182  }
183 
184  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
185  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
186  << ": finished forming 'm_xs_asterisks_standard'"
187  << std::endl;
188  }
189 
190  //***********************************************************************
191  // Standardize 'ts_asterisks_original'
192  //***********************************************************************
193  // Form 'm_tSeq_original', compute 'm_tSeq_mins' and 'm_tSeq_ranges', form 'm_tSeq_standard'
194  for (unsigned int i = 0; i < m_paper_m; ++i) {
195  m_tSeq_original.setPositionValues(i,*(simulStorage.ts_asterisks_original()[i]));
196  }
197  m_tSeq_mins = m_tSeq_original.subMinPlain();
198  m_tSeq_maxs = m_tSeq_original.subMaxPlain();
199  m_tSeq_ranges = m_tSeq_maxs - m_tSeq_mins;
200  S_V tmpTVec(m_p_t_space.zeroVector());
201  for (unsigned int i = 0; i < m_paper_m; ++i) {
202  m_tSeq_original.getPositionValues(i,tmpTVec);
203  tmpTVec -= m_tSeq_mins;
204  for (unsigned int j = 0; j < m_paper_p_t; ++j) {
205  if (m_tSeq_ranges[j] == 0.) { // CSRI, 2013-Aug-06, with Laura
206  tmpTVec[j] = 0.5;
207  }
208  else {
209  tmpTVec[j] /= m_tSeq_ranges[j];
210  }
211  }
212  m_tSeq_standard.setPositionValues(i,tmpTVec);
213  }
214 
215  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
216  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
217  << ": finished setting 'm_tSeq_original'"
218  << ", computing mins of 'm_tSeq_original'"
219  << ", computing ranges of 'm_tSeq_original'"
220  << ", computing 'm_tSeq_standard'"
221  << std::endl;
222  }
223 
224  for (unsigned int i = 0; i < m_paper_m; ++i) {
225  m_tSeq_standard.getPositionValues(i,tmpTVec);
226  m_ts_asterisks_standard[i] = new S_V(tmpTVec); // to be deleted on destructor
227  }
228 
229  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 99)) {
230  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor():";
231  for (unsigned int i = 0; i < m_paper_m; ++i) {
232  m_tSeq_original.getPositionValues(i,tmpTVec);
233  *m_env.subDisplayFile() << "\n m_tSeq_original[" << i << "] = " << tmpTVec;
234  }
235  *m_env.subDisplayFile() << "\n m_tSeq_original_mins = " << m_tSeq_mins
236  << "\n m_tSeq_original_maxs = " << m_tSeq_maxs
237  << "\n m_tSeq_original_ranges = " << m_tSeq_ranges
238  << std::endl;
239  for (unsigned int i = 0; i < m_paper_m; ++i) {
240  m_tSeq_standard.getPositionValues(i,tmpTVec);
241  *m_env.subDisplayFile() << "\n m_tSeq_standard[" << i << "] = " << tmpTVec;
242  }
243  m_tSeq_mins = m_tSeq_standard.subMinPlain();
244  m_tSeq_maxs = m_tSeq_standard.subMaxPlain();
245  m_tSeq_ranges = m_tSeq_maxs - m_tSeq_mins;
246  *m_env.subDisplayFile() << "\n m_tSeq_standard_mins = " << m_tSeq_mins
247  << "\n m_tSeq_standard_maxs = " << m_tSeq_maxs
248  << "\n m_tSeq_standard_ranges = " << m_tSeq_ranges
249  << std::endl;
250  }
251 
252  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
253  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
254  << ": finished forming 'm_ts_asterisks_standard'"
255  << std::endl;
256  }
257 
258  //***********************************************************************
259  // Form 'etaVec_transformed' and 'etaMat_transformed' matrix
260  //***********************************************************************
261  // Form 'm_etaSeq_original', compute 'm_etaSeq_original_mean' and 'm_etaSeq_original_std', form 'm_etaSeq_transformed'
262  for (unsigned int i = 0; i < m_paper_m; ++i) {
263  m_etaSeq_original.setPositionValues(i,simulStorage.outputVec_original(i));
264  }
265  m_etaSeq_original_mean = m_etaSeq_original.subMeanPlain();
266  m_etaSeq_original.subSampleStd(0,m_paper_m,m_etaSeq_original_mean,m_etaSeq_original_std);
267  Q_V tmpEtaVec(m_n_eta_space.zeroVector());
268  for (unsigned int i = 0; i < m_paper_m; ++i) {
269  m_etaSeq_original.getPositionValues(i,tmpEtaVec);
270  tmpEtaVec -= m_etaSeq_original_mean;
271  //for (unsigned int j = 0; j < m_paper_n_eta; ++j) {
272  // tmpEtaVec[j] /= m_etaSeq_original_std[j];
273  //}
274  m_etaSeq_transformed.setPositionValues(i,tmpEtaVec);
275  }
276 
277 #ifdef _GPMSA_CODE_TREATS_SIMULATION_VECTORS_IN_CHUNKS
278  unsigned int cumulativeSize = 0;
279  for (unsigned int chunkId = 0; chunkId < m_etaSeq_chunkStds.size(); ++chunkId) {
280  m_etaSeq_chunkMeans[chunkId] = 0.;
281  for (unsigned int i = 0; i < m_paper_m; ++i) {
282  m_etaSeq_transformed.getPositionValues(i,tmpEtaVec);
283  for (unsigned int j = cumulativeSize; j < (cumulativeSize + simulStorage.chunkSizes()[chunkId]); ++j) {
284  m_etaSeq_chunkMeans[chunkId] += tmpEtaVec[j];
285  }
286  }
287  m_etaSeq_chunkMeans[chunkId] /= ((double) (m_paper_m * simulStorage.chunkSizes()[chunkId]));
288  cumulativeSize += simulStorage.chunkSizes()[chunkId];
289  }
290  queso_require_equal_to_msg(cumulativeSize, m_paper_n_eta, "inconsistency in 'cumulativeSize' (1)");
291 
292  cumulativeSize = 0;
293  for (unsigned int chunkId = 0; chunkId < m_etaSeq_chunkStds.size(); ++chunkId) {
294  m_etaSeq_chunkStds[chunkId] = 0.;
295  for (unsigned int i = 0; i < m_paper_m; ++i) {
296  m_etaSeq_transformed.getPositionValues(i,tmpEtaVec);
297  for (unsigned int j = cumulativeSize; j < (cumulativeSize + simulStorage.chunkSizes()[chunkId]); ++j) {
298  double diff = tmpEtaVec[j] - m_etaSeq_chunkMeans[chunkId];
299  m_etaSeq_chunkStds[chunkId] += diff * diff;
300  }
301  }
302  m_etaSeq_chunkStds[chunkId] = sqrt( m_etaSeq_chunkStds[chunkId] / ((double) (m_paper_m * simulStorage.chunkSizes()[chunkId] - 1)) );
303  cumulativeSize += simulStorage.chunkSizes()[chunkId];
304  }
305  queso_require_equal_to_msg(cumulativeSize, m_paper_n_eta, "inconsistency in 'cumulativeSize' (2)");
306 
307  for (unsigned int i = 0; i < m_paper_m; ++i) {
308  m_etaSeq_transformed.getPositionValues(i,tmpEtaVec);
309  cumulativeSize = 0;
310  for (unsigned int chunkId = 0; chunkId < m_etaSeq_chunkStds.size(); ++chunkId) {
311  for (unsigned int j = cumulativeSize; j < (cumulativeSize + simulStorage.chunkSizes()[chunkId]); ++j) {
312  tmpEtaVec[j] /= m_etaSeq_chunkStds[chunkId];
313  }
314  cumulativeSize += simulStorage.chunkSizes()[chunkId];
315  }
316  queso_require_equal_to_msg(cumulativeSize, m_paper_n_eta, "inconsistency in 'cumulativeSize' (3)");
317  m_etaSeq_transformed.setPositionValues(i,tmpEtaVec);
318  }
319 #else
320  m_etaSeq_allMean = 0.;
321  for (unsigned int i = 0; i < m_paper_m; ++i) {
322  m_etaSeq_transformed.getPositionValues(i,tmpEtaVec);
323  for (unsigned int j = 0; j < m_paper_n_eta; ++j) {
324  m_etaSeq_allMean += tmpEtaVec[j];
325  }
326  }
327  m_etaSeq_allMean /= ((double) (m_paper_m * m_paper_n_eta));
328 
329  m_etaSeq_allStd = 0.;
330  for (unsigned int i = 0; i < m_paper_m; ++i) {
331  m_etaSeq_transformed.getPositionValues(i,tmpEtaVec);
332  for (unsigned int j = 0; j < m_paper_n_eta; ++j) {
333  double diff = tmpEtaVec[j] - m_etaSeq_allMean;
334  m_etaSeq_allStd += diff * diff;
335  }
336  }
337  m_etaSeq_allStd = sqrt( m_etaSeq_allStd / ((double) (m_paper_m * m_paper_n_eta - 1)) );
338 
339  for (unsigned int i = 0; i < m_paper_m; ++i) {
340  m_etaSeq_transformed.getPositionValues(i,tmpEtaVec);
341  tmpEtaVec /= m_etaSeq_allStd;
342  m_etaSeq_transformed.setPositionValues(i,tmpEtaVec);
343  }
344 #endif
345 
346  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
347  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
348  << ": finished setting 'm_etaSeq_original'"
349  << ", computing means of 'm_etaSeq_original'"
350  << ", computing sample stds of 'm_etaSeq_original'"
351 #ifdef _GPMSA_CODE_TREATS_SIMULATION_VECTORS_IN_CHUNKS
352  << ", computing 'm_etaSeq_chunkMeans'"
353  << ", computing 'm_etaSeq_chunkStds'"
354 #else
355  << ", computing 'm_etaSeq_allmean' = " << m_etaSeq_allMean
356  << ", computing 'm_etaSeq_allStd' = " << m_etaSeq_allStd
357 #endif
358  << ", computing 'm_etaSeq_transformed'"
359  << std::endl;
360  }
361 
362 #ifdef _GPMSA_CODE_TREATS_SIMULATION_VECTORS_IN_CHUNKS
363  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
364  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
365  << ": 'm_etaSeq_chunkMeans' =";
366  for (unsigned int chunkId = 0; chunkId < m_etaSeq_chunkStds.size(); ++chunkId) {
367  *m_env.subDisplayFile() << " " << m_etaSeq_chunkMeans[chunkId];
368  }
369  *m_env.subDisplayFile() << std::endl
370  << ": 'm_etaSeq_chunkStds' =";
371  for (unsigned int chunkId = 0; chunkId < m_etaSeq_chunkStds.size(); ++chunkId) {
372  *m_env.subDisplayFile() << " " << m_etaSeq_chunkStds[chunkId];
373  }
374  *m_env.subDisplayFile() << std::endl;
375  }
376 #endif
377 
378  // Form 'etaVec_transformed' and 'etaMat_transformed' matrix
379  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
380  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
381  << ": before forming 'm_etaVec_transformed' and 'm_etaMat_transformed"
382  << "\n m_paper_p_eta = " << m_paper_p_eta
383  << "\n m_paper_m = " << m_paper_m
384  << "\n m_etaSeq_original.subSequenceSize() = " << m_etaSeq_original.subSequenceSize()
385  << "\n m_etaSeq_transformed.subSequenceSize() = " << m_etaSeq_transformed.subSequenceSize()
386  << "\n m_etaMat_transformed.numRowsLocal() = " << m_etaMat_transformed.numRowsLocal()
387  << "\n m_etaMat_transformed.numCols() = " << m_etaMat_transformed.numCols()
388  << std::endl;
389  }
390 
391  unsigned int cumulativeEtaVecPosition = 0;
392  for (unsigned int i = 0; i < m_paper_m; ++i) {
393  m_etaSeq_transformed.getPositionValues(i,tmpEtaVec);
394  m_etaVec_transformed.cwSet(cumulativeEtaVecPosition,tmpEtaVec);
395  cumulativeEtaVecPosition += tmpEtaVec.sizeLocal();
396  m_etaMat_transformed.setColumn(i,tmpEtaVec);
397  }
398  queso_require_equal_to_msg(cumulativeEtaVecPosition, m_etaVec_transformed.sizeLocal(), "inconsistency in the assemble of 'm_etaVec_transformed'");
399 
400  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 99)) {
401  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()";
402  for (unsigned int i = 0; i < m_paper_m; ++i) {
403  m_etaSeq_original.getPositionValues(i,tmpEtaVec);
404  *m_env.subDisplayFile() << "\n m_etaSeq_original[" << i << "] = " << tmpEtaVec;
405  }
406  *m_env.subDisplayFile() << "\n m_etaSeq_original_mean = " << m_etaSeq_original_mean
407  << "\n m_etaSeq_original_std = " << m_etaSeq_original_std
408  << std::endl;
409  for (unsigned int i = 0; i < m_paper_m; ++i) {
410  m_etaSeq_transformed.getPositionValues(i,tmpEtaVec);
411  *m_env.subDisplayFile() << "\n m_etaSeq_transformed[" << i << "] = " << tmpEtaVec;
412  }
413  m_etaSeq_transformed_mean = m_etaSeq_transformed.subMeanPlain();
414  m_etaSeq_transformed.subSampleStd(0,m_paper_m,m_etaSeq_transformed_mean,m_etaSeq_transformed_std);
415  *m_env.subDisplayFile() << "\n m_etaSeq_transformed_mean = " << m_etaSeq_transformed_mean
416  << "\n m_etaSeq_transformed_std = " << m_etaSeq_transformed_std
417  << std::endl;
418  }
419 
420  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
421  m_etaMat_transformed.setPrintHorizontally(false);
422  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
423  << ": finished forming 'm_etaVec_transformed' and 'm_etaMat_transformed'"
424  << "\n (key-debug) m_etaMat_transformed =\n" << m_etaMat_transformed
425  << std::endl;
426  }
427 
428  //***********************************************************************
429  // Compute SVD of 'etaMat_transformed' matrix
430  //***********************************************************************
431  if (m_etaMat_transformed.numRowsGlobal() >= m_etaMat_transformed.numCols()) {
432  //***********************************************************************
433  // Case 'etaMat.numRows >= etaMat.numCols': perform svd on 'etaMat_transformed'
434  //***********************************************************************
435  Q_M svdU_mat (m_env, m_n_eta_space.map(), m_paper_m);
436  Q_V svdS_vec (m_m_space.zeroVector());
437  Q_M svdVt_mat(m_m_space.zeroVector());
438  //std::cout << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
439  // << ": before m_etaMat_transformed.svd()"
440  // << std::endl;
441  int iRC = m_etaMat_transformed.svd(svdU_mat, svdS_vec, svdVt_mat);
442  //std::cout << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
443  // << ": after m_etaMat_transformed.svd(), iRC = " << iRC
444  // << std::endl;
445  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
446  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
447  << "\n iRC from m_etaMat_transformed.svd(1) = " << iRC
448  << std::endl;
449  }
450  queso_require_equal_to_msg(iRC, 0, "svd(1) failed");
451 
452  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
453  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
454  << ": finished performing svd(1) on 'm_etaMat_transformed'"
455  << "\n svdU_mat =\n" << svdU_mat
456  << "\n svdS_vec =\n" << svdS_vec
457  << "\n svdVt_mat =\n" << svdVt_mat
458  << std::endl;
459  }
460 
461  //***********************************************************************
462  // Case 'etaMat.numRows >= etaMat.numCols': determine 'm_paper_p_eta'
463  // [U,S,V]=svd(ysimStd,0);
464  // lam=diag(S).^2/sum(diag(S).^2);
465  // lam=cumsum(lam);
466  // pu=sum(lam<pcpct)+1;
467  // Ksim=U(:,1:pu)*S(1:pu,1:pu)./sqrt(m);
468  //***********************************************************************
469  m_paper_p_eta = computePEta(svdS_vec);
470  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
471  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
472  << ", m_optionsObj->m_cdfThresholdForPEta = " << m_optionsObj->m_cdfThresholdForPEta
473  << ", m_optionsObj->m_zeroRelativeSingularValue = " << m_optionsObj->m_zeroRelativeSingularValue
474  << ": m_paper_p_eta(1) = " << m_paper_p_eta
475  << std::endl;
476  }
477 
478  m_p_eta_space = new VectorSpace<Q_V,Q_M>(m_env, "p_eta_", m_paper_p_eta, NULL); // to be deleted on destructor
479  m_Kmat_eta = new Q_M (m_env, m_n_eta_space.map(), m_paper_p_eta); // to be deleted on destructor
480 
481  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 2)) {
482  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
483  << ": key-debug"
484  << ", m_Kmat_eta just created (situation 1; not yet populated)"
485  << ", numRowsLocal = " << m_Kmat_eta->numRowsLocal()
486  << ", numCols = " << m_Kmat_eta->numCols()
487  << std::endl;
488  }
489 
490  //***********************************************************************
491  // Case 'etaMat.numRows >= etaMat.numCols': form 'Kmat_eta' matrix
492  //***********************************************************************
493  Q_M matU_tmp(m_env, m_n_eta_space.map(), m_paper_p_eta);
494  Q_V vecU_tmp(m_n_eta_space.zeroVector());
495  for (unsigned int i = 0; i < m_paper_p_eta; ++i) {
496  svdU_mat.getColumn(i, vecU_tmp);
497  matU_tmp.setColumn(i, vecU_tmp);
498  }
499 
500  Q_M matS_tmp(m_p_eta_space->zeroVector());
501  for (unsigned int i = 0; i < m_paper_p_eta; ++i) {
502  matS_tmp(i,i) = svdS_vec[i];
503  }
504 
505  (*m_Kmat_eta) = matU_tmp * matS_tmp;
506 
507  if (m_env.identifyingString() == "CASLexample") {
508  // IMPORTANT: temporary, just to match with CASL example of the Matlab version of GPMSA
509  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
510  *m_env.subDisplayFile() << "IMPORTANT In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
511  << ": multiplying first and third columns of 'm_Kmat_eta' by '-1' in order to match the results of Matlab GPMSA"
512  << std::endl;
513  }
514  for (unsigned int i = 0; i < m_Kmat_eta->numRowsLocal(); ++i) {
515  (*m_Kmat_eta)(i,0) *= -1.;
516  }
517  for (unsigned int i = 0; i < m_Kmat_eta->numRowsLocal(); ++i) {
518  (*m_Kmat_eta)(i,2) *= -1.;
519  }
520  }
521  }
522  else {
523  //***********************************************************************
524  // Case 'etaMat.numRows < etaMat.numCols': perform svd on 'etaMat_transformed'
525  //***********************************************************************
526  Q_M svdU_mat (m_env, m_m_space.map(), m_paper_n_eta);
527  Q_V svdS_vec (m_n_eta_space.zeroVector());
528  Q_M svdVt_mat(m_n_eta_space.zeroVector());
529  Q_M etaMat_transfored_transpose(svdU_mat);
530  etaMat_transfored_transpose.fillWithTranspose(0,0,m_etaMat_transformed,true,true);
531  int iRC = etaMat_transfored_transpose.svd(svdU_mat, svdS_vec, svdVt_mat);
532 
533  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
534  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
535  << "\n iRC from m_etaMat_transformed.svd(2) = " << iRC
536  << std::endl;
537  }
538  queso_require_equal_to_msg(iRC, 0, "svd(2) failed");
539 
540  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
541  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
542  << ": finished performing svd(2) on 'm_etaMat_transformed'"
543  << "\n svdU_mat =\n" << svdU_mat
544  << "\n svdS_vec =\n" << svdS_vec
545  << "\n svdVt_mat =\n" << svdVt_mat
546  << std::endl;
547  }
548 
549  //***********************************************************************
550  // Case 'etaMat.numRows < etaMat.numCols': determine 'm_paper_p_eta'
551  // [U,S,V]=svd(ysimStd,0);
552  // lam=diag(S).^2/sum(diag(S).^2);
553  // lam=cumsum(lam);
554  // pu=sum(lam<pcpct)+1;
555  // Ksim=U(:,1:pu)*S(1:pu,1:pu)./sqrt(m);
556  //***********************************************************************
557  m_paper_p_eta = computePEta(svdS_vec);
558  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
559  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
560  << ", m_optionsObj->m_cdfThresholdForPEta = " << m_optionsObj->m_cdfThresholdForPEta
561  << ", m_optionsObj->m_zeroRelativeSingularValue = " << m_optionsObj->m_zeroRelativeSingularValue
562  << ": m_paper_p_eta(2) = " << m_paper_p_eta
563  << std::endl;
564  }
565 
566  m_p_eta_space = new VectorSpace<Q_V,Q_M>(m_env, "p_eta_", m_paper_p_eta, NULL); // to be deleted on destructor
567  m_Kmat_eta = new Q_M (m_env, m_n_eta_space.map(), m_paper_p_eta); // to be deleted on destructor
568 
569  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 2)) {
570  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
571  << ": key-debug"
572  << ", m_Kmat_eta just created (situation 2; not yet populated)"
573  << ", numRowsLocal = " << m_Kmat_eta->numRowsLocal()
574  << ", numCols = " << m_Kmat_eta->numCols()
575  << std::endl;
576  }
577 
578  //***********************************************************************
579  // Case 'etaMat.numRows < etaMat.numCols': form 'Kmat_eta' matrix
580  //***********************************************************************
581  Q_M svdV_mat(svdVt_mat.transpose());
582  Q_M matV_tmp(m_env, m_n_eta_space.map(), m_paper_p_eta);
583  Q_V vecV_tmp(m_n_eta_space.zeroVector());
584  for (unsigned int i = 0; i < m_paper_p_eta; ++i) {
585  svdV_mat.getColumn(i, vecV_tmp);
586  matV_tmp.setColumn(i, vecV_tmp);
587  }
588 
589  Q_M matS_tmp(m_p_eta_space->zeroVector());
590  for (unsigned int i = 0; i < m_paper_p_eta; ++i) {
591  matS_tmp(i,i) = svdS_vec[i];
592  }
593 
594  (*m_Kmat_eta) = matV_tmp * matS_tmp;
595 
596  if (m_env.identifyingString() == "towerExampleToMatchGPMSA") {
597  // IMPORTANT: temporary, just to match with tower example of the Matlab version of GPMSA
598  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
599  *m_env.subDisplayFile() << "IMPORTANT In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
600  << ": multiplying first column of 'm_Kmat_eta' by '-1' in order to match the results of Matlab GPMSA"
601  << std::endl;
602  }
603  for (unsigned int i = 0; i < m_Kmat_eta->numRowsLocal(); ++i) {
604  (*m_Kmat_eta)(i,0) *= -1.;
605  }
606  }
607  }
608  (*m_Kmat_eta) /= sqrt(m_paper_m); // Scaling (check)
609 
610  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
611  m_Kmat_eta->setPrintHorizontally(false);
612  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
613  << ": finished forming 'm_Kmat_eta'"
614  << "\n (key-debug) *m_Kmat_eta =\n" << *m_Kmat_eta
615  << std::endl;
616  }
617 
618  std::set<unsigned int> tmpSet;
619  tmpSet.insert(m_env.subId());
620  if (m_optionsObj->m_dataOutputAllowedSet.find(m_env.subId()) != m_optionsObj->m_dataOutputAllowedSet.end()) {
621  m_Kmat_eta->subWriteContents("Kmat_eta",
622  "mat_K_eta",
623  "m",
624  tmpSet);
625  }
626 
627 
628  m_kvec_is.resize(m_paper_p_eta, (Q_V*) NULL); // to be deleted on destructor
629  m_Kmat_is.resize(m_paper_p_eta, (Q_M*) NULL); // to be deleted on destructor
630  m_Kmat = new Q_M(m_env, m_eta_space.map(), m_paper_m*m_paper_p_eta); // to be deleted on destructor
631  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 2)) {
632  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
633  << ": key-debug"
634  << ", m_Kmat just created (not yet populated)"
635  << ", numRowsLocal = " << m_Kmat->numRowsLocal()
636  << ", numCols = " << m_Kmat->numCols()
637  << std::endl;
638  }
639 
640  //***********************************************************************
641  // Form 'kvec_is'
642  //***********************************************************************
643  for (unsigned int i = 0; i < m_paper_p_eta; ++i) {
644  m_kvec_is[i] = new Q_V(m_n_eta_space.zeroVector());
645  m_Kmat_eta->getColumn(i,*(m_kvec_is[i]));
646  }
647 
648  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
649  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
650  << ": finished computing 'm_kvec_is'"
651  << std::endl;
652  }
653 
654  //***********************************************************************
655  // Compute 'Kmat_is' (via tensor product)
656  //***********************************************************************
657  for (unsigned int i = 0; i < m_paper_p_eta; ++i) {
658  m_Kmat_is[i] = new Q_M(m_env, m_eta_space.map(), m_paper_m);
659  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 2)) {
660  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
661  << ": key-debug"
662  << ", before 'm_Kmat_is[i]->fillWithTensorProduct()"
663  << ", m_Kmat_is[" << i << "]->numRowsLocal() = " << m_Kmat_is[i]->numRowsLocal()
664  << ", m_Kmat_is[" << i << "]->numCols() = " << m_Kmat_is[i]->numCols()
665  << ", m_m_Imat.numRowsLocal() = " << m_m_Imat.numRowsLocal()
666  << ", m_m_Imat.numCols() = " << m_m_Imat.numCols()
667  << ", m_kvec_is[" << i << "]->sizeLocal() = " << m_kvec_is[i]->sizeLocal()
668  << std::endl;
669  }
670  m_Kmat_is[i]->fillWithTensorProduct(0,0,m_m_Imat,*(m_kvec_is[i]),true,true);
671  }
672 
673  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
674  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
675  << ": finished computing 'm_Kmat_is'"
676  << std::endl;
677  }
678 
679  //***********************************************************************
680  // Form 'Kmat' matrix
681  //***********************************************************************
682  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 3)) {
683  *m_env.subDisplayFile() << "In SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
684  << ": before 'm_Kmat->fillWithBlocksDiagonal()"
685  << "\n m_Kmat->numRowsLocal() = " << m_Kmat->numRowsLocal()
686  << "\n m_Kmat->numCols() = " << m_Kmat->numCols()
687  << std::endl;
688  }
689  m_Kmat->fillWithBlocksHorizontally(0,0,m_Kmat_is,true,true);
690 
691  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 2)) {
692  *m_env.subDisplayFile() << "Leaving SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::constructor()"
693  << ": prefix = " << m_optionsObj->m_prefix
694  << std::endl;
695  }
696 }
697 
698 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
699 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::~SimulationModel()
700 {
701  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 2)) {
702  *m_env.subDisplayFile() << "Entering SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::destructor()..."
703  << std::endl;
704  }
705 
706  delete m_Kmat; // to be deleted on destructor
707  for (unsigned int i = 0; i < m_paper_p_eta; ++i) {
708  delete m_Kmat_is[i]; // to be deleted on destructor
709  delete m_kvec_is[i]; // to be deleted on destructor
710  }
711  delete m_Kmat_eta; // to be deleted on destructor
712  delete m_p_eta_space; // to be deleted on destructor
713  for (unsigned int i = 0; i < m_paper_m; ++i) {
714  delete m_xs_asterisks_standard[i]; // to be deleted on destructor
715  delete m_ts_asterisks_standard[i]; // to be deleted on destructor
716  }
717  if (m_optionsObj) delete m_optionsObj;
718 
719  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 2)) {
720  *m_env.subDisplayFile() << "Leaving SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::destructor()"
721  << std::endl;
722  }
723 }
724 
725 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
726 unsigned int
727 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::computePEta(const Q_V& svdS_vec)
728 {
729  unsigned int result = 0;
730  if (m_optionsObj->m_cdfThresholdForPEta <= 0.) {
731  result = m_optionsObj->m_p_eta;
732  }
733  else if (m_optionsObj->m_cdfThresholdForPEta >= 1.) {
734  for (unsigned int i = 0; i < svdS_vec.sizeLocal(); ++i) {
735  if ((svdS_vec[i]/svdS_vec[0]) > m_optionsObj->m_zeroRelativeSingularValue) {
736  result++;
737  }
738  }
739  }
740  else {
741  Q_V auxVec(svdS_vec);
742  double auxSum = 0.;
743  for (unsigned int i = 0; i < auxVec.sizeLocal(); ++i) {
744  double auxTerm = auxVec[i];
745  auxVec[i] = auxTerm * auxTerm;
746  auxSum += auxVec[i];
747  }
748  for (unsigned int i = 0; i < auxVec.sizeLocal(); ++i) {
749  auxVec[i] /= auxSum;
750  }
751  Q_V auxCumSum(m_m_space.zeroVector());
752  auxCumSum[0] = auxVec[0];
753  for (unsigned int i = 1; i < auxVec.sizeLocal(); ++i) {
754  auxCumSum[i] = auxCumSum[i-1] + auxVec[i];
755  }
756  for (unsigned int i = 0; i < auxVec.sizeLocal(); ++i) {
757  if (auxCumSum[i] < m_optionsObj->m_cdfThresholdForPEta) {
758  result++;
759  }
760  }
761  result += 1;
762  }
763 
764  return result;
765 }
766 
767 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
768 unsigned int
769 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::numBasis() const
770 {
771  return m_paper_p_eta;
772 }
773 
774 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
775 const std::vector<const S_V* >&
776 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::xs_asterisks_standard() const
777 {
778  return m_xs_asterisks_standard;
779 }
780 
781 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
782 const S_V&
783 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::xSeq_original_mins() const
784 {
785  return m_xSeq_original_mins;
786 }
787 
788 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
789 const S_V&
790 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::xSeq_original_ranges() const
791 {
792  return m_xSeq_original_ranges;
793 }
794 
795 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
796 const std::vector<const P_V* >&
797 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::ts_asterisks_standard() const
798 {
799  return m_ts_asterisks_standard;
800 }
801 
802 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
803 const Q_V&
804 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::etaSeq_original_mean() const
805 {
806  return m_etaSeq_original_mean;
807 }
808 
809 #ifdef _GPMSA_CODE_TREATS_SIMULATION_VECTORS_IN_CHUNKS
810 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
811 double
812 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::etaSeq_chunkStd(unsigned int chunkId) const
813 {
814  queso_require_less_equal_msg(chunkId, m_etaSeq_chunkStds.size(), "'chunkId' is too big");
815 
816  return m_etaSeq_chunkStds[chunkId];
817 }
818 #else
819 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
820 double
821 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::etaSeq_allStd() const
822 {
823  return m_etaSeq_allStd;
824 }
825 #endif
826 
827 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
828 const Q_V&
829 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::etaVec_transformed(const std::string& debugString) const
830 {
831  //if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 99)) {
832  // *m_env.subDisplayFile() << "Entering SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::etaVec_transformed()"
833  // << ": debugString = " << debugString
834  // << ", m_etaVec_transformed = " << m_etaVec_transformed
835  // << std::endl;
836  //}
837 
838  return m_etaVec_transformed;
839 }
840 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
841 const Q_V&
842 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::basisVec(unsigned int basisId) const
843 {
844  return *(m_kvec_is[basisId]);
845 }
846 
847 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
848 const Q_M&
849 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::Kmat_eta() const
850 {
851  return *m_Kmat_eta;
852 }
853 
854 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
855 const Q_M&
856 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::Kmat() const
857 {
858  return *m_Kmat;
859 }
860 
861 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
862 const SimulationModelOptions&
863 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::optionsObj() const
864 {
865  // This method returns the old SimulationModelOptions object. We're using
866  // the new SmOptionsValues so this method is deprecated.
867  queso_deprecated();
868  return *m_simulationModelOptions;
869 }
870 
871 template<class S_V,class S_M,class P_V,class P_M,class Q_V,class Q_M>
872 void
873 SimulationModel<S_V,S_M,P_V,P_M,Q_V,Q_M>::print(std::ostream& os) const
874 {
875  return;
876 }
877 
878 } // End namespace QUESO
879 
880 template class QUESO::SimulationModel<QUESO::GslVector, QUESO::GslMatrix, QUESO::GslVector, QUESO::GslMatrix, QUESO::GslVector, QUESO::GslMatrix>;
QUESO::BaseEnvironment::displayVerbosity
unsigned int displayVerbosity() const
Definition: Environment.C:396
QUESO::VectorSpace::dimLocal
unsigned int dimLocal() const
Definition: VectorSpace.C:170
QUESO::SimulationModel::xSeq_original_ranges
const S_V & xSeq_original_ranges() const
Definition: SimulationModel.C:790
QUESO::SimulationModel::computePEta
unsigned int computePEta(const Q_V &svdS_vec)
Definition: SimulationModel.C:727
QUESO::SimulationModel::m_etaSeq_original
SequenceOfVectors< Q_V, Q_M > m_etaSeq_original
Definition: SimulationModel.h:103
QUESO::SimulationModel::optionsObj
const SimulationModelOptions & optionsObj() const
Definition: SimulationModel.C:863
QUESO::SimulationModel::m_m_space
VectorSpace< Q_V, Q_M > m_m_space
Definition: SimulationModel.h:120
QUESO::SimulationModel::Kmat_eta
const Q_M & Kmat_eta() const
Definition: SimulationModel.C:849
QUESO::SimulationModel::m_xSeq_standard
SequenceOfVectors< S_V, S_M > m_xSeq_standard
Definition: SimulationModel.h:88
QUESO::VectorSpace::map
const Map & map() const
Map.
Definition: VectorSpace.C:157
QUESO::SimulationModel::m_p_x_space
VectorSpace< S_V, S_M > m_p_x_space
Definition: SimulationModel.h:83
QUESO::SimulationModel::print
void print(std::ostream &os) const
Definition: SimulationModel.C:873
QUESO::SimulationStorage::xs_asterisks_original
const std::vector< const S_V * > & xs_asterisks_original() const
Definition: SimulationStorage.C:134
QUESO::SimulationModel::etaSeq_original_mean
const Q_V & etaSeq_original_mean() const
Definition: SimulationModel.C:804
QUESO::SmOptionsValues::m_dataOutputAllowedSet
std::set< unsigned int > m_dataOutputAllowedSet
Definition: SimulationModelOptions.h:70
QUESO::SimulationModel::m_etaSeq_transformed
SequenceOfVectors< Q_V, Q_M > m_etaSeq_transformed
Definition: SimulationModel.h:113
QUESO::SimulationModel::m_env
const BaseEnvironment & m_env
Definition: SimulationModel.h:74
queso_require_less_equal_msg
#define queso_require_less_equal_msg(expr1, expr2, msg)
Definition: asserts.h:89
QUESO::SimulationModel::m_xs_asterisks_standard
std::vector< const S_V * > m_xs_asterisks_standard
Definition: SimulationModel.h:92
QUESO::SimulationModel::m_eta_space
VectorSpace< Q_V, Q_M > m_eta_space
Definition: SimulationModel.h:116
QUESO::SimulationModel::numBasis
unsigned int numBasis() const
Definition: SimulationModel.C:769
QUESO::BaseVectorSequence::subMaxPlain
const V & subMaxPlain() const
Finds the maximum value of the sub-sequence.
Definition: VectorSequence.C:172
QUESO::VectorSpace::zeroVector
const V & zeroVector() const
Returns a vector filled with zeros.
Definition: VectorSpace.C:189
QUESO::SequenceOfVectors::subSampleStd
void subSampleStd(unsigned int initialPos, unsigned int numPos, const V &meanVec, V &stdVec) const
Finds the sample standard deviation of the sub-sequence, considering numPos positions starting at pos...
Definition: SequenceOfVectors.C:601
queso_require_equal_to_msg
#define queso_require_equal_to_msg(expr1, expr2, msg)
Definition: asserts.h:85
QUESO::SimulationModel::basisVec
const Q_V & basisVec(unsigned int basisId) const
Definition: SimulationModel.C:842
QUESO::SimulationModel::etaSeq_allStd
double etaSeq_allStd() const
Definition: SimulationModel.C:821
QUESO::SimulationModel::Kmat
const Q_M & Kmat() const
Definition: SimulationModel.C:856
QUESO::BaseEnvironment::optionsInputFileName
std::string optionsInputFileName() const
Access to the attribute m_optionsInputFileName, which stores the name of the input file passed by the...
Definition: Environment.C:307
QUESO::SequenceOfVectors::setPositionValues
void setPositionValues(unsigned int posId, const V &vec)
Set the values in vec at position posId of the sequence.
Definition: SequenceOfVectors.C:206
QUESO::SimulationModel::m_n_eta_space
VectorSpace< Q_V, Q_M > m_n_eta_space
Definition: SimulationModel.h:102
QUESO::BaseEnvironment::subDisplayFile
std::ofstream * subDisplayFile() const
Access function for m_subDisplayFile (displays file on stream).
Definition: Environment.C:274
QUESO::BaseVectorSequence::subMinPlain
const V & subMinPlain() const
Finds the minimum value of the sub-sequence.
Definition: VectorSequence.C:146
QUESO::SimulationStorage::ts_asterisks_original
const std::vector< const P_V * > & ts_asterisks_original() const
Definition: SimulationStorage.C:141
QUESO::BaseEnvironment::subId
unsigned int subId() const
Access function to the number of each sub-environment Id: m_subId.
Definition: Environment.C:295
QUESO::SimulationModel::ts_asterisks_standard
const std::vector< const P_V * > & ts_asterisks_standard() const
Definition: SimulationModel.C:797
QUESO::SequenceOfVectors::subSequenceSize
unsigned int subSequenceSize() const
Size of the sub-sequence of vectors.
Definition: SequenceOfVectors.C:103
QUESO::SimulationModel::m_Kmat_is
std::vector< Q_M * > m_Kmat_is
Definition: SimulationModel.h:128
QUESO::SimulationModel::m_p_eta_space
VectorSpace< Q_V, Q_M > * m_p_eta_space
Definition: SimulationModel.h:125
QUESO::SimulationModel::m_optionsObj
const SmOptionsValues * m_optionsObj
Definition: SimulationModel.h:75
QUESO::SimulationModel::m_kvec_is
std::vector< Q_V * > m_kvec_is
Definition: SimulationModel.h:127
QUESO::SimulationModel::etaVec_transformed
const Q_V & etaVec_transformed(const std::string &debugString) const
Definition: SimulationModel.C:829
queso_deprecated
#define queso_deprecated()
Definition: Defines.h:120
QUESO::SimulationModel::xs_asterisks_standard
const std::vector< const S_V * > & xs_asterisks_standard() const
Definition: SimulationModel.C:776
QUESO::SequenceOfVectors::getPositionValues
void getPositionValues(unsigned int posId, V &vec) const
Gets the values of the sequence at position posId and stores them at vec.
Definition: SequenceOfVectors.C:188
QUESO::SimulationModel::m_tSeq_standard
SequenceOfVectors< P_V, P_M > m_tSeq_standard
Definition: SimulationModel.h:99
QUESO::SimulationModel::m_paper_n_eta
unsigned int m_paper_n_eta
Definition: SimulationModel.h:81
QUESO::SimulationModel::SimulationModel
SimulationModel(const char *prefix, const SmOptionsValues *alternativeOptionsValues, const SimulationStorage< S_V, S_M, P_V, P_M, Q_V, Q_M > &simulStorage)
Definition: SimulationModel.C:32
QUESO::SimulationModel::m_simulationModelOptions
SimulationModelOptions * m_simulationModelOptions
Definition: SimulationModel.h:76
QUESO::SimulationModel::m_p_t_space
VectorSpace< P_V, P_M > m_p_t_space
Definition: SimulationModel.h:94
QUESO::SimulationStorage::outputVec_original
const Q_V & outputVec_original(unsigned int simulationId) const
Definition: SimulationStorage.C:191
QUESO::SimulationModel::xSeq_original_mins
const S_V & xSeq_original_mins() const
Definition: SimulationModel.C:783
QUESO::SimulationModel::m_xSeq_original
SequenceOfVectors< S_V, S_M > m_xSeq_original
Definition: SimulationModel.h:84
QUESO::BaseVectorSequence::subMeanPlain
const V & subMeanPlain() const
Finds the mean value of the sub-sequence.
Definition: VectorSequence.C:198
QUESO::SimulationModel::m_ts_asterisks_standard
std::vector< const P_V * > m_ts_asterisks_standard
Definition: SimulationModel.h:100
QUESO::SimulationModel::m_tSeq_original
SequenceOfVectors< P_V, P_M > m_tSeq_original
Definition: SimulationModel.h:95
Generated on Thu Jun 11 2015 13:52:32 for queso-0.53.0 by  doxygen 1.8.5