queso-0.57.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
SequenceOfVectors.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-2017 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/SequenceOfVectors.h>
26 #include <queso/GslVector.h>
27 #include <queso/GslMatrix.h>
28 
29 namespace QUESO {
30 
31 // Default constructor -----------------------------
32 template <class V, class M>
33 SequenceOfVectors<V,M>::SequenceOfVectors(
34  const VectorSpace<V,M>& vectorSpace,
35  unsigned int subSequenceSize,
36  const std::string& name)
37  :
38  BaseVectorSequence<V,M>(vectorSpace,subSequenceSize,name),
39  m_seq (subSequenceSize,NULL)
40 #ifdef UQ_CODE_HAS_MONITORS
41  ,
42  m_subMeanMonitorPosSeq (NULL),
43  m_subMeanVecSeq (NULL),
44  m_subMeanCltStdSeq (NULL),
45  m_subMeanInter0MonitorPosSeq (NULL),
46  m_subMeanInter0Mean (NULL),
47  m_subMeanInter0Clt95 (NULL),
48  m_subMeanInter0Empirical90 (NULL),
49  m_subMeanInter0Min (NULL),
50  m_subMeanInter0Max (NULL),
51  m_unifiedMeanMonitorPosSeq (NULL),
52  m_unifiedMeanVecSeq (NULL),
53  m_unifiedMeanCltStdSeq (NULL)
54 #endif
55 {
56  //if (m_env.subDisplayFile()) {
57  // *m_env.subDisplayFile() << "Entering SequenceOfVectors<V,M>::constructor()"
58  // << std::endl;
59  //}
60 
61  //if (m_env.subDisplayFile()) {
62  // *m_env.subDisplayFile() << "Leaving SequenceOfVectors<V,M>::constructor()"
63  // << std::endl;
64  //}
65 }
66 // Destructor ---------------------------------------
67 template <class V, class M>
68 SequenceOfVectors<V,M>::~SequenceOfVectors()
69 {
70 #ifdef UQ_CODE_HAS_MONITORS
71  if (m_subMeanMonitorPosSeq) delete m_subMeanMonitorPosSeq;
72  if (m_subMeanVecSeq ) delete m_subMeanVecSeq;
73  if (m_subMeanCltStdSeq ) delete m_subMeanCltStdSeq;
74 
75  if (m_subMeanInter0MonitorPosSeq) delete m_subMeanInter0MonitorPosSeq;
76  if (m_subMeanInter0Mean ) delete m_subMeanInter0Mean;
77  if (m_subMeanInter0Clt95 ) delete m_subMeanInter0Clt95;
78  if (m_subMeanInter0Empirical90 ) delete m_subMeanInter0Empirical90;
79  if (m_subMeanInter0Min ) delete m_subMeanInter0Min;
80  if (m_subMeanInter0Max ) delete m_subMeanInter0Max;
81 
82  if (m_unifiedMeanMonitorPosSeq) delete m_unifiedMeanMonitorPosSeq;
83  if (m_unifiedMeanVecSeq ) delete m_unifiedMeanVecSeq;
84  if (m_unifiedMeanCltStdSeq ) delete m_unifiedMeanCltStdSeq;
85 #endif
86 
87  for (unsigned int i = 0; i < (unsigned int) m_seq.size(); ++i) {
88  if (m_seq[i]) delete m_seq[i];
89  }
90 }
91 // Set methods --------------------------------------
92 template <class V, class M>
93 SequenceOfVectors<V,M>&
94 SequenceOfVectors<V,M>::operator= (const SequenceOfVectors<V,M>& rhs)
95 {
96  this->copy(rhs);
97  return *this;
98 }
99 
100 // Sequence methods ---------------------------------
101 template <class V, class M>
102 unsigned int
103 SequenceOfVectors<V,M>::subSequenceSize() const
104 {
105  return m_seq.size();
106 }
107 //---------------------------------------------------
108 template <class V, class M>
109 void
110 SequenceOfVectors<V,M>::resizeSequence(unsigned int newSubSequenceSize)
111 {
112  if (newSubSequenceSize != this->subSequenceSize()) {
113  if (newSubSequenceSize < this->subSequenceSize()) {
114  this->resetValues(newSubSequenceSize,this->subSequenceSize()-newSubSequenceSize);
115  }
116  m_seq.resize(newSubSequenceSize,NULL);
117  std::vector<const V*>(m_seq).swap(m_seq);
118  BaseVectorSequence<V,M>::deleteStoredVectors();
119  }
120 
121  return;
122 }
123 //---------------------------------------------------
124 template <class V, class M>
125 void
126 SequenceOfVectors<V,M>::resetValues(unsigned int initialPos, unsigned int numPos)
127 {
128  bool bRC = ((initialPos < this->subSequenceSize()) &&
129  (0 < numPos ) &&
130  ((initialPos+numPos) <= this->subSequenceSize()));
131  if ((bRC == false) && (m_env.subDisplayFile())) {
132  *m_env.subDisplayFile() << "In SequenceOfVectors<V,M>::resetValues()"
133  << ", initialPos = " << initialPos
134  << ", this->subSequenceSize() = " << this->subSequenceSize()
135  << ", numPos = " << numPos
136  << std::endl;
137  }
138  queso_require_msg(bRC, "invalid input data");
139 
140  for (unsigned int j = 0; j < numPos; ++j) {
141  if (m_seq[initialPos+j] != NULL) {
142  delete m_seq[initialPos+j];
143  m_seq[initialPos+j] = NULL;
144  }
145  }
146 
147  BaseVectorSequence<V,M>::deleteStoredVectors();
148 
149  return;
150 }
151 //---------------------------------------------------
152 template <class V, class M>
153 void
154 SequenceOfVectors<V,M>::erasePositions(unsigned int initialPos, unsigned int numPos)
155 {
156  bool bRC = ((initialPos < this->subSequenceSize()) &&
157  (0 < numPos ) &&
158  ((initialPos+numPos) <= this->subSequenceSize()));
159  queso_require_msg(bRC, "invalid input data");
160 
161  for (unsigned int j = 0; j < numPos; ++j) {
162  if (m_seq[initialPos+j] != NULL) {
163  delete m_seq[initialPos+j];
164  m_seq[initialPos+j] = NULL;
165  }
166  }
167 
168  seqVectorPositionIteratorTypedef posIteratorBegin = m_seq.begin();
169  if (initialPos < this->subSequenceSize()) std::advance(posIteratorBegin,initialPos);
170  else posIteratorBegin = m_seq.end();
171 
172  unsigned int posEnd = initialPos + numPos;
173  seqVectorPositionIteratorTypedef posIteratorEnd = m_seq.begin();
174  if (posEnd < this->subSequenceSize()) std::advance(posIteratorEnd,posEnd);
175  else posIteratorEnd = m_seq.end();
176 
177  unsigned int oldSubSequenceSize = this->subSequenceSize();
178  m_seq.erase(posIteratorBegin,posIteratorEnd);
179  queso_require_equal_to_msg((oldSubSequenceSize - numPos), this->subSequenceSize(), "(oldSubSequenceSize - numPos) != this->subSequenceSize()");
180 
181  BaseVectorSequence<V,M>::deleteStoredVectors();
182 
183  return;
184 }
185 //---------------------------------------------------
186 template <class V, class M>
187 void
188 SequenceOfVectors<V,M>::getPositionValues(unsigned int posId, V& vec) const
189 {
190  queso_require_less_msg(posId, this->subSequenceSize(), "posId > subSequenceSize()");
191 
192  queso_require_msg(m_seq[posId], "posId is NULL");
193 
194  //if (posId == 0) { // mox
195  // std::cout << "In SequenceOfVectors<V,M>::getPositionValues(): m_seq[0] = " << m_seq[0] << ", *(m_seq[0]) = " << *(m_seq[0])
196  // << std::endl;
197  //}
198 
199  vec = *(m_seq[posId]); // *(const_cast<V*>(m_seq[posId])); // prudenci 2010-06-17 mox
200 
201  return;
202 }
203 //---------------------------------------------------
204 template <class V, class M>
205 void
206 SequenceOfVectors<V,M>::setPositionValues(unsigned int posId, const V& vec)
207 {
208  queso_require_less_msg(posId, this->subSequenceSize(), "posId > subSequenceSize()");
209 
210  queso_require_equal_to_msg(vec.sizeLocal(), m_vectorSpace.zeroVector().sizeLocal(), "invalid vec");
211 
212  if (m_seq[posId] != NULL) delete m_seq[posId];
213  m_seq[posId] = new V(vec);
214 
215  //if (posId == 0) { // mox
216  // std::cout << "In SequenceOfVectors<V,M>::setPositionValues(): m_seq[0] = " << m_seq[0] << ", *(m_seq[0]) = " << *(m_seq[0])
217  // << std::endl;
218  //}
219 
220  BaseVectorSequence<V,M>::deleteStoredVectors();
221 
222  return;
223 }
224 //---------------------------------------------------
225 template <class V, class M>
226 void
227 SequenceOfVectors<V,M>::subUniformlySampledCdf(
228  const V& numEvaluationPointsVec,
229  ArrayOfOneDGrids <V,M>& cdfGrids,
230  ArrayOfOneDTables<V,M>& cdfValues) const
231 {
232  V minDomainValues(m_vectorSpace.zeroVector());
233  V maxDomainValues(m_vectorSpace.zeroVector());
234 
235  ScalarSequence<double> data(m_env,0,"");
236 
237  unsigned int numParams = this->vectorSizeLocal();
238  for (unsigned int i = 0; i < numParams; ++i) {
239  this->extractScalarSeq(0, // initialPos
240  1, // spacing
241  subSequenceSize(), // numPos
242  i,
243  data);
244 
245  std::vector<double> aCdf(0);
246  data.subUniformlySampledCdf((unsigned int) numEvaluationPointsVec[i],
247  minDomainValues[i],
248  maxDomainValues[i],
249  aCdf);
250  cdfValues.setOneDTable(i,aCdf);
251  }
252 
253  cdfGrids.setUniformGrids(numEvaluationPointsVec,
254  minDomainValues,
255  maxDomainValues);
256 
257  return;
258 }
259 //---------------------------------------------------
260 template <class V, class M>
261 void
262 SequenceOfVectors<V,M>::unifiedUniformlySampledCdf(
263  const V& numEvaluationPointsVec,
264  ArrayOfOneDGrids <V,M>& unifiedCdfGrids,
265  ArrayOfOneDTables<V,M>& unifiedCdfValues) const
266 {
267  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 10)) {
268  *m_env.subDisplayFile() << "Entering SequenceOfVectors<V,M>::unifiedUniformlySampledCdf()"
269  << std::endl;
270  }
271 
272  V unifiedMinDomainValues(m_vectorSpace.zeroVector());
273  V unifiedMaxDomainValues(m_vectorSpace.zeroVector());
274 
275  ScalarSequence<double> data(m_env,0,"");
276 
277  unsigned int numParams = this->vectorSizeLocal();
278  for (unsigned int i = 0; i < numParams; ++i) {
279  this->extractScalarSeq(0, // initialPos
280  1, // spacing
281  subSequenceSize(), // numPos
282  i,
283  data);
284 
285  std::vector<double> aCdf(0);
286  data.unifiedUniformlySampledCdf(m_vectorSpace.numOfProcsForStorage() == 1,
287  (unsigned int) numEvaluationPointsVec[i],
288  unifiedMinDomainValues[i],
289  unifiedMaxDomainValues[i],
290  aCdf);
291  unifiedCdfValues.setOneDTable(i,aCdf);
292  }
293 
294  unifiedCdfGrids.setUniformGrids(numEvaluationPointsVec,
295  unifiedMinDomainValues,
296  unifiedMaxDomainValues);
297 
298  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 10)) {
299  *m_env.subDisplayFile() << "Leaving SequenceOfVectors<V,M>::unifiedUniformlySampledCdf()"
300  << std::endl;
301  }
302 
303  return;
304 }
305 //---------------------------------------------------
306 template <class V, class M>
307 void
308 SequenceOfVectors<V,M>::subMeanExtra(
309  unsigned int initialPos,
310  unsigned int numPos,
311  V& meanVec) const
312 {
313  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 5)) {
314  *m_env.subDisplayFile() << "Entering SequenceOfVectors<V,M>::subMeanExtra()"
315  << ": initialPos = " << initialPos
316  << ", numPos = " << numPos
317  << ", sub sequence size = " << this->subSequenceSize()
318  << std::endl;
319  }
320 
321  bool bRC = ((initialPos < this->subSequenceSize()) &&
322  (0 < numPos ) &&
323  ((initialPos+numPos) <= this->subSequenceSize()) &&
324  (this->vectorSizeLocal() == meanVec.sizeLocal() ));
325  if ((bRC == false) && (m_env.subDisplayFile())) {
326  *m_env.subDisplayFile() << "In SequenceOfVectors<V,M>::subMeanExtra()"
327  << ", initialPos = " << initialPos
328  << ", this->subSequenceSize() = " << this->subSequenceSize()
329  << ", numPos = " << numPos
330  << ", this->vectorSizeLocal() = " << this->vectorSizeLocal()
331  << ", meanVec.sizeLocal() = " << meanVec.sizeLocal()
332  << std::endl;
333  }
334  queso_require_msg(bRC, "invalid input data");
335 
336  ScalarSequence<double> data(m_env,0,"");
337 
338  unsigned int numParams = this->vectorSizeLocal();
339  for (unsigned int i = 0; i < numParams; ++i) {
340  this->extractScalarSeq(initialPos,
341  1, // spacing
342  numPos,
343  i,
344  data);
345  meanVec[i] = data.subMeanExtra(0,
346  numPos);
347  }
348 
349  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 5)) {
350  *m_env.subDisplayFile() << "Leaving SequenceOfVectors<V,M>::subMeanExtra()"
351  << ": initialPos = " << initialPos
352  << ", numPos = " << numPos
353  << ", sub sequence size = " << this->subSequenceSize()
354  << ", meanVec = " << meanVec
355  << std::endl;
356  }
357 
358  return;
359 }
360 //---------------------------------------------------
361 template <class V, class M>
362 void
363 SequenceOfVectors<V,M>::unifiedMeanExtra(
364  unsigned int initialPos,
365  unsigned int numPos,
366  V& unifiedMeanVec) const
367 {
368  unsigned int tmpUnif = this->unifiedSequenceSize();
369  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 5)) {
370  *m_env.subDisplayFile() << "Entering SequenceOfVectors<V,M>::unifiedMeanExtra()"
371  << ": initialPos = " << initialPos
372  << ", numPos = " << numPos
373  << ", sub sequence size = " << this->subSequenceSize()
374  << ", unified sequence size = " << tmpUnif
375  << std::endl;
376  }
377 
378  bool bRC = ((initialPos < this->subSequenceSize() ) &&
379  (0 < numPos ) &&
380  ((initialPos+numPos) <= this->subSequenceSize() ) &&
381  (this->vectorSizeLocal() == unifiedMeanVec.sizeLocal()));
382  if ((bRC == false) && (m_env.subDisplayFile())) {
383  *m_env.subDisplayFile() << "In SequenceOfVectors<V,M>::unifiedMeanExtra()"
384  << ", initialPos = " << initialPos
385  << ", this->subSequenceSize() = " << this->subSequenceSize()
386  << ", numPos = " << numPos
387  << ", this->vectorSizeLocal() = " << this->vectorSizeLocal()
388  << ", unifiedMeanVec.sizeLocal() = " << unifiedMeanVec.sizeLocal()
389  << std::endl;
390  }
391  queso_require_msg(bRC, "invalid input data");
392 
393  ScalarSequence<double> data(m_env,0,"");
394 
395  unsigned int numParams = this->vectorSizeLocal();
396  for (unsigned int i = 0; i < numParams; ++i) {
397  this->extractScalarSeq(initialPos,
398  1, // spacing
399  numPos,
400  i,
401  data);
402  unifiedMeanVec[i] = data.unifiedMeanExtra(m_vectorSpace.numOfProcsForStorage() == 1,
403  0,
404  numPos);
405  }
406 
407  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 5)) {
408  *m_env.subDisplayFile() << "Leaving SequenceOfVectors<V,M>::unifiedMeanExtra()"
409  << ": initialPos = " << initialPos
410  << ", numPos = " << numPos
411  << ", sub sequence size = " << this->subSequenceSize()
412  << ", unified sequence size = " << tmpUnif
413  << ", unifiedMeanVec = " << unifiedMeanVec
414  << std::endl;
415  }
416 
417  return;
418 }
419 //---------------------------------------------------
420 template <class V, class M>
421 void
422 SequenceOfVectors<V,M>::subMedianExtra(
423  unsigned int initialPos,
424  unsigned int numPos,
425  V& medianVec) const
426 {
427  if (this->subSequenceSize() == 0) return;
428 
429  bool bRC = ((initialPos < this->subSequenceSize()) &&
430  (0 < numPos ) &&
431  ((initialPos+numPos) <= this->subSequenceSize()));
432  if (bRC == false) {
433  std::cerr << "In SequenceOfVectors<V,M>::subMedianExtra()"
434  << ": ERROR at fullRank " << m_env.fullRank()
435  << ", initialPos = " << initialPos
436  << ", numPos = " << numPos
437  << ", this->subSequenceSize() = " << this->subSequenceSize()
438  << std::endl;
439  if (m_env.subDisplayFile()) {
440  *m_env.subDisplayFile() << "In SequenceOfVectors<V,M>::subMedianExtra()"
441  << ": ERROR at fullRank " << m_env.fullRank()
442  << ", initialPos = " << initialPos
443  << ", numPos = " << numPos
444  << ", this->subSequenceSize() = " << this->subSequenceSize()
445  << std::endl;
446  }
447  }
448  queso_require_msg(bRC, "invalid input data");
449 
450  ScalarSequence<double> data(m_env,0,"");
451 
452  unsigned int numParams = this->vectorSizeLocal();
453  for (unsigned int i = 0; i < numParams; ++i) {
454  this->extractScalarSeq(initialPos,
455  1, // spacing
456  numPos,
457  i,
458  data);
459  medianVec[i] = data.subMedianExtra(0,
460  numPos);
461  }
462 
463  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 5)) {
464  *m_env.subDisplayFile() << "Leaving SequenceOfVectors<V,M>::subMedianExtra()"
465  << ": initialPos = " << initialPos
466  << ", numPos = " << numPos
467  << ", sub sequence size = " << this->subSequenceSize()
468  << ", medianVec = " << medianVec
469  << std::endl;
470  }
471 
472  return;
473 }
474 //---------------------------------------------------
475 template <class V, class M>
476 void
477 SequenceOfVectors<V,M>::unifiedMedianExtra(
478  unsigned int initialPos,
479  unsigned int numPos,
480  V& unifiedMedianVec) const
481 {
482  unsigned int tmpUnif = this->unifiedSequenceSize();
483  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 5)) {
484  *m_env.subDisplayFile() << "Entering SequenceOfVectors<V,M>::unifiedMedianExtra()"
485  << ": initialPos = " << initialPos
486  << ", numPos = " << numPos
487  << ", sub sequence size = " << this->subSequenceSize()
488  << ", unified sequence size = " << tmpUnif
489  << std::endl;
490  }
491 
492  bool bRC = ((initialPos < this->subSequenceSize() ) &&
493  (0 < numPos ) &&
494  ((initialPos+numPos) <= this->subSequenceSize() ) &&
495  (this->vectorSizeLocal() == unifiedMedianVec.sizeLocal()));
496  if ((bRC == false) && (m_env.subDisplayFile())) {
497  *m_env.subDisplayFile() << "In SequenceOfVectors<V,M>::unifiedMedianExtra()"
498  << ", initialPos = " << initialPos
499  << ", this->subSequenceSize() = " << this->subSequenceSize()
500  << ", numPos = " << numPos
501  << ", this->vectorSizeLocal() = " << this->vectorSizeLocal()
502  << ", unifiedMedianVec.sizeLocal() = " << unifiedMedianVec.sizeLocal()
503  << std::endl;
504  }
505  queso_require_msg(bRC, "invalid input data");
506 
507  ScalarSequence<double> data(m_env,0,"");
508 
509  unsigned int numParams = this->vectorSizeLocal();
510  for (unsigned int i = 0; i < numParams; ++i) {
511  this->extractScalarSeq(initialPos,
512  1, // spacing
513  numPos,
514  i,
515  data);
516  unifiedMedianVec[i] = data.unifiedMedianExtra(m_vectorSpace.numOfProcsForStorage() == 1,
517  0,
518  numPos);
519  }
520 
521  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 5)) {
522  *m_env.subDisplayFile() << "Leaving SequenceOfVectors<V,M>::unifiedMedianExtra()"
523  << ": initialPos = " << initialPos
524  << ", numPos = " << numPos
525  << ", sub sequence size = " << this->subSequenceSize()
526  << ", unified sequence size = " << tmpUnif
527  << ", unifiedMedianVec = " << unifiedMedianVec
528  << std::endl;
529  }
530 
531  return;
532 }
533 //---------------------------------------------------
534 template <class V, class M>
535 void
536 SequenceOfVectors<V,M>::subSampleVarianceExtra(
537  unsigned int initialPos,
538  unsigned int numPos,
539  const V& meanVec,
540  V& samVec) const
541 {
542  bool bRC = ((initialPos < this->subSequenceSize()) &&
543  (0 < numPos ) &&
544  ((initialPos+numPos) <= this->subSequenceSize()) &&
545  (this->vectorSizeLocal() == meanVec.sizeLocal() ) &&
546  (this->vectorSizeLocal() == samVec.sizeLocal() ));
547  queso_require_msg(bRC, "invalid input data");
548 
549  ScalarSequence<double> data(m_env,0,"");
550 
551  unsigned int numParams = this->vectorSizeLocal();
552  for (unsigned int i = 0; i < numParams; ++i) {
553  this->extractScalarSeq(initialPos,
554  1, // spacing
555  numPos,
556  i,
557  data);
558  samVec[i] = data.subSampleVarianceExtra(0,
559  numPos,
560  meanVec[i]);
561  }
562 
563  return;
564 }
565 //---------------------------------------------------
566 template <class V, class M>
567 void
568 SequenceOfVectors<V,M>::unifiedSampleVarianceExtra(
569  unsigned int initialPos,
570  unsigned int numPos,
571  const V& unifiedMeanVec,
572  V& unifiedSamVec) const
573 {
574  bool bRC = ((initialPos < this->subSequenceSize() ) &&
575  (0 < numPos ) &&
576  ((initialPos+numPos) <= this->subSequenceSize() ) &&
577  (this->vectorSizeLocal() == unifiedMeanVec.sizeLocal()) &&
578  (this->vectorSizeLocal() == unifiedSamVec.sizeLocal() ));
579  queso_require_msg(bRC, "invalid input data");
580 
581  ScalarSequence<double> data(m_env,0,"");
582 
583  unsigned int numParams = this->vectorSizeLocal();
584  for (unsigned int i = 0; i < numParams; ++i) {
585  this->extractScalarSeq(initialPos,
586  1, // spacing
587  numPos,
588  i,
589  data);
590  unifiedSamVec[i] = data.unifiedSampleVarianceExtra(m_vectorSpace.numOfProcsForStorage() == 1,
591  0,
592  numPos,
593  unifiedMeanVec[i]);
594  }
595 
596  return;
597 }
598 //---------------------------------------------------
599 template <class V, class M>
600 void
601 SequenceOfVectors<V,M>::subSampleStd(
602  unsigned int initialPos,
603  unsigned int numPos,
604  const V& meanVec,
605  V& stdvec) const
606 {
607  bool bRC = ((initialPos < this->subSequenceSize()) &&
608  (0 < numPos ) &&
609  ((initialPos+numPos) <= this->subSequenceSize()) &&
610  (this->vectorSizeLocal() == meanVec.sizeLocal() ) &&
611  (this->vectorSizeLocal() == stdvec.sizeLocal() ));
612  queso_require_msg(bRC, "invalid input data");
613 
614  ScalarSequence<double> data(m_env,0,"");
615 
616  unsigned int numParams = this->vectorSizeLocal();
617  for (unsigned int i = 0; i < numParams; ++i) {
618  this->extractScalarSeq(initialPos,
619  1, // spacing
620  numPos,
621  i,
622  data);
623  stdvec[i] = data.subSampleStd(0,
624  numPos,
625  meanVec[i]);
626  }
627 
628  return;
629 }
630 //---------------------------------------------------
631 template <class V, class M>
632 void
633 SequenceOfVectors<V,M>::unifiedSampleStd(
634  unsigned int initialPos,
635  unsigned int numPos,
636  const V& unifiedMeanVec,
637  V& unifiedStdVec) const
638 {
639  bool bRC = ((initialPos < this->subSequenceSize() ) &&
640  (0 < numPos ) &&
641  ((initialPos+numPos) <= this->subSequenceSize() ) &&
642  (this->vectorSizeLocal() == unifiedMeanVec.sizeLocal()) &&
643  (this->vectorSizeLocal() == unifiedStdVec.sizeLocal() ));
644  queso_require_msg(bRC, "invalid input data");
645 
646  ScalarSequence<double> data(m_env,0,"");
647 
648  unsigned int numParams = this->vectorSizeLocal();
649  for (unsigned int i = 0; i < numParams; ++i) {
650  this->extractScalarSeq(initialPos,
651  1, // spacing
652  numPos,
653  i,
654  data);
655  unifiedStdVec[i] = data.unifiedSampleStd(m_vectorSpace.numOfProcsForStorage() == 1,
656  0,
657  numPos,
658  unifiedMeanVec[i]);
659  }
660 
661  return;
662 }
663 //---------------------------------------------------
664 template <class V, class M>
665 void
666 SequenceOfVectors<V,M>::subPopulationVariance(
667  unsigned int initialPos,
668  unsigned int numPos,
669  const V& meanVec,
670  V& popVec) const
671 {
672  bool bRC = ((initialPos < this->subSequenceSize()) &&
673  (0 < numPos ) &&
674  ((initialPos+numPos) <= this->subSequenceSize()) &&
675  (this->vectorSizeLocal() == meanVec.sizeLocal() ) &&
676  (this->vectorSizeLocal() == popVec.sizeLocal() ));
677  queso_require_msg(bRC, "invalid input data");
678 
679  ScalarSequence<double> data(m_env,0,"");
680 
681  unsigned int numParams = this->vectorSizeLocal();
682  for (unsigned int i = 0; i < numParams; ++i) {
683  this->extractScalarSeq(initialPos,
684  1, // spacing
685  numPos,
686  i,
687  data);
688  popVec[i] = data.subPopulationVariance(0,
689  numPos,
690  meanVec[i]);
691  }
692 
693  return;
694 }
695 //---------------------------------------------------
696 template <class V, class M>
697 void
698 SequenceOfVectors<V,M>::unifiedPopulationVariance(
699  unsigned int initialPos,
700  unsigned int numPos,
701  const V& unifiedMeanVec,
702  V& unifiedPopVec) const
703 {
704  bool bRC = ((initialPos < this->subSequenceSize() ) &&
705  (0 < numPos ) &&
706  ((initialPos+numPos) <= this->subSequenceSize() ) &&
707  (this->vectorSizeLocal() == unifiedMeanVec.sizeLocal()) &&
708  (this->vectorSizeLocal() == unifiedPopVec.sizeLocal() ));
709  queso_require_msg(bRC, "invalid input data");
710 
711  ScalarSequence<double> data(m_env,0,"");
712 
713  unsigned int numParams = this->vectorSizeLocal();
714  for (unsigned int i = 0; i < numParams; ++i) {
715  this->extractScalarSeq(initialPos,
716  1, // spacing
717  numPos,
718  i,
719  data);
720  unifiedPopVec[i] = data.unifiedPopulationVariance(m_vectorSpace.numOfProcsForStorage() == 1,
721  0,
722  numPos,
723  unifiedMeanVec[i]);
724  }
725 
726  return;
727 }
728 //---------------------------------------------------
729 template <class V, class M>
730 void
731 SequenceOfVectors<V,M>::autoCovariance(
732  unsigned int initialPos,
733  unsigned int numPos,
734  const V& meanVec,
735  unsigned int lag,
736  V& covVec) const
737 {
738  bool bRC = ((initialPos < this->subSequenceSize()) &&
739  (0 < numPos ) &&
740  ((initialPos+numPos) <= this->subSequenceSize()) &&
741  (this->vectorSizeLocal() == meanVec.sizeLocal() ) &&
742  (lag < numPos ) && // lag should not be too large
743  (this->vectorSizeLocal() == covVec.sizeLocal() ));
744  queso_require_msg(bRC, "invalid input data");
745 
746  ScalarSequence<double> data(m_env,0,"");
747 
748  unsigned int numParams = this->vectorSizeLocal();
749  for (unsigned int i = 0; i < numParams; ++i) {
750  this->extractScalarSeq(initialPos,
751  1, // spacing
752  numPos,
753  i,
754  data);
755  covVec[i] = data.autoCovariance(0,
756  numPos,
757  meanVec[i],
758  lag);
759  }
760 
761  return;
762 }
763 //---------------------------------------------------
764 template <class V, class M>
765 void
766 SequenceOfVectors<V,M>::autoCorrViaDef(
767  unsigned int initialPos,
768  unsigned int numPos,
769  unsigned int lag,
770  V& corrVec) const
771 {
772  bool bRC = ((initialPos < this->subSequenceSize()) &&
773  (0 < numPos ) &&
774  ((initialPos+numPos) <= this->subSequenceSize()) &&
775  (lag < numPos ) && // lag should not be too large
776  (this->vectorSizeLocal() == corrVec.sizeLocal() ));
777  queso_require_msg(bRC, "invalid input data");
778 
779  ScalarSequence<double> data(m_env,0,"");
780 
781  unsigned int numParams = this->vectorSizeLocal();
782  for (unsigned int i = 0; i < numParams; ++i) {
783  this->extractScalarSeq(initialPos,
784  1, // spacing
785  numPos,
786  i,
787  data);
788  corrVec[i] = data.autoCorrViaDef(0,
789  numPos,
790  lag);
791  }
792 
793  return;
794 }
795 //---------------------------------------------------
796 template <class V, class M>
797 void
798 SequenceOfVectors<V,M>::autoCorrViaFft(
799  unsigned int initialPos,
800  unsigned int numPos,
801  const std::vector<unsigned int>& lags,
802  std::vector<V*>& corrVecs) const
803 {
804  bool bRC = ((initialPos < this->subSequenceSize()) &&
805  (0 < numPos ) &&
806  ((initialPos+numPos) <= this->subSequenceSize()) &&
807  (0 < lags.size() ) &&
808  (lags[lags.size()-1] < numPos )); // lag should not be too large
809  queso_require_msg(bRC, "invalid input data");
810 
811  for (unsigned int j = lags.size(); j < corrVecs.size(); ++j) {
812  if (corrVecs[j] != NULL) {
813  delete corrVecs[j];
814  corrVecs[j] = NULL;
815  }
816  }
817  corrVecs.resize(lags.size(),NULL);
818  for (unsigned int j = 0; j < corrVecs.size(); ++j) {
819  if (corrVecs[j] == NULL) corrVecs[j] = new V(m_vectorSpace.zeroVector());
820  }
821 
822  ScalarSequence<double> data(m_env,0,"");
823  unsigned int maxLag = lags[lags.size()-1];
824  std::vector<double> autoCorrs(maxLag+1,0.); // Yes, +1
825 
826  unsigned int numParams = this->vectorSizeLocal();
827  for (unsigned int i = 0; i < numParams; ++i) {
828  this->extractScalarSeq(initialPos,
829  1, // spacing
830  numPos,
831  i,
832  data);
833 
834  //if (m_env.subDisplayFile()) {
835  // *m_env.subDisplayFile() << "In SequenceOfVectors<V,M>::autoCorrViaFft()"
836  // << ": about to call data.autoCorrViaFft() for paramId = " << i
837  // << ", with numPos = " << numPos
838  // << ", maxLag = " << maxLag
839  // << ", autoCorrs.size() = " << autoCorrs.size()
840  // << std::endl;
841  //}
842  data.autoCorrViaFft(0,
843  numPos,
844  maxLag,
845  autoCorrs);
846 
847  for (unsigned int j = 0; j < lags.size(); ++j) {
848  (*(corrVecs[j]))[i] = autoCorrs[lags[j]];
849  }
850  }
851 
852  return;
853 }
854 //---------------------------------------------------
855 template <class V, class M>
856 void
857 SequenceOfVectors<V,M>::autoCorrViaFft(
858  unsigned int initialPos,
859  unsigned int numPos,
860  unsigned int numSum,
861  V& autoCorrsSumVec) const
862 {
863  bool bRC = ((initialPos < this->subSequenceSize()) &&
864  (0 < numPos ) &&
865  ((initialPos+numPos) <= this->subSequenceSize()) &&
866  (0 < numSum ) &&
867  (numSum <= numPos ) &&
868  (autoCorrsSumVec.sizeLocal() == this->vectorSizeLocal()));
869  queso_require_msg(bRC, "invalid input data");
870 
871  ScalarSequence<double> data(m_env,0,"");
872 
873  unsigned int numParams = this->vectorSizeLocal();
874  for (unsigned int i = 0; i < numParams; ++i) {
875  this->extractScalarSeq(initialPos,
876  1, // spacing
877  numPos,
878  i,
879  data);
880 
881  data.autoCorrViaFft(0,
882  numPos,
883  numSum,
884  autoCorrsSumVec[i]);
885  }
886 
887  return;
888 }
889 //---------------------------------------------------
890 template <class V, class M>
891 void
892 SequenceOfVectors<V,M>::subMinMaxExtra(
893  unsigned int initialPos,
894  unsigned int numPos,
895  V& minVec,
896  V& maxVec) const
897 {
898  bool bRC = ((0 < numPos ) &&
899  ((initialPos+numPos) <= this->subSequenceSize()) &&
900  (this->vectorSizeLocal() == minVec.sizeLocal() ) &&
901  (this->vectorSizeLocal() == maxVec.sizeLocal() ));
902  queso_require_msg(bRC, "invalid input data");
903 
904  //unsigned int numPos = this->subSequenceSize() - initialPos;
905  unsigned int numParams = this->vectorSizeLocal();
906  ScalarSequence<double> data(m_env,0,"");
907 
908  for (unsigned int i = 0; i < numParams; ++i) {
909  this->extractScalarSeq(initialPos,
910  1, // spacing
911  numPos,
912  i,
913  data);
914  data.subMinMaxExtra(0,numPos,minVec[i],maxVec[i]);
915  }
916 
917  return;
918 }
919 //---------------------------------------------------
920 template <class V, class M>
921 void
922 SequenceOfVectors<V,M>::unifiedMinMaxExtra(
923  unsigned int initialPos,
924  unsigned int numPos,
925  V& unifiedMinVec,
926  V& unifiedMaxVec) const
927 {
928  bool bRC = ((0 < numPos ) &&
929  ((initialPos+numPos) <= this->subSequenceSize() ) &&
930  (this->vectorSizeLocal() == unifiedMinVec.sizeLocal()) &&
931  (this->vectorSizeLocal() == unifiedMaxVec.sizeLocal()));
932  queso_require_msg(bRC, "invalid input data");
933 
934  //unsigned int numPos = this->subSequenceSize() - initialPos;
935  unsigned int numParams = this->vectorSizeLocal();
936  ScalarSequence<double> data(m_env,0,"");
937 
938  for (unsigned int i = 0; i < numParams; ++i) {
939  this->extractScalarSeq(initialPos,
940  1, // spacing
941  numPos,
942  i,
943  data);
944  data.unifiedMinMaxExtra(m_vectorSpace.numOfProcsForStorage() == 1,
945  0,
946  numPos,
947  unifiedMinVec[i],
948  unifiedMaxVec[i]);
949  }
950 
951  return;
952 }
953 //---------------------------------------------------
954 template <class V, class M>
955 void
956 SequenceOfVectors<V,M>::subHistogram(
957  unsigned int initialPos,
958  const V& minVec,
959  const V& maxVec,
960  std::vector<V*>& centersForAllBins,
961  std::vector<V*>& quanttsForAllBins) const
962 {
963  bool bRC = ((initialPos < this->subSequenceSize() ) &&
964  (this->vectorSizeLocal() == minVec.sizeLocal() ) &&
965  (this->vectorSizeLocal() == maxVec.sizeLocal() ) &&
966  (0 < centersForAllBins.size()) &&
967  (centersForAllBins.size() == quanttsForAllBins.size()));
968  queso_require_msg(bRC, "invalid input data");
969 
970  for (unsigned int j = 0; j < quanttsForAllBins.size(); ++j) {
971  centersForAllBins[j] = new V(m_vectorSpace.zeroVector());
972  quanttsForAllBins [j] = new V(m_vectorSpace.zeroVector());
973  }
974 
975  unsigned int dataSize = this->subSequenceSize() - initialPos;
976  unsigned int numParams = this->vectorSizeLocal();
977  for (unsigned int i = 0; i < numParams; ++i) {
978  ScalarSequence<double> data(m_env,dataSize,"");
979  for (unsigned int j = 0; j < dataSize; ++j) {
980  data[j] = (*(m_seq[initialPos+j]))[i];
981  }
982 
983  std::vector<double > centers(centersForAllBins.size(),0.);
984  std::vector<unsigned int> quantts(quanttsForAllBins.size(), 0 );
985  data.subHistogram(0,
986  minVec[i],
987  maxVec[i],
988  centers,
989  quantts);
990 
991  for (unsigned int j = 0; j < quantts.size(); ++j) {
992  (*(centersForAllBins[j]))[i] = centers[j];
993  (*(quanttsForAllBins[j]))[i] = (double) quantts[j];
994  }
995  }
996 
997  return;
998 }
999 //---------------------------------------------------
1000 template <class V, class M>
1001 void
1002 SequenceOfVectors<V,M>::unifiedHistogram(
1003  unsigned int initialPos,
1004  const V& unifiedMinVec,
1005  const V& unifiedMaxVec,
1006  std::vector<V*>& unifiedCentersForAllBins,
1007  std::vector<V*>& unifiedQuanttsForAllBins) const
1008 {
1009  bool bRC = ((initialPos < this->subSequenceSize() ) &&
1010  (this->vectorSizeLocal() == unifiedMinVec.sizeLocal() ) &&
1011  (this->vectorSizeLocal() == unifiedMaxVec.sizeLocal() ) &&
1012  (0 < unifiedCentersForAllBins.size()) &&
1013  (unifiedCentersForAllBins.size() == unifiedQuanttsForAllBins.size()));
1014  queso_require_msg(bRC, "invalid input data");
1015 
1016  for (unsigned int j = 0; j < unifiedQuanttsForAllBins.size(); ++j) {
1017  unifiedCentersForAllBins[j] = new V(m_vectorSpace.zeroVector());
1018  unifiedQuanttsForAllBins [j] = new V(m_vectorSpace.zeroVector());
1019  }
1020 
1021  unsigned int dataSize = this->subSequenceSize() - initialPos;
1022  unsigned int numParams = this->vectorSizeLocal();
1023  for (unsigned int i = 0; i < numParams; ++i) {
1024  ScalarSequence<double> data(m_env,dataSize,"");
1025  for (unsigned int j = 0; j < dataSize; ++j) {
1026  data[j] = (*(m_seq[initialPos+j]))[i];
1027  }
1028 
1029  std::vector<double > unifiedCenters(unifiedCentersForAllBins.size(),0.);
1030  std::vector<unsigned int> unifiedQuantts(unifiedQuanttsForAllBins.size(), 0 );
1031  data.unifiedHistogram(m_vectorSpace.numOfProcsForStorage() == 1,
1032  0,
1033  unifiedMinVec[i],
1034  unifiedMaxVec[i],
1035  unifiedCenters,
1036  unifiedQuantts);
1037 
1038  for (unsigned int j = 0; j < unifiedQuantts.size(); ++j) {
1039  (*(unifiedCentersForAllBins[j]))[i] = unifiedCenters[j];
1040  (*(unifiedQuanttsForAllBins[j]))[i] = (double) unifiedQuantts[j];
1041  }
1042  }
1043 
1044  return;
1045 }
1046 //---------------------------------------------------
1047 template <class V, class M>
1048 void
1049 SequenceOfVectors<V,M>::subInterQuantileRange(
1050  unsigned int initialPos,
1051  V& iqrVec) const
1052 {
1053  bool bRC = ((initialPos < this->subSequenceSize()) &&
1054  (this->vectorSizeLocal() == iqrVec.sizeLocal() ));
1055  queso_require_msg(bRC, "invalid input data");
1056 
1057  unsigned int numPos = this->subSequenceSize() - initialPos;
1058  ScalarSequence<double> data(m_env,0,"");
1059 
1060  unsigned int numParams = this->vectorSizeLocal();
1061  for (unsigned int i = 0; i < numParams; ++i) {
1062  this->extractScalarSeq(initialPos,
1063  1, // spacing
1064  numPos,
1065  i,
1066  data);
1067  iqrVec[i] = data.subInterQuantileRange(0);
1068  }
1069 
1070  return;
1071 }
1072 //---------------------------------------------------
1073 template <class V, class M>
1074 void
1075 SequenceOfVectors<V,M>::unifiedInterQuantileRange(
1076  unsigned int initialPos,
1077  V& unifiedIqrVec) const
1078 {
1079  bool bRC = ((initialPos < this->subSequenceSize() ) &&
1080  (this->vectorSizeLocal() == unifiedIqrVec.sizeLocal()));
1081  queso_require_msg(bRC, "invalid input data");
1082 
1083  unsigned int numPos = this->subSequenceSize() - initialPos;
1084  ScalarSequence<double> data(m_env,0,"");
1085 
1086  unsigned int numParams = this->vectorSizeLocal();
1087  for (unsigned int i = 0; i < numParams; ++i) {
1088  this->extractScalarSeq(initialPos,
1089  1, // spacing
1090  numPos,
1091  i,
1092  data);
1093  unifiedIqrVec[i] = data.unifiedInterQuantileRange(m_vectorSpace.numOfProcsForStorage() == 1,
1094  0);
1095  }
1096 
1097  return;
1098 }
1099 //---------------------------------------------------
1100 template <class V, class M>
1101 void
1102 SequenceOfVectors<V,M>::subScalesForKde(
1103  unsigned int initialPos,
1104  const V& iqrVec,
1105  unsigned int kdeDimension,
1106  V& scaleVec) const
1107 {
1108  bool bRC = ((initialPos < this->subSequenceSize()) &&
1109  (this->vectorSizeLocal() == iqrVec.sizeLocal() ) &&
1110  (this->vectorSizeLocal() == scaleVec.sizeLocal() ));
1111  queso_require_msg(bRC, "invalid input data");
1112 
1113  unsigned int numPos = this->subSequenceSize() - initialPos;
1114  ScalarSequence<double> data(m_env,0,"");
1115 
1116  unsigned int numParams = this->vectorSizeLocal();
1117  for (unsigned int i = 0; i < numParams; ++i) {
1118  this->extractScalarSeq(initialPos,
1119  1, // spacing
1120  numPos,
1121  i,
1122  data);
1123  scaleVec[i] = data.subScaleForKde(0,
1124  iqrVec[i],
1125  kdeDimension);
1126  }
1127 
1128  return;
1129 }
1130 //---------------------------------------------------
1131 template <class V, class M>
1132 void
1133 SequenceOfVectors<V,M>::unifiedScalesForKde(
1134  unsigned int initialPos,
1135  const V& unifiedIqrVec,
1136  unsigned int kdeDimension,
1137  V& unifiedScaleVec) const
1138 {
1139  bool bRC = ((initialPos < this->subSequenceSize() ) &&
1140  (this->vectorSizeLocal() == unifiedIqrVec.sizeLocal() ) &&
1141  (this->vectorSizeLocal() == unifiedScaleVec.sizeLocal()));
1142  queso_require_msg(bRC, "invalid input data");
1143 
1144  unsigned int numPos = this->subSequenceSize() - initialPos;
1145  ScalarSequence<double> data(m_env,0,"");
1146 
1147  unsigned int numParams = this->vectorSizeLocal();
1148  for (unsigned int i = 0; i < numParams; ++i) {
1149  this->extractScalarSeq(initialPos,
1150  1, // spacing
1151  numPos,
1152  i,
1153  data);
1154  unifiedScaleVec[i] = data.unifiedScaleForKde(m_vectorSpace.numOfProcsForStorage() == 1,
1155  0,
1156  unifiedIqrVec[i],
1157  kdeDimension);
1158  }
1159 
1160  return;
1161 }
1162 //---------------------------------------------------
1163 template <class V, class M>
1164 void
1165 SequenceOfVectors<V,M>::subGaussian1dKde(
1166  unsigned int initialPos,
1167  const V& scaleVec,
1168  const std::vector<V*>& evalParamVecs,
1169  std::vector<V*>& densityVecs) const
1170 {
1171  bool bRC = ((initialPos < this->subSequenceSize()) &&
1172  (this->vectorSizeLocal() == scaleVec.sizeLocal() ) &&
1173  (0 < evalParamVecs.size() ) &&
1174  (evalParamVecs.size() == densityVecs.size() ));
1175  queso_require_msg(bRC, "invalid input data");
1176 
1177  unsigned int numPos = this->subSequenceSize() - initialPos;
1178  ScalarSequence<double> data(m_env,0,"");
1179 
1180  unsigned int numEvals = evalParamVecs.size();
1181  for (unsigned int j = 0; j < numEvals; ++j) {
1182  densityVecs[j] = new V(m_vectorSpace.zeroVector());
1183  }
1184  std::vector<double> evalParams(numEvals,0.);
1185  std::vector<double> densities (numEvals,0.);
1186 
1187  unsigned int numParams = this->vectorSizeLocal();
1188  for (unsigned int i = 0; i < numParams; ++i) {
1189  this->extractScalarSeq(initialPos,
1190  1, // spacing
1191  numPos,
1192  i,
1193  data);
1194 
1195  for (unsigned int j = 0; j < numEvals; ++j) {
1196  evalParams[j] = (*evalParamVecs[j])[i];
1197  }
1198 
1199  data.subGaussian1dKde(0,
1200  scaleVec[i],
1201  evalParams,
1202  densities);
1203 
1204  for (unsigned int j = 0; j < numEvals; ++j) {
1205  (*densityVecs[j])[i] = densities[j];
1206  }
1207  }
1208 
1209  return;
1210 }
1211 //---------------------------------------------------
1212 template <class V, class M>
1213 void
1214 SequenceOfVectors<V,M>::unifiedGaussian1dKde(
1215  unsigned int initialPos,
1216  const V& unifiedScaleVec,
1217  const std::vector<V*>& unifiedEvalParamVecs,
1218  std::vector<V*>& unifiedDensityVecs) const
1219 {
1220  bool bRC = ((initialPos < this->subSequenceSize() ) &&
1221  (this->vectorSizeLocal() == unifiedScaleVec.sizeLocal()) &&
1222  (0 < unifiedEvalParamVecs.size()) &&
1223  (unifiedEvalParamVecs.size() == unifiedDensityVecs.size() ));
1224  queso_require_msg(bRC, "invalid input data");
1225 
1226  unsigned int numPos = this->subSequenceSize() - initialPos;
1227  ScalarSequence<double> data(m_env,0,"");
1228 
1229  unsigned int numEvals = unifiedEvalParamVecs.size();
1230  for (unsigned int j = 0; j < numEvals; ++j) {
1231  unifiedDensityVecs[j] = new V(m_vectorSpace.zeroVector());
1232  }
1233  std::vector<double> unifiedEvalParams(numEvals,0.);
1234  std::vector<double> unifiedDensities (numEvals,0.);
1235 
1236  unsigned int numParams = this->vectorSizeLocal();
1237  for (unsigned int i = 0; i < numParams; ++i) {
1238  this->extractScalarSeq(initialPos,
1239  1, // spacing
1240  numPos,
1241  i,
1242  data);
1243 
1244  for (unsigned int j = 0; j < numEvals; ++j) {
1245  unifiedEvalParams[j] = (*unifiedEvalParamVecs[j])[i];
1246  }
1247 
1248  data.unifiedGaussian1dKde(m_vectorSpace.numOfProcsForStorage() == 1,
1249  0,
1250  unifiedScaleVec[i],
1251  unifiedEvalParams,
1252  unifiedDensities);
1253 
1254  for (unsigned int j = 0; j < numEvals; ++j) {
1255  (*unifiedDensityVecs[j])[i] = unifiedDensities[j];
1256  }
1257  }
1258 
1259  return;
1260 }
1261 //---------------------------------------------------
1262 template <class V, class M>
1263 void
1264 SequenceOfVectors<V,M>::subWriteContents(
1265  unsigned int initialPos,
1266  unsigned int numPos,
1267  const std::string& fileName,
1268  const std::string& fileType,
1269  const std::set<unsigned int>& allowedSubEnvIds) const
1270 {
1271  queso_require_greater_equal_msg(m_env.subRank(), 0, "unexpected subRank");
1272 
1273  FilePtrSetStruct filePtrSet;
1274  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 99)) {
1275  *m_env.subDisplayFile() << "In SequenceOfVectors<V,M>::subWriteContents()"
1276  << ": about to try to open file '" << fileName << "." << fileType
1277  << "'"
1278  << ", initialPos = " << initialPos
1279  << ", numPos = " << numPos
1280  << std::endl;
1281  }
1282  if (m_env.openOutputFile(fileName,
1283  fileType,
1284  allowedSubEnvIds,
1285  false, // A 'true' causes problems when the user chooses (via options
1286  // in the input file) to use just one file for all outputs.
1287  filePtrSet)) {
1288  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 99)) {
1289  *m_env.subDisplayFile() << "In SequenceOfVectors<V,M>::subWriteContents()"
1290  << ": successfully opened file '" << fileName << "." << fileType
1291  << "'"
1292  << std::endl;
1293  }
1294  this->subWriteContents(initialPos,
1295  numPos,
1296  filePtrSet,
1297  fileType);
1298  m_env.closeFile(filePtrSet,fileType);
1299  }
1300  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 99)) {
1301  *m_env.subDisplayFile() << "In SequenceOfVectors<V,M>::subWriteContents()"
1302  << ": before Barrier()"
1303  << std::endl;
1304  }
1305  m_env.subComm().Barrier();
1306 
1307  return;
1308 }
1309 //---------------------------------------------------
1310 template <class V, class M>
1311 void
1312 SequenceOfVectors<V,M>::subWriteContents(
1313  unsigned int initialPos,
1314  unsigned int numPos,
1315  FilePtrSetStruct& filePtrSet,
1316  const std::string& fileType) const
1317 {
1318  if (fileType == UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT ||
1319  fileType == UQ_FILE_EXTENSION_FOR_TXT_FORMAT) {
1320  queso_require_msg(filePtrSet.ofsVar, "filePtrSet.ofsVar should not be NULL");
1321  this->subWriteContents(initialPos,
1322  numPos,
1323  *filePtrSet.ofsVar,
1324  fileType);
1325  }
1326 #ifdef QUESO_HAS_HDF5
1327  else if (fileType == UQ_FILE_EXTENSION_FOR_HDF_FORMAT) {
1328 
1329  // Check the file is still legit
1330  queso_require_greater_equal_msg(
1331  filePtrSet.h5Var,
1332  0,
1333  "filePtrSet.h5Var should not be non-negative");
1334 
1335  // Sanity check extent
1336  queso_require_less_equal_msg(
1337  (initialPos+numPos),
1338  this->subSequenceSize(),
1339  "invalid routine input parameters");
1340 
1341  unsigned int numParams = m_vectorSpace.dimLocal();
1342  unsigned int chainSize = this->subSequenceSize();
1343  hsize_t dims[2] = { chainSize, numParams };
1344 
1345  // Create dataspace
1346  hid_t dataspace_id = H5Screate_simple(2, dims, dims);
1347 
1348  // Sanity check dataspace
1349  queso_require_greater_equal_msg(
1350  dataspace_id,
1351  0,
1352  "error creating dataspace with id: " << dataspace_id);
1353 
1354  // Create dataset
1355  hid_t dataset_id = H5Dcreate(filePtrSet.h5Var,
1356  "data",
1357  H5T_IEEE_F64LE,
1358  dataspace_id,
1359  H5P_DEFAULT,
1360  H5P_DEFAULT,
1361  H5P_DEFAULT);
1362 
1363  // Sanity check dataset
1364  queso_require_greater_equal_msg(
1365  dataset_id,
1366  0,
1367  "error creating dataset with id: " << dataset_id);
1368 
1369  // This is so egregiously awfully badly terribly horrific I want to die.
1370  //
1371  // And, of course, if any of the subsequent H5* sanity check fail we
1372  // throw an exception and leak a metric fuckton of memory.
1373  double * data = (double *)malloc(numParams * chainSize * sizeof(double));
1374 
1375  for (unsigned int i = 0; i < chainSize; i++) {
1376  V tmpVec(*(m_seq[i]));
1377  for (unsigned int j = 0; j < numParams; j++) {
1378  data[numParams*i+j] = tmpVec[j];
1379  }
1380  }
1381 
1382  // Write the dataset
1383  herr_t status = H5Dwrite(
1384  dataset_id,
1385  H5T_NATIVE_DOUBLE, // The type in memory
1386  H5S_ALL, // The dataspace in memory
1387  dataspace_id, // The file dataspace
1388  H5P_DEFAULT, // Xfer property list
1389  data);
1390 
1391  // Sanity check the write
1392  queso_require_greater_equal_msg(
1393  status,
1394  0,
1395  "error writing dataset to file with id: " << filePtrSet.h5Var);
1396 
1397  // Clean up
1398  free(data);
1399  H5Dclose(dataset_id);
1400  H5Sclose(dataspace_id);
1401 
1402  // Should we close the file too? It's unclear.
1403 
1404  }
1405 #endif
1406  else {
1407  queso_error_msg("invalid file type");
1408  }
1409 
1410  return;
1411 }
1412 //---------------------------------------------------
1413 template <class V, class M>
1414 void
1415 SequenceOfVectors<V,M>::subWriteContents(
1416  unsigned int initialPos,
1417  unsigned int numPos,
1418  std::ofstream& ofs,
1419  const std::string& fileType) const
1420 {
1421  queso_require_less_equal_msg((initialPos+numPos), this->subSequenceSize(), "invalid routine input parameters");
1422 
1423  if (initialPos == 0) {
1424  if (fileType == UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT) {
1425  // Need sub matlab header here since this is subWriteContents
1426  this->writeSubMatlabHeader(ofs,
1427  this->subSequenceSize(),
1428  this->vectorSizeLocal());
1429  }
1430  else if (fileType == UQ_FILE_EXTENSION_FOR_TXT_FORMAT) {
1431  this->writeTxtHeader(ofs,
1432  this->subSequenceSize(),
1433  this->vectorSizeLocal());
1434  }
1435  }
1436 
1437  for (unsigned int j = initialPos; j < initialPos+numPos; ++j) {
1438  bool savedVectorPrintScientific = m_seq[j]->getPrintScientific();
1439  bool savedVectorPrintState = m_seq[j]->getPrintHorizontally();
1440  m_seq[j]->setPrintScientific (true);
1441  m_seq[j]->setPrintHorizontally(true);
1442 
1443  ofs << *(m_seq[j])
1444  << std::endl;
1445 
1446  m_seq[j]->setPrintHorizontally(savedVectorPrintState);
1447  m_seq[j]->setPrintScientific (savedVectorPrintScientific);
1448  }
1449 
1450  // Write Matlab-specific ending if desired
1451  if ((fileType == UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT) &&
1452  ((initialPos + numPos) == this->subSequenceSize())) {
1453  ofs << "];\n";
1454  }
1455 }
1456 
1457 template <class V, class M>
1458 void
1459 SequenceOfVectors<V,M>::writeSubMatlabHeader(std::ofstream & ofs,
1460  double sequenceSize, double vectorSizeLocal) const
1461 {
1462  ofs << m_name << "_sub" << m_env.subIdString() << " = zeros(" << sequenceSize
1463  << "," << vectorSizeLocal
1464  << ");"
1465  << std::endl;
1466  ofs << m_name << "_sub" << m_env.subIdString() << " = [";
1467 }
1468 
1469 template <class V, class M>
1470 void
1471 SequenceOfVectors<V,M>::writeUnifiedMatlabHeader(std::ofstream & ofs,
1472  double sequenceSize, double vectorSizeLocal) const
1473 {
1474  ofs << m_name << "_unified" << " = zeros(" << sequenceSize
1475  << "," << vectorSizeLocal
1476  << ");"
1477  << std::endl;
1478  ofs<< m_name << "_unified" << " = [";
1479 }
1480 
1481 template <class V, class M>
1482 void
1483 SequenceOfVectors<V,M>::writeTxtHeader(std::ofstream & ofs,
1484  double sequenceSize, double vectorSizeLocal) const
1485 {
1486  ofs << sequenceSize << " " << vectorSizeLocal
1487  << std::endl;
1488 }
1489 
1490 //---------------------------------------------------
1491 template <class V, class M>
1492 void
1493 SequenceOfVectors<V,M>::unifiedWriteContents(
1494  const std::string& fileName,
1495  const std::string& inputFileType) const
1496 {
1497  std::string fileType(inputFileType);
1498 #ifdef QUESO_HAS_HDF5
1499  // Do nothing
1500 #else
1501  if (fileType == UQ_FILE_EXTENSION_FOR_HDF_FORMAT) {
1502  if (m_env.subDisplayFile()) {
1503  *m_env.subDisplayFile() << "WARNING in SequenceOfVectors<V,M>::unifiedWriteContents()"
1504  << ": file format '" << UQ_FILE_EXTENSION_FOR_HDF_FORMAT
1505  << "' has been requested, but this QUESO library has not been built with 'hdf5'"
1506  << ". Code will therefore process the file format '" << UQ_FILE_EXTENSION_FOR_HDF_FORMAT
1507  << "' instead..."
1508  << std::endl;
1509  }
1510  if (m_env.subRank() == 0) {
1511  std::cerr << "WARNING in SequenceOfVectors<V,M>::unifiedWriteContents()"
1512  << ": file format '" << UQ_FILE_EXTENSION_FOR_HDF_FORMAT
1513  << "' has been requested, but this QUESO library has not been built with 'hdf5'"
1514  << ". Code will therefore process the file format '" << UQ_FILE_EXTENSION_FOR_HDF_FORMAT
1515  << "' instead..."
1516  << std::endl;
1517  }
1518  fileType = UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT;
1519  }
1520 #endif
1521 
1522  // All processors in 'fullComm' should call this routine...
1523 
1524  //m_env.fullComm().Barrier(); // Dangerous to barrier on fullComm ... // prudenci-2011-01-17
1525  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 10)) {
1526  *m_env.subDisplayFile() << "Entering SequenceOfVectors<V,M>::unifiedWriteContents()"
1527  << ": worldRank " << m_env.worldRank()
1528  << ", fullRank " << m_env.fullRank()
1529  << ", subEnvironment " << m_env.subId()
1530  << ", subRank " << m_env.subRank()
1531  << ", inter0Rank " << m_env.inter0Rank()
1532  //<< ", m_env.inter0Comm().NumProc() = " << m_env.inter0Comm().NumProc()
1533  << ", fileName = " << fileName
1534  << std::endl;
1535  }
1536 
1537  if (m_env.inter0Rank() >= 0) {
1538  for (unsigned int r = 0; r < (unsigned int) m_env.inter0Comm().NumProc(); ++r) {
1539  if (m_env.inter0Rank() == (int) r) {
1540  // My turn
1541  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 10)) {
1542  *m_env.subDisplayFile() << "In SequenceOfVectors<V,M>::unifiedWriteContents()"
1543  << ": worldRank " << m_env.worldRank()
1544  << ", fullRank " << m_env.fullRank()
1545  << ", subEnvironment " << m_env.subId()
1546  << ", subRank " << m_env.subRank()
1547  << ", inter0Rank " << m_env.inter0Rank()
1548  //<< ", m_env.inter0Comm().NumProc() = " << m_env.inter0Comm().NumProc()
1549  << ", fileName = " << fileName
1550  << ", about to open file for r = " << r
1551  << std::endl;
1552  }
1553 
1554  // bool writeOver = (r == 0);
1555  bool writeOver = false; // A 'true' causes problems when the user chooses (via options
1556  // in the input file) to use just one file for all outputs.
1557  FilePtrSetStruct unifiedFilePtrSet;
1558  if (m_env.openUnifiedOutputFile(fileName,
1559  fileType, // "m or hdf"
1560  writeOver,
1561  unifiedFilePtrSet)) {
1562  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 10)) { // 2013-02-23
1563  *m_env.subDisplayFile() << "In SequenceOfVectors<V,M>::unifiedWriteContents()"
1564  << ": worldRank " << m_env.worldRank()
1565  << ", fullRank " << m_env.fullRank()
1566  << ", subEnvironment " << m_env.subId()
1567  << ", subRank " << m_env.subRank()
1568  << ", inter0Rank " << m_env.inter0Rank()
1569  //<< ", m_env.inter0Comm().NumProc() = " << m_env.inter0Comm().NumProc()
1570  << ", fileName = " << fileName
1571  << ", just opened file for r = " << r
1572  << std::endl;
1573  }
1574 
1575  unsigned int chainSize = this->subSequenceSize();
1576  if ((fileType == UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT) ||
1577  (fileType == UQ_FILE_EXTENSION_FOR_TXT_FORMAT)) {
1578  if (r == 0) {
1579  if (fileType == UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT) {
1580  // Need unified matlab header here since this is unifiedWriteContents
1581  writeUnifiedMatlabHeader(*unifiedFilePtrSet.ofsVar,
1582  this->subSequenceSize()*m_env.inter0Comm().NumProc(),
1583  this->vectorSizeLocal());
1584  }
1585  else { // If we get here it's definitely a txt file not matlab
1586  writeTxtHeader(*unifiedFilePtrSet.ofsVar,
1587  this->subSequenceSize()*m_env.inter0Comm().NumProc(),
1588  this->vectorSizeLocal());
1589  }
1590  }
1591 
1592  for (unsigned int j = 0; j < chainSize; ++j) { // 2013-02-23
1593  //std::cout << "In SequenceOfVectors<V,M>::unifiedWriteContents(): m_seq[" << j << "] = " << m_seq[j]
1594  // << std::endl;
1595  //std::cout << "In SequenceOfVectors<V,M>::unifiedWriteContents(): &(m_seq[" << j << "].map()) = " << &(m_seq[j]->map())
1596  // << std::endl;
1597  //std::cout << "In SequenceOfVectors<V,M>::unifiedWriteContents(): (m_seq[" << j << "].map().NumMyElements = " << m_seq[j]->map().NumMyElements()
1598  // << std::endl;
1599  //V tmpVec(*(m_seq[j]));
1600  //std::cout << "*(m_seq[" << j << "]) = " << tmpVec
1601  // << std::endl;
1602  //std::cout << "*(m_seq[" << j << "]) = " << *(m_seq[j])
1603  // << std::endl;
1604 
1605  bool savedVectorPrintScientific = m_seq[j]->getPrintScientific();
1606  bool savedVectorPrintState = m_seq[j]->getPrintHorizontally();
1607  m_seq[j]->setPrintScientific (true);
1608  m_seq[j]->setPrintHorizontally(true);
1609 
1610  *unifiedFilePtrSet.ofsVar << *(m_seq[j])
1611  << std::endl;
1612 
1613  m_seq[j]->setPrintHorizontally(savedVectorPrintState);
1614  m_seq[j]->setPrintScientific (savedVectorPrintScientific);
1615  }
1616  }
1617 #ifdef QUESO_HAS_HDF5
1618  else if (fileType == UQ_FILE_EXTENSION_FOR_HDF_FORMAT) {
1619  unsigned int numParams = m_vectorSpace.dimLocal();
1620  if (r == 0) {
1621  hid_t datatype = H5Tcopy(H5T_NATIVE_DOUBLE);
1622  //std::cout << "In SequenceOfVectors<V,M>::unifiedWriteContents(): h5 case, data type created" << std::endl;
1623  hsize_t dimsf[2];
1624  dimsf[0] = chainSize;
1625  dimsf[1] = numParams;
1626  hid_t dataspace = H5Screate_simple(2, dimsf, NULL); // HDF5_rank = 2
1627  //std::cout << "In SequenceOfVectors<V,M>::unifiedWriteContents(): h5 case, data space created" << std::endl;
1628  hid_t dataset = H5Dcreate2(unifiedFilePtrSet.h5Var,
1629  "data",
1630  datatype,
1631  dataspace,
1632  H5P_DEFAULT, // Link creation property list
1633  H5P_DEFAULT, // Dataset creation property list
1634  H5P_DEFAULT); // Dataset access property list
1635  //std::cout << "In SequenceOfVectors<V,M>::unifiedWriteContents(): h5 case, data set created" << std::endl;
1636 
1637  struct timeval timevalBegin;
1638  int iRC = UQ_OK_RC;
1639  iRC = gettimeofday(&timevalBegin,NULL);
1640  if (iRC) {}; // just to remover compiler warning
1641 
1642  double * data;
1643  data = (double *)malloc(numParams * chainSize * sizeof(double));
1644 
1645  for (unsigned int i = 0; i < chainSize; ++i) {
1646  V tmpVec(*(m_seq[i]));
1647  for (unsigned int j = 0; j < numParams; ++j) {
1648  data[numParams*i+j] = tmpVec[j];
1649  }
1650  }
1651 
1652  herr_t status;
1653  status = H5Dwrite(dataset,
1654  H5T_NATIVE_DOUBLE,
1655  H5S_ALL,
1656  H5S_ALL,
1657  H5P_DEFAULT,
1658  data);
1659  if (status) {}; // just to remover compiler warning
1660  //std::cout << "In SequenceOfVectors<V,M>::unifiedWriteContents(): h5 case, data written" << std::endl;
1661 
1662  double writeTime = MiscGetEllapsedSeconds(&timevalBegin);
1663  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 2)) {
1664  *m_env.subDisplayFile() << "In SequenceOfVectors<V,M>::unifiedWriteContents()"
1665  << ": worldRank " << m_env.worldRank()
1666  << ", fullRank " << m_env.fullRank()
1667  << ", subEnvironment " << m_env.subId()
1668  << ", subRank " << m_env.subRank()
1669  << ", inter0Rank " << m_env.inter0Rank()
1670  << ", fileName = " << fileName
1671  << ", numParams = " << numParams
1672  << ", chainSize = " << chainSize
1673  << ", writeTime = " << writeTime << " seconds"
1674  << std::endl;
1675  }
1676 
1677  H5Dclose(dataset);
1678  //std::cout << "In SequenceOfVectors<V,M>::unifiedWriteContents(): h5 case, data set closed" << std::endl;
1679  H5Sclose(dataspace);
1680  //std::cout << "In SequenceOfVectors<V,M>::unifiedWriteContents(): h5 case, data space closed" << std::endl;
1681  H5Tclose(datatype);
1682  //std::cout << "In SequenceOfVectors<V,M>::unifiedWriteContents(): h5 case, data type closed" << std::endl;
1683  free(data);
1684  }
1685  else {
1686  queso_error_msg("hdf file type not supported for multiple sub-environments yet");
1687  }
1688  }
1689 #endif
1690  else {
1691  queso_error_msg("invalid file type");
1692  }
1693 
1694  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 10)) {
1695  *m_env.subDisplayFile() << "In SequenceOfVectors<V,M>::unifiedWriteContents()"
1696  << ": worldRank " << m_env.worldRank()
1697  << ", fullRank " << m_env.fullRank()
1698  << ", subEnvironment " << m_env.subId()
1699  << ", subRank " << m_env.subRank()
1700  << ", inter0Rank " << m_env.inter0Rank()
1701  //<< ", m_env.inter0Comm().NumProc() = " << m_env.inter0Comm().NumProc()
1702  << ", fileName = " << fileName
1703  << ", about to close file for r = " << r
1704  << std::endl;
1705  }
1706 
1707  m_env.closeFile(unifiedFilePtrSet,fileType);
1708 
1709  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 10)) {
1710  *m_env.subDisplayFile() << "In SequenceOfVectors<V,M>::unifiedWriteContents()"
1711  << ": worldRank " << m_env.worldRank()
1712  << ", fullRank " << m_env.fullRank()
1713  << ", subEnvironment " << m_env.subId()
1714  << ", subRank " << m_env.subRank()
1715  << ", inter0Rank " << m_env.inter0Rank()
1716  //<< ", m_env.inter0Comm().NumProc() = " << m_env.inter0Comm().NumProc()
1717  << ", fileName = " << fileName
1718  << ", just closed file for r = " << r
1719  << std::endl;
1720  }
1721  } // if (m_env.openUnifiedOutputFile())
1722  } // if (m_env.inter0Rank() == (int) r)
1723  m_env.inter0Comm().Barrier();
1724  } // for r
1725 
1726  if (m_env.inter0Rank() == 0) {
1727  if ((fileType == UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT) ||
1728  (fileType == UQ_FILE_EXTENSION_FOR_TXT_FORMAT)) {
1729  FilePtrSetStruct unifiedFilePtrSet;
1730  if (m_env.openUnifiedOutputFile(fileName,
1731  fileType,
1732  false, // Yes, 'writeOver = false' in order to close the array for matlab
1733  unifiedFilePtrSet)) {
1734 
1735  if (fileType == UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT) {
1736  *unifiedFilePtrSet.ofsVar << "];\n";
1737  }
1738 
1739  m_env.closeFile(unifiedFilePtrSet,fileType);
1740  }
1741  }
1742  else if (fileType == UQ_FILE_EXTENSION_FOR_HDF_FORMAT) {
1743  // Do nothing
1744  }
1745  else {
1746  queso_error_msg("invalid file type");
1747  }
1748  }
1749  } // if (m_env.inter0Rank() >= 0)
1750 
1751  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 10)) {
1752  *m_env.subDisplayFile() << "Leaving SequenceOfVectors<V,M>::unifiedWriteContents()"
1753  << ", fileName = " << fileName
1754  << std::endl;
1755  }
1756  //m_env.fullComm().Barrier(); // Dangerous to barrier on fullComm ... // prudenci-2011-01-17
1757 
1758  return;
1759 }
1760 //---------------------------------------------------
1761 template <class V, class M>
1762 void
1763 SequenceOfVectors<V,M>::unifiedReadContents(
1764  const std::string& fileName,
1765  const std::string& inputFileType,
1766  const unsigned int subReadSize)
1767 {
1768  std::string fileType(inputFileType);
1769 #ifdef QUESO_HAS_HDF5
1770  // Do nothing
1771 #else
1772  if (fileType == UQ_FILE_EXTENSION_FOR_HDF_FORMAT) {
1773  if (m_env.subDisplayFile()) {
1774  *m_env.subDisplayFile() << "WARNING in SequenceOfVectors<V,M>::unifiedReadContents()"
1775  << ": file format '" << UQ_FILE_EXTENSION_FOR_HDF_FORMAT
1776  << "' has been requested, but this QUESO library has not been built with 'hdf5'"
1777  << ". Code will therefore process the file format '" << UQ_FILE_EXTENSION_FOR_HDF_FORMAT
1778  << "' instead..."
1779  << std::endl;
1780  }
1781  if (m_env.subRank() == 0) {
1782  std::cerr << "WARNING in SequenceOfVectors<V,M>::unifiedReadContents()"
1783  << ": file format '" << UQ_FILE_EXTENSION_FOR_HDF_FORMAT
1784  << "' has been requested, but this QUESO library has not been built with 'hdf5'"
1785  << ". Code will therefore process the file format '" << UQ_FILE_EXTENSION_FOR_HDF_FORMAT
1786  << "' instead..."
1787  << std::endl;
1788  }
1789  fileType = UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT;
1790  }
1791 #endif
1792 
1793  //m_env.fullComm().Barrier(); // Dangerous to barrier on fullComm ...
1794  if (m_env.subDisplayFile()) {
1795  *m_env.subDisplayFile() << "Entering SequenceOfVectors<V,M>::unifiedReadContents()"
1796  << ": worldRank " << m_env.worldRank()
1797  << ", fullRank " << m_env.fullRank()
1798  << ", subEnvironment " << m_env.subId()
1799  << ", subRank " << m_env.subRank()
1800  << ", inter0Rank " << m_env.inter0Rank()
1801  //<< ", m_env.inter0Comm().NumProc() = " << m_env.inter0Comm().NumProc()
1802  << ", fileName = " << fileName
1803  << ", subReadSize = " << subReadSize
1804  //<< ", unifiedReadSize = " << unifiedReadSize
1805  << std::endl;
1806  }
1807 
1808  this->resizeSequence(subReadSize);
1809 
1810  if (m_env.inter0Rank() >= 0) {
1811  double unifiedReadSize = subReadSize*m_env.inter0Comm().NumProc();
1812 
1813  // In the logic below, the id of a line' begins with value 0 (zero)
1814  unsigned int idOfMyFirstLine = 1 + m_env.inter0Rank()*subReadSize;
1815  unsigned int idOfMyLastLine = (1 + m_env.inter0Rank())*subReadSize;
1816  unsigned int numParams = this->vectorSizeLocal();
1817 
1818  for (unsigned int r = 0; r < (unsigned int) m_env.inter0Comm().NumProc(); ++r) { // "m or hdf"
1819  if (m_env.inter0Rank() == (int) r) {
1820  // My turn
1821  FilePtrSetStruct unifiedFilePtrSet;
1822  if (m_env.openUnifiedInputFile(fileName,
1823  fileType,
1824  unifiedFilePtrSet)) {
1825  if ((fileType == UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT) ||
1826  (fileType == UQ_FILE_EXTENSION_FOR_TXT_FORMAT)) {
1827  if (r == 0) {
1828  // Read number of chain positions in the file by taking care of the first line,
1829  // which resembles something like 'variable_name = zeros(n_positions,m_params);'
1830  std::string tmpString;
1831 
1832  // Read 'variable name' string
1833  *unifiedFilePtrSet.ifsVar >> tmpString;
1834  //std::cout << "Just read '" << tmpString << "'" << std::endl;
1835 
1836  // Read '=' sign
1837  *unifiedFilePtrSet.ifsVar >> tmpString;
1838  //std::cout << "Just read '" << tmpString << "'" << std::endl;
1839  queso_require_equal_to_msg(tmpString, std::string("="), std::string("string should be the '=' sign"));
1840 
1841  // Read 'zeros(n_positions,n_params)' string
1842  *unifiedFilePtrSet.ifsVar >> tmpString;
1843  //std::cout << "Just read '" << tmpString << "'" << std::endl;
1844  unsigned int posInTmpString = 6;
1845 
1846  // Isolate 'n_positions' in a string
1847  //char nPositionsString[tmpString.size()-posInTmpString+1]; // avoid compiler warning
1848  std::string nPositionsString((size_t) (tmpString.size()-posInTmpString+1),' ');
1849  unsigned int posInPositionsString = 0;
1850  do {
1851  queso_require_less_msg(posInTmpString, tmpString.size(), "symbol ',' not found in first line of file");
1852  nPositionsString[posInPositionsString++] = tmpString[posInTmpString++];
1853  } while (tmpString[posInTmpString] != ',');
1854  nPositionsString[posInPositionsString] = '\0';
1855 
1856  // Isolate 'n_params' in a string
1857  posInTmpString++; // Avoid reading ',' char
1858  //char nParamsString[tmpString.size()-posInTmpString+1]; // avoid compiler warning
1859  std::string nParamsString((size_t) (tmpString.size()-posInTmpString+1),' ');
1860  unsigned int posInParamsString = 0;
1861  do {
1862  queso_require_less_msg(posInTmpString, tmpString.size(), "symbol ')' not found in first line of file");
1863  nParamsString[posInParamsString++] = tmpString[posInTmpString++];
1864  } while (tmpString[posInTmpString] != ')');
1865  nParamsString[posInParamsString] = '\0';
1866 
1867  // Convert 'n_positions' and 'n_params' strings to numbers
1868  unsigned int sizeOfChainInFile = (unsigned int) strtod(nPositionsString.c_str(),NULL);
1869  unsigned int numParamsInFile = (unsigned int) strtod(nParamsString.c_str(), NULL);
1870  if (m_env.subDisplayFile()) {
1871  *m_env.subDisplayFile() << "In SequenceOfVectors<V,M>::unifiedReadContents()"
1872  << ": worldRank " << m_env.worldRank()
1873  << ", fullRank " << m_env.fullRank()
1874  << ", sizeOfChainInFile = " << sizeOfChainInFile
1875  << ", numParamsInFile = " << numParamsInFile
1876  << std::endl;
1877  }
1878 
1879  // Check if [size of chain in file] >= [requested unified sequence size]
1880  queso_require_greater_equal_msg(sizeOfChainInFile, unifiedReadSize, "size of chain in file is not big enough");
1881 
1882  // Check if [num params in file] == [num params in current chain]
1883  queso_require_equal_to_msg(numParamsInFile, numParams, "number of parameters of chain in file is different than number of parameters in this chain object");
1884  } // if (r == 0)
1885 
1886  // Code common to any core in 'inter0Comm', including core of rank 0
1887  unsigned int maxCharsPerLine = 64*numParams; // Up to about 60 characters to represent each parameter value
1888 
1889  unsigned int lineId = 0;
1890  while (lineId < idOfMyFirstLine) {
1891  unifiedFilePtrSet.ifsVar->ignore(maxCharsPerLine,'\n');
1892  lineId++;
1893  };
1894 
1895  if (r == 0) {
1896  // Take care of initial part of the first data line,
1897  // which resembles something like 'variable_name = [value1 value2 ...'
1898  std::string tmpString;
1899 
1900  // Read 'variable name' string
1901  *unifiedFilePtrSet.ifsVar >> tmpString;
1902  //std::cout << "Core 0 just read '" << tmpString << "'" << std::endl;
1903 
1904  // Read '=' sign
1905  *unifiedFilePtrSet.ifsVar >> tmpString;
1906  //std::cout << "Core 0 just read '" << tmpString << "'" << std::endl;
1907  queso_require_equal_to_msg(tmpString, std::string("="), std::string("in core 0, string should be the '=' sign"));
1908 
1909  // Take into account the ' [' portion
1910  std::streampos tmpPos = unifiedFilePtrSet.ifsVar->tellg();
1911  unifiedFilePtrSet.ifsVar->seekg(tmpPos+(std::streampos)2);
1912  }
1913 
1914  V tmpVec(m_vectorSpace.zeroVector());
1915  while (lineId <= idOfMyLastLine) {
1916  for (unsigned int i = 0; i < numParams; ++i) {
1917  *unifiedFilePtrSet.ifsVar >> tmpVec[i];
1918  }
1919  this->setPositionValues(lineId - idOfMyFirstLine, tmpVec);
1920  lineId++;
1921  };
1922  }
1923 #ifdef QUESO_HAS_HDF5
1924  else if (fileType == UQ_FILE_EXTENSION_FOR_HDF_FORMAT) {
1925  if (r == 0) {
1926  hid_t dataset = H5Dopen2(unifiedFilePtrSet.h5Var,
1927  "data",
1928  H5P_DEFAULT); // Dataset access property list
1929  hid_t datatype = H5Dget_type(dataset);
1930  H5T_class_t t_class = H5Tget_class(datatype);
1931  queso_require_equal_to_msg(t_class, H5T_FLOAT, "t_class is not H5T_DOUBLE");
1932  hid_t dataspace = H5Dget_space(dataset);
1933  int rank = H5Sget_simple_extent_ndims(dataspace);
1934  queso_require_equal_to_msg(rank, 2, "hdf rank is not 2");
1935  hsize_t dims_in[2];
1936  int status_n;
1937  status_n = H5Sget_simple_extent_dims(dataspace, dims_in, NULL);
1938  if (status_n) {}; // just to remover compiler warning
1939  //std::cout << "In SequenceOfVectors<V,M>::unifiedReadContents()"
1940  // << ": dims_in[0] = " << dims_in[0]
1941  // << ", dims_in[1] = " << dims_in[1]
1942  // << std::endl;
1943  queso_require_equal_to_msg(dims_in[0], numParams, "dims_in[0] is not equal to 'numParams'");
1944  queso_require_greater_equal_msg(dims_in[1], subReadSize, "dims_in[1] is smaller that requested 'subReadSize'");
1945 
1946  struct timeval timevalBegin;
1947  int iRC = UQ_OK_RC;
1948  iRC = gettimeofday(&timevalBegin,NULL);
1949  if (iRC) {}; // just to remover compiler warning
1950 
1951  unsigned int chainSizeIn = (unsigned int) dims_in[1];
1952  //double* dataIn[numParams]; // avoid compiler warning
1953  std::vector<double*> dataIn((size_t) numParams,NULL);
1954  dataIn[0] = (double*) malloc(numParams*chainSizeIn*sizeof(double));
1955  for (unsigned int i = 1; i < numParams; ++i) { // Yes, from '1'
1956  dataIn[i] = dataIn[i-1] + chainSizeIn; // Yes, just 'chainSizeIn', not 'chainSizeIn*sizeof(double)'
1957  }
1958  //std::cout << "In SequenceOfVectors<V,M>::unifiedReadContents(): h5 case, memory allocated" << std::endl;
1959  herr_t status;
1960  status = H5Dread(dataset,
1961  H5T_NATIVE_DOUBLE,
1962  H5S_ALL,
1963  dataspace,
1964  H5P_DEFAULT,
1965  dataIn[0]);
1966  if (status) {}; // just to remover compiler warning
1967  //std::cout << "In SequenceOfVectors<V,M>::unifiedReadContents(): h5 case, data read" << std::endl;
1968  V tmpVec(m_vectorSpace.zeroVector());
1969  for (unsigned int j = 0; j < subReadSize; ++j) { // Yes, 'subReadSize', not 'chainSizeIn'
1970  for (unsigned int i = 0; i < numParams; ++i) {
1971  tmpVec[i] = dataIn[i][j];
1972  }
1973  this->setPositionValues(j, tmpVec);
1974  }
1975 
1976  double readTime = MiscGetEllapsedSeconds(&timevalBegin);
1977  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 2)) {
1978  *m_env.subDisplayFile() << "In SequenceOfVectors<V,M>::unifiedReadContents()"
1979  << ": worldRank " << m_env.worldRank()
1980  << ", fullRank " << m_env.fullRank()
1981  << ", subEnvironment " << m_env.subId()
1982  << ", subRank " << m_env.subRank()
1983  << ", inter0Rank " << m_env.inter0Rank()
1984  << ", fileName = " << fileName
1985  << ", numParams = " << numParams
1986  << ", chainSizeIn = " << chainSizeIn
1987  << ", subReadSize = " << subReadSize
1988  << ", readTime = " << readTime << " seconds"
1989  << std::endl;
1990  }
1991 
1992  H5Sclose(dataspace);
1993  H5Tclose(datatype);
1994  H5Dclose(dataset);
1995  //free(dataIn[0]); // related to the changes above for compiler warning
1996  for (unsigned int tmpIndex = 0; tmpIndex < dataIn.size(); tmpIndex++) {
1997  free (dataIn[tmpIndex]);
1998  }
1999  }
2000  else {
2001  queso_error_msg("hdf file type not supported for multiple sub-environments yet");
2002  }
2003  }
2004 #endif
2005  else {
2006  queso_error_msg("invalid file type");
2007  }
2008  m_env.closeFile(unifiedFilePtrSet,fileType);
2009  } // if (m_env.openUnifiedInputFile())
2010  } // if (m_env.inter0Rank() == (int) r)
2011  m_env.inter0Comm().Barrier();
2012  } // for r
2013  } // if (m_env.inter0Rank() >= 0)
2014  else {
2015  V tmpVec(m_vectorSpace.zeroVector());
2016  for (unsigned int i = 1; i < subReadSize; ++i) {
2017  this->setPositionValues(i,tmpVec);
2018  }
2019  }
2020 
2021  if (m_env.subDisplayFile()) {
2022  *m_env.subDisplayFile() << "Leaving SequenceOfVectors<V,M>::unifiedReadContents()"
2023  << ", fileName = " << fileName
2024  << std::endl;
2025  }
2026  //m_env.fullComm().Barrier(); // Dangerous to barrier on fullComm ...
2027 
2028  return;
2029 }
2030 //---------------------------------------------------
2031 template <class V, class M>
2032 void
2033 SequenceOfVectors<V,M>::select(const std::vector<unsigned int>& /* idsOfUniquePositions */)
2034 {
2035  queso_error_msg("Code is not complete yet");
2036 
2037  return;
2038 }
2039 //---------------------------------------------------
2040 template <class V, class M>
2041 void
2042 SequenceOfVectors<V,M>::filter(
2043  unsigned int initialPos,
2044  unsigned int spacing)
2045 {
2046  if (m_env.subDisplayFile()) {
2047  *m_env.subDisplayFile() << "Entering SequenceOfVectors<V,M>::filter()"
2048  << ": initialPos = " << initialPos
2049  << ", spacing = " << spacing
2050  << ", subSequenceSize = " << this->subSequenceSize()
2051  << std::endl;
2052  }
2053 
2054  unsigned int i = 0;
2055  unsigned int j = initialPos;
2056  unsigned int originalSubSequenceSize = this->subSequenceSize();
2057  while (j < originalSubSequenceSize) {
2058  if (i != j) {
2059  //*m_env.subDisplayFile() << i << "--" << j << " ";
2060  delete m_seq[i];
2061  m_seq[i] = new V(*(m_seq[j]));
2062  }
2063  i++;
2064  j += spacing;
2065  }
2066 
2067  this->resetValues(i,originalSubSequenceSize-i);
2068  this->resizeSequence(i);
2069 
2070  if (m_env.subDisplayFile()) {
2071  *m_env.subDisplayFile() << "Leaving SequenceOfVectors<V,M>::filter()"
2072  << ": initialPos = " << initialPos
2073  << ", spacing = " << spacing
2074  << ", subSequenceSize = " << this->subSequenceSize()
2075  << std::endl;
2076  }
2077 
2078  return;
2079 }
2080 //---------------------------------------------------
2081 template <class V, class M>
2082 double
2083 SequenceOfVectors<V,M>::estimateConvBrooksGelman(
2084  unsigned int initialPos,
2085  unsigned int numPos) const
2086 {
2087  // This method requires *at least* two sequences. Error if there is only one.
2088  queso_require_greater_equal_msg(m_env.numSubEnvironments(), 2, "At least two sequences required for Brooks-Gelman convergence test.");
2089 
2090  // TODO: Need special case for 1-dimensional parameter space.
2091 
2092  // Initialize with garbage to give the user a clue something is funky.
2093  double convMeasure = -1.0;
2094 
2095  // We only do the work on the subenvironment where the sequence data
2096  // is stashed.
2097  if( m_env.inter0Rank() >= 0 )
2098  {
2099  // Sanity Checks
2100 
2101  // REMEMBER: \psi is a *vector* of parameters
2102  // Get quantities we will use several times
2103  V psi_j_dot = m_vectorSpace.zeroVector();
2104  V psi_dot_dot = m_vectorSpace.zeroVector();
2105  V work = m_vectorSpace.zeroVector();
2106 
2107  // m = number of chains > 1
2108  // n = number of steps for which we are computing the metric
2109  int m = m_env.numSubEnvironments();
2110  int n = numPos;
2111 
2112  this->subMeanExtra ( initialPos, numPos, psi_j_dot );
2113  this->unifiedMeanExtra( initialPos, numPos, psi_dot_dot );
2114 
2115 #if 0
2116  std::cout << "psi_j_dot = " << psi_j_dot << std::endl;
2117  std::cout << "psi_dot_dot = " << psi_dot_dot << std::endl;
2118 #endif
2119 
2120  /* W = \frac{1}{m*(n-1)}*\sum_{j=1}^m \sum{t=1}^n
2121  (\psi_{jt} - \overline{\psi_{j\cdot}})*(\psi_{jt} - \overline{\psi_{j\cdot}})^T
2122  This corresponds to the "within-sequence" covariance matrix. */
2123  M* W_local = m_vectorSpace.newDiagMatrix( m_vectorSpace.zeroVector() );
2124  M* W = m_vectorSpace.newDiagMatrix( m_vectorSpace.zeroVector() );
2125  V psi_j_t = m_vectorSpace.zeroVector();
2126 
2127  // Sum within the chain
2128  for( unsigned int t = initialPos; t < initialPos+numPos; ++t )
2129  {
2130  psi_j_t = *(m_seq[t]);
2131 
2132  work = psi_j_t - psi_j_dot;
2133 
2134  (*W_local) += matrixProduct( work, work );
2135  }
2136 
2137  // Now do the sum over the chains
2138  // W will be available on all inter0 processors
2139  W_local->mpiSum( m_env.inter0Comm(), (*W) );
2140 
2141  (*W) = 1.0/(double(m)*(double(n)-1.0)) * (*W);
2142 
2143 #if 0
2144  std::cout << "n, m = " << n << ", " << m << std::endl;
2145  std::cout << "W_local = " << *W_local << std::endl;
2146  std::cout << "W = " << *W << std::endl;
2147 #endif
2148 
2149  // Need to delete pointers to temporary covariance matrices
2150  delete W_local;
2151 
2152  /* B/n = \frac{1}{m-1}\sum_{j=1}^m
2153  (\overline{\psi_{j\cdot}} - \overline{\psi_{\cdot \cdot}})*
2154  (\overline{\psi_{j\cdot}} - \overline{\psi_{\cdot \cdot}})^T
2155  This corresponds to the "between-sequence" covariance matrix. */
2156  M* B_over_n_local = m_vectorSpace.newDiagMatrix( m_vectorSpace.zeroVector() );
2157  M* B_over_n = m_vectorSpace.newDiagMatrix( m_vectorSpace.zeroVector() );
2158 
2159  work = psi_j_dot - psi_dot_dot;
2160  (*B_over_n_local) = matrixProduct( work, work );
2161 
2162  B_over_n_local->mpiSum( m_env.inter0Comm(), (*B_over_n) );
2163 
2164  // Need to delete pointers to temporary covariance matrices
2165  delete B_over_n_local;
2166 
2167  (*B_over_n) = 1.0/(double(m)-1.0) * (*B_over_n);
2168 
2169 #if 0
2170  std::cout << "B_over_n = " << *B_over_n << std::endl;
2171 #endif
2172 
2173 
2174  /* R_p = (n-1)/n + (m+1)/m * \lambda
2175  \lambda = largest eigenvalue of W^{-1}*B/n */
2176  M* A = m_vectorSpace.newDiagMatrix( m_vectorSpace.zeroVector() );
2177 
2178  W->invertMultiply( *B_over_n, *A );
2179 
2180 #if 0
2181  std::cout << "A = " << *A << std::endl;
2182  std::cout.flush();
2183 #endif
2184  // Need to delete pointers to temporary covariance matrices
2185  delete W;
2186  delete B_over_n;
2187 
2188  double eigenValue;
2189  V eigenVector = m_vectorSpace.zeroVector();
2190 
2191  A->largestEigen( eigenValue, eigenVector );
2192 
2193  // Need to delete pointers to temporary covariance matrices
2194  delete A;
2195 
2196  // Now, finally compute the final convMeasure
2197  convMeasure = (double(n)-1.0)/double(n) + (double(m)+1.0)/double(m)*eigenValue;
2198 
2199  } // End of check on inter0Rank
2200 
2201  //TODO: Do we need a Barrier here?
2202  //TODO: Error checking on MPI.
2203  //m_env.fullComm().Barrier(); // Dangerous to barrier on fullComm ...
2204 
2205  return convMeasure;
2206 }
2207 //---------------------------------------------------
2208 template <class V, class M>
2209 void
2210 SequenceOfVectors<V,M>::extractScalarSeq(
2211  unsigned int initialPos,
2212  unsigned int spacing,
2213  unsigned int numPos,
2214  unsigned int paramId,
2215  ScalarSequence<double>& scalarSeq) const
2216 {
2217  scalarSeq.resizeSequence(numPos);
2218  if (spacing == 1) {
2219  for (unsigned int j = 0; j < numPos; ++j) {
2220  scalarSeq[j] = (*(m_seq[initialPos+j ]))[paramId];
2221  }
2222  }
2223  else {
2224  for (unsigned int j = 0; j < numPos; ++j) {
2225  scalarSeq[j] = (*(m_seq[initialPos+j*spacing]))[paramId];
2226  }
2227  }
2228 
2229  return;
2230 }
2231 // Private methods ------------------------------------
2232 template <class V, class M>
2233 void
2234 SequenceOfVectors<V,M>::copy(const SequenceOfVectors<V,M>& src)
2235 {
2236  BaseVectorSequence<V,M>::copy(src);
2237  for (unsigned int i = 0; i < (unsigned int) m_seq.size(); ++i) {
2238  if (m_seq[i]) {
2239  delete m_seq[i];
2240  m_seq[i] = NULL;
2241  }
2242  }
2243  m_seq.resize(src.subSequenceSize(),NULL);
2244  for (unsigned int i = 0; i < m_seq.size(); ++i) {
2245  m_seq[i] = new V(*(src.m_seq[i]));
2246  }
2247 
2248  return;
2249 }
2250 //---------------------------------------------------
2251 template <class V, class M>
2252 void
2253 SequenceOfVectors<V,M>::extractRawData(
2254  unsigned int initialPos,
2255  unsigned int spacing,
2256  unsigned int numPos,
2257  unsigned int paramId,
2258  std::vector<double>& rawData) const
2259 {
2260  rawData.resize(numPos);
2261  if (spacing == 1) {
2262  for (unsigned int j = 0; j < numPos; ++j) {
2263  rawData[j] = (*(m_seq[initialPos+j ]))[paramId];
2264  }
2265  }
2266  else {
2267  for (unsigned int j = 0; j < numPos; ++j) {
2268  rawData[j] = (*(m_seq[initialPos+j*spacing]))[paramId];
2269  }
2270  }
2271 
2272  return;
2273 }
2274 
2275 // --------------------------------------------------
2276 // Methods conditionally available ------------------
2277 // --------------------------------------------------
2278 // --------------------------------------------------
2279 #ifdef UQ_SEQ_VEC_USES_OPERATOR
2280 template <class V, class M>
2281 const V*
2282 SequenceOfVectors<V,M>::operator[](unsigned int posId) const
2283 {
2284  queso_require_less_msg(posId, this->subSequenceSize(), "posId > subSequenceSize()");
2285 
2286  return (const V*) (m_seq[posId]);
2287 }
2288 // --------------------------------------------------
2289 template <class V, class M>
2290 const V*&
2291 SequenceOfVectors<V,M>::operator[](unsigned int posId)
2292 {
2293  queso_require_less_msg(posId, this->subSequenceSize(), "posId > subSequenceSize()");
2294 
2295  return m_seq[posId];
2296 }
2297 #endif
2298 
2299 // --------------------------------------------------
2300 // --------------------------------------------------
2301 // --------------------------------------------------
2302 
2303 #ifdef UQ_ALSO_COMPUTE_MDFS_WITHOUT_KDE
2304 template <class V, class M>
2305 void
2306 SequenceOfVectors<V,M>::subUniformlySampledMdf(
2307  const V& numEvaluationPointsVec,
2308  ArrayOfOneDGrids <V,M>& mdfGrids,
2309  ArrayOfOneDTables<V,M>& mdfValues) const
2310 {
2311  V minDomainValues(m_vectorSpace.zeroVector());
2312  V maxDomainValues(m_vectorSpace.zeroVector());
2313 
2314  ScalarSequence<double> data(m_env,0,"");
2315 
2316  unsigned int numParams = this->vectorSizeLocal();
2317  for (unsigned int i = 0; i < numParams; ++i) {
2318  this->extractScalarSeq(0, // initialPos
2319  1, // spacing
2320  subSequenceSize(), // numPos
2321  i,
2322  data);
2323 
2324  std::vector<double> aMdf(0);
2325  data.subUniformlySampledMdf((unsigned int) numEvaluationPointsVec[i],
2326  minDomainValues[i],
2327  maxDomainValues[i],
2328  aMdf);
2329  mdfValues.setOneDTable(i,aMdf);
2330  }
2331 
2332  mdfGrids.setUniformGrids(numEvaluationPointsVec,
2333  minDomainValues,
2334  maxDomainValues);
2335 
2336  return;
2337 }
2338 #endif
2339 
2340 // --------------------------------------------------
2341 // --------------------------------------------------
2342 // --------------------------------------------------
2343 
2344 #ifdef QUESO_COMPUTES_EXTRA_POST_PROCESSING_STATISTICS
2345 template <class V, class M>
2346 void
2347 SequenceOfVectors<V,M>::subMeanCltStd(
2348  unsigned int initialPos,
2349  unsigned int numPos,
2350  const V& meanVec,
2351  V& stdVec) const
2352 {
2353  bool bRC = ((initialPos < this->subSequenceSize()) &&
2354  (0 < numPos ) &&
2355  ((initialPos+numPos) <= this->subSequenceSize()) &&
2356  (this->vectorSizeLocal() == meanVec.sizeLocal() ) &&
2357  (this->vectorSizeLocal() == stdVec.sizeLocal() ));
2358  queso_require_msg(bRC, "invalid input data");
2359 
2360  ScalarSequence<double> data(m_env,0,"");
2361 
2362  unsigned int numParams = this->vectorSizeLocal();
2363  for (unsigned int i = 0; i < numParams; ++i) {
2364  this->extractScalarSeq(initialPos,
2365  1, // spacing
2366  numPos,
2367  i,
2368  data);
2369  stdVec[i] = data.subMeanCltStd(0,
2370  numPos,
2371  meanVec[i]);
2372  }
2373 
2374  return;
2375 }
2376 
2377 template <class V, class M>
2378 void
2379 SequenceOfVectors<V,M>::unifiedMeanCltStd(
2380  unsigned int initialPos,
2381  unsigned int numPos,
2382  const V& unifiedMeanVec,
2383  V& unifiedSamVec) const
2384 {
2385  bool bRC = ((initialPos < this->subSequenceSize() ) &&
2386  (0 < numPos ) &&
2387  ((initialPos+numPos) <= this->subSequenceSize() ) &&
2388  (this->vectorSizeLocal() == unifiedMeanVec.sizeLocal()) &&
2389  (this->vectorSizeLocal() == unifiedSamVec.sizeLocal() ));
2390  queso_require_msg(bRC, "invalid input data");
2391 
2392  ScalarSequence<double> data(m_env,0,"");
2393 
2394  unsigned int numParams = this->vectorSizeLocal();
2395  for (unsigned int i = 0; i < numParams; ++i) {
2396  this->extractScalarSeq(initialPos,
2397  1, // spacing
2398  numPos,
2399  i,
2400  data);
2401  unifiedSamVec[i] = data.unifiedMeanCltStd(m_vectorSpace.numOfProcsForStorage() == 1,
2402  0,
2403  numPos,
2404  unifiedMeanVec[i]);
2405  }
2406 
2407  return;
2408 }
2409 //---------------------------------------------------
2410 template <class V, class M>
2411 void
2412 SequenceOfVectors<V,M>::bmm(
2413  unsigned int initialPos,
2414  unsigned int batchLength,
2415  V& bmmVec) const
2416 {
2417  bool bRC = ((initialPos < this->subSequenceSize() ) &&
2418  (batchLength < (this->subSequenceSize()-initialPos)) &&
2419  (this->vectorSizeLocal() == bmmVec.sizeLocal() ));
2420  queso_require_msg(bRC, "invalid input data");
2421 
2422  ScalarSequence<double> data(m_env,0,"");
2423 
2424  unsigned int numParams = this->vectorSizeLocal();
2425  for (unsigned int i = 0; i < numParams; ++i) {
2426  this->extractScalarSeq(initialPos,
2427  1, // spacing
2428  this->subSequenceSize()-initialPos,
2429  i,
2430  data);
2431  bmmVec[i] = data.bmm(0,
2432  batchLength);
2433  }
2434 
2435  return;
2436 }
2437 //---------------------------------------------------
2438 template <class V, class M>
2439 void
2440 SequenceOfVectors<V,M>::fftForward(
2441  unsigned int initialPos,
2442  unsigned int fftSize,
2443  unsigned int paramId,
2444  std::vector<std::complex<double> >& fftResult) const
2445 {
2446  bool bRC = ((initialPos < this->subSequenceSize()) &&
2447  (paramId < this->vectorSizeLocal()) &&
2448  (0 < fftSize ) &&
2449  ((initialPos+fftSize) <= this->subSequenceSize()) &&
2450  (fftSize < this->subSequenceSize()));
2451  queso_require_msg(bRC, "invalid input data");
2452 
2453  std::vector<double> rawData(fftSize,0.);
2454  this->extractRawData(initialPos,
2455  1, // spacing
2456  fftSize,
2457  paramId,
2458  rawData);
2459 
2460  m_fftObj->forward(rawData,fftSize,fftResult);
2461 
2462  return;
2463 }
2464 //---------------------------------------------------
2465 template <class V, class M>
2466 void
2467 SequenceOfVectors<V,M>::psd(
2468  unsigned int initialPos,
2469  unsigned int numBlocks,
2470  double hopSizeRatio,
2471  unsigned int paramId,
2472  std::vector<double>& psdResult) const
2473 {
2474  bool bRC = ((initialPos < this->subSequenceSize()) &&
2475  (paramId < this->vectorSizeLocal()));
2476  queso_require_msg(bRC, "invalid input data");
2477 
2478  ScalarSequence<double> data(m_env,0,"");
2479 
2480  this->extractScalarSeq(initialPos,
2481  1, // spacing
2482  this->subSequenceSize()-initialPos,
2483  paramId,
2484  data);
2485  data.psd(0,
2486  numBlocks,
2487  hopSizeRatio,
2488  psdResult);
2489 
2490  return;
2491 }
2492 //---------------------------------------------------
2493 template <class V, class M>
2494 void
2495 SequenceOfVectors<V,M>::psdAtZero(
2496  unsigned int initialPos,
2497  unsigned int numBlocks,
2498  double hopSizeRatio,
2499  V& psdVec) const
2500 {
2501  bool bRC = ((initialPos < this->subSequenceSize()) &&
2502  (this->vectorSizeLocal() == psdVec.sizeLocal()));
2503  queso_require_msg(bRC, "invalid input data");
2504 
2505  ScalarSequence<double> data(m_env,0,"");
2506  std::vector<double> psdResult(0,0.); // size will be determined by 'data.psd()'
2507 
2508  unsigned int numParams = this->vectorSizeLocal();
2509  for (unsigned int i = 0; i < numParams; ++i) {
2510  this->extractScalarSeq(initialPos,
2511  1, // spacing
2512  this->subSequenceSize()-initialPos,
2513  i,
2514  data);
2515  data.psd(0,
2516  numBlocks,
2517  hopSizeRatio,
2518  psdResult);
2519  psdVec[i] = psdResult[0];
2520  //*m_env.subDisplayFile() << "psdResult[0] = " << psdResult[0] << std::endl;
2521  }
2522 
2523  return;
2524 }
2525 //---------------------------------------------------
2526 template <class V, class M>
2527 void
2528 SequenceOfVectors<V,M>::geweke(
2529  unsigned int initialPos,
2530  double ratioNa,
2531  double ratioNb,
2532  V& gewVec) const
2533 {
2534  bool bRC = ((initialPos < this->subSequenceSize()) &&
2535  (this->vectorSizeLocal() == gewVec.sizeLocal() ));
2536  queso_require_msg(bRC, "invalid input data");
2537 
2538  unsigned int numPos = this->subSequenceSize() - initialPos;
2539  ScalarSequence<double> data(m_env,0,"");
2540 
2541  unsigned int numParams = this->vectorSizeLocal();
2542  for (unsigned int i = 0; i < numParams; ++i) {
2543  this->extractScalarSeq(initialPos,
2544  1, // spacing
2545  numPos,
2546  i,
2547  data);
2548  gewVec[i] = data.geweke(0,
2549  ratioNa,
2550  ratioNb);
2551  }
2552 
2553  return;
2554 }
2555 //---------------------------------------------------
2556 template <class V, class M>
2557 void
2558 SequenceOfVectors<V,M>::meanStacc(
2559  unsigned int initialPos,
2560  V& meanStaccVec) const
2561 {
2562  bool bRC = ((initialPos < this->subSequenceSize() ) &&
2563  (this->vectorSizeLocal() == meanStaccVec.sizeLocal()));
2564  queso_require_msg(bRC, "invalid input data");
2565 
2566  unsigned int numPos = this->subSequenceSize() - initialPos;
2567  ScalarSequence<double> data(m_env,0,"");
2568 
2569  unsigned int numParams = this->vectorSizeLocal();
2570  for (unsigned int i = 0; i < numParams; ++i) {
2571  this->extractScalarSeq(initialPos,
2572  1, // spacing
2573  numPos,
2574  i,
2575  data);
2576  meanStaccVec[i] = data.meanStacc(0);
2577  }
2578 
2579  return;
2580 }
2581 //---------------------------------------------------
2582 template <class V, class M>
2583 void
2584 SequenceOfVectors<V,M>::subCdfPercentageRange(
2585  unsigned int initialPos,
2586  unsigned int numPos,
2587  double range, // \in [0,1]
2588  V& lowerVec,
2589  V& upperVec) const
2590 {
2591  bool bRC = ((0 < numPos ) &&
2592  ((initialPos+numPos) <= this->subSequenceSize()) &&
2593  (this->vectorSizeLocal() == lowerVec.sizeLocal() ) &&
2594  (this->vectorSizeLocal() == upperVec.sizeLocal() ));
2595  queso_require_msg(bRC, "invalid input data");
2596 
2597  unsigned int numParams = this->vectorSizeLocal();
2598  ScalarSequence<double> data(m_env,0,"");
2599 
2600  for (unsigned int i = 0; i < numParams; ++i) {
2601  this->extractScalarSeq(initialPos,
2602  1, // spacing
2603  numPos,
2604  i,
2605  data);
2606  data.subCdfPercentageRange(0,
2607  numPos,
2608  range,
2609  lowerVec[i],
2610  upperVec[i]);
2611  }
2612 
2613  return;
2614 }
2615 //---------------------------------------------------
2616 template <class V, class M>
2617 void
2618 SequenceOfVectors<V,M>::unifiedCdfPercentageRange(
2619  unsigned int initialPos,
2620  unsigned int numPos,
2621  double range, // \in [0,1]
2622  V& lowerVec,
2623  V& upperVec) const
2624 {
2625  bool bRC = ((0 < numPos ) &&
2626  ((initialPos+numPos) <= this->subSequenceSize()) &&
2627  (this->vectorSizeLocal() == lowerVec.sizeLocal() ) &&
2628  (this->vectorSizeLocal() == upperVec.sizeLocal() ));
2629  queso_require_msg(bRC, "invalid input data");
2630 
2631  unsigned int numParams = this->vectorSizeLocal();
2632  ScalarSequence<double> data(m_env,0,"");
2633 
2634  for (unsigned int i = 0; i < numParams; ++i) {
2635  this->extractScalarSeq(initialPos,
2636  1, // spacing
2637  numPos,
2638  i,
2639  data);
2640  data.unifiedCdfPercentageRange(m_vectorSpace.numOfProcsForStorage() == 1,
2641  0,
2642  numPos,
2643  range,
2644  lowerVec[i],
2645  upperVec[i]);
2646  }
2647 
2648  return;
2649 }
2650 //---------------------------------------------------
2651 template <class V, class M>
2652 void
2653 SequenceOfVectors<V,M>::subCdfStacc(
2654  unsigned int initialPos,
2655  std::vector<V*>& cdfStaccVecs,
2656  std::vector<V*>& cdfStaccVecsUp,
2657  std::vector<V*>& cdfStaccVecsLow,
2658  std::vector<V*>& sortedDataVecs) const
2659 {
2660  bool bRC = (initialPos < this->subSequenceSize());
2661  queso_require_msg(bRC, "invalid input data");
2662 
2663  unsigned int numPos = this->subSequenceSize() - initialPos;
2664  unsigned int numEvals = numPos;
2665  for (unsigned int j = 0; j < numEvals; ++j) {
2666  cdfStaccVecs [j] = new V(m_vectorSpace.zeroVector());
2667  cdfStaccVecsUp [j] = new V(m_vectorSpace.zeroVector());
2668  cdfStaccVecsLow[j] = new V(m_vectorSpace.zeroVector());
2669  sortedDataVecs [j] = new V(m_vectorSpace.zeroVector());
2670  }
2671  std::vector<double> cdfStaccs (numEvals,0.);
2672  std::vector<double> cdfStaccsup (numEvals,0.);
2673  std::vector<double> cdfStaccslow(numEvals,0.);
2674 
2675  ScalarSequence<double> data (m_env,0,"");
2676  ScalarSequence<double> sortedData(m_env,0,"");
2677  unsigned int numParams = this->vectorSizeLocal();
2678  for (unsigned int i = 0; i < numParams; ++i) {
2679  this->extractScalarSeq(initialPos,
2680  1, // spacing
2681  numPos,
2682  i,
2683  data);
2684  //std::cout << "x-data" << data<< std::endl;
2685  data.subSort(initialPos,sortedData);
2686  data.subCdfStacc(initialPos,
2687  cdfStaccs,
2688  cdfStaccsup,
2689  cdfStaccslow,
2690  sortedData);
2691 
2692  for (unsigned int j = 0; j < numEvals; ++j) {
2693  (*sortedDataVecs [j])[i] = sortedData [j];
2694  (*cdfStaccVecs [j])[i] = cdfStaccs [j];
2695  (*cdfStaccVecsUp [j])[i] = cdfStaccsup [j];
2696  (*cdfStaccVecsLow[j])[i] = cdfStaccslow[j];
2697  }
2698  }
2699 
2700  return;
2701 }
2702 //---------------------------------------------------
2703 template <class V, class M>
2704 void
2705 SequenceOfVectors<V,M>::subCdfStacc(
2706  unsigned int initialPos,
2707  const std::vector<V*>& evalPositionsVecs,
2708  std::vector<V*>& cdfStaccVecs) const
2709 {
2710  bool bRC = ((initialPos < this->subSequenceSize() ) &&
2711  (0 < evalPositionsVecs.size()) &&
2712  (evalPositionsVecs.size() == cdfStaccVecs.size() ));
2713  queso_require_msg(bRC, "invalid input data");
2714 
2715  unsigned int numPos = this->subSequenceSize() - initialPos;
2716  ScalarSequence<double> data(m_env,0,"");
2717 
2718  unsigned int numEvals = evalPositionsVecs.size();
2719  for (unsigned int j = 0; j < numEvals; ++j) {
2720  cdfStaccVecs[j] = new V(m_vectorSpace.zeroVector());
2721  }
2722  std::vector<double> evalPositions(numEvals,0.);
2723  std::vector<double> cdfStaccs (numEvals,0.);
2724 
2725  unsigned int numParams = this->vectorSizeLocal();
2726  for (unsigned int i = 0; i < numParams; ++i) {
2727  this->extractScalarSeq(initialPos,
2728  1, // spacing
2729  numPos,
2730  i,
2731  data);
2732 
2733  for (unsigned int j = 0; j < numEvals; ++j) {
2734  evalPositions[j] = (*evalPositionsVecs[j])[i];
2735  }
2736 
2737  data.subCdfStacc(0,
2738  evalPositions,
2739  cdfStaccs);
2740 
2741  for (unsigned int j = 0; j < numEvals; ++j) {
2742  (*cdfStaccVecs[j])[i] = cdfStaccs[j];
2743  }
2744  }
2745 
2746  return;
2747 }
2748 #endif // #ifdef QUESO_COMPUTES_EXTRA_POST_PROCESSING_STATISTICS
2749 
2750 // --------------------------------------------------
2751 // --------------------------------------------------
2752 // --------------------------------------------------
2753 
2754 #ifdef UQ_CODE_HAS_MONITORS
2755 template <class V, class M>
2756 void
2757 SequenceOfVectors<V,M>::subMeanMonitorAlloc(unsigned int numberOfMonitorPositions)
2758 {
2759  m_subMeanMonitorPosSeq = new ScalarSequence<double>(m_env, numberOfMonitorPositions,(m_name+"_subMeanMonitorPosSeq").c_str());
2760  m_subMeanVecSeq = new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+"_subMeanVecSeq").c_str() );
2761  m_subMeanCltStdSeq = new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+"_subMeanCltStdSeq").c_str() );
2762 
2763  return;
2764 }
2765 // --------------------------------------------------
2766 template <class V, class M>
2767 void
2768 SequenceOfVectors<V,M>::subMeanInter0MonitorAlloc(unsigned int numberOfMonitorPositions)
2769 {
2770  m_subMeanInter0MonitorPosSeq = new ScalarSequence<double>(m_env, numberOfMonitorPositions,(m_name+"_subMeanInter0MonitorPosSeq").c_str() );
2771  m_subMeanInter0Mean = new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+"_subMeanInter0MeanSeq").c_str() );
2772  m_subMeanInter0Clt95 = new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+"_subMeanInter0Clt95Seq").c_str() );
2773  m_subMeanInter0Empirical90 = new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+"_subMeanInter0Empirical90Seq").c_str());
2774  m_subMeanInter0Min = new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+"_subMeanInter0MinSeq").c_str() );
2775  m_subMeanInter0Max = new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+"_subMeanInter0MaxSeq").c_str() );
2776 
2777  return;
2778 }
2779 // --------------------------------------------------
2780 template <class V, class M>
2781 void
2782 SequenceOfVectors<V,M>::unifiedMeanMonitorAlloc(unsigned int numberOfMonitorPositions)
2783 {
2784  m_unifiedMeanMonitorPosSeq = new ScalarSequence<double>(m_env, numberOfMonitorPositions,(m_name+"_unifiedMeanMonitorPosSeq").c_str());
2785  m_unifiedMeanVecSeq = new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+"_unifiedMeanVecSeq").c_str() );
2786  m_unifiedMeanCltStdSeq = new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+"_unifiedMeanCltStdSeq").c_str() );
2787 
2788  return;
2789 }
2790 // --------------------------------------------------
2791 template <class V, class M>
2792 void
2793 SequenceOfVectors<V,M>::subMeanMonitorRun(unsigned int monitorPosition,
2794  V& subMeanVec,
2795  V& subMeanCltStd)
2796 {
2797  this->subMeanExtra(0,
2798  monitorPosition,
2799  subMeanVec);
2800 
2801  this->subMeanCltStd(0,
2802  monitorPosition,
2803  subMeanVec,
2804  subMeanCltStd);
2805 
2806  return;
2807 }
2808 // --------------------------------------------------
2809 template <class V, class M>
2810 void
2811 SequenceOfVectors<V,M>::subMeanInter0MonitorRun(unsigned int monitorPosition,
2812  V& subMeanInter0Mean,
2813  V& subMeanInter0Clt95,
2814  V& subMeanInter0Empirical90,
2815  V& subMeanInter0Min,
2816  V& subMeanInter0Max)
2817 {
2818  V subMeanVec(m_vectorSpace.zeroVector());
2819  this->subMeanExtra(0,
2820  monitorPosition,
2821  subMeanVec);
2822 
2823  subMeanVec.mpiAllReduce(RawValue_MPI_SUM,m_env.inter0Comm(),subMeanInter0Mean);
2824  subMeanInter0Mean /= ((double) m_env.inter0Comm().NumProc());
2825 
2826  V subMeanInter0CltVariance = subMeanVec-subMeanInter0Mean;
2827  subMeanInter0CltVariance *= subMeanInter0CltVariance;
2828  subMeanInter0CltVariance.mpiAllReduce(RawValue_MPI_SUM,m_env.inter0Comm(),subMeanInter0Clt95);
2829  subMeanInter0Clt95 /= ((double) (m_env.inter0Comm().NumProc()-1));
2830  subMeanInter0Clt95 /= ((double) (m_env.inter0Comm().NumProc()-1));
2831  subMeanInter0Clt95.cwSqrt();
2832  subMeanInter0Clt95 *= 3.;
2833 
2834  V subMeanInter0Quantile5(m_vectorSpace.zeroVector());
2835  subMeanVec.mpiAllQuantile(.05,m_env.inter0Comm(),subMeanInter0Quantile5);
2836  V subMeanInter0Quantile95(m_vectorSpace.zeroVector());
2837  subMeanVec.mpiAllQuantile(.95,m_env.inter0Comm(),subMeanInter0Quantile95);
2838  subMeanInter0Empirical90 = subMeanInter0Quantile95 - subMeanInter0Quantile5;
2839 
2840  subMeanVec.mpiAllReduce(RawValue_MPI_MIN,m_env.inter0Comm(),subMeanInter0Min);
2841 
2842  subMeanVec.mpiAllReduce(RawValue_MPI_MAX,m_env.inter0Comm(),subMeanInter0Max);
2843 
2844  return;
2845 }
2846 // --------------------------------------------------
2847 template <class V, class M>
2848 void
2849 SequenceOfVectors<V,M>::unifiedMeanMonitorRun(unsigned int monitorPosition,
2850  V& unifiedMeanVec,
2851  V& unifiedMeanCltStd)
2852 {
2853  this->unifiedMeanExtra(0,
2854  monitorPosition,
2855  unifiedMeanVec);
2856 
2857  this->unifiedMeanCltStd(0,
2858  monitorPosition,
2859  unifiedMeanVec,
2860  unifiedMeanCltStd);
2861  return;
2862 }
2863 // --------------------------------------------------
2864 template <class V, class M>
2865 void
2866 SequenceOfVectors<V,M>::subMeanMonitorStore(unsigned int i,
2867  unsigned int monitorPosition,
2868  const V& subMeanVec,
2869  const V& subMeanCltStd)
2870 {
2871  (*m_subMeanMonitorPosSeq)[i] = monitorPosition;
2872  m_subMeanVecSeq->setPositionValues(i,subMeanVec);
2873  m_subMeanCltStdSeq->setPositionValues(i,subMeanCltStd);
2874 
2875  return;
2876 }
2877 // --------------------------------------------------
2878 template <class V, class M>
2879 void
2880 SequenceOfVectors<V,M>::subMeanInter0MonitorStore(unsigned int i,
2881  unsigned int monitorPosition,
2882  const V& subMeanInter0Mean,
2883  const V& subMeanInter0Clt95,
2884  const V& subMeanInter0Empirical90,
2885  const V& subMeanInter0Min,
2886  const V& subMeanInter0Max)
2887 {
2888  (*m_subMeanInter0MonitorPosSeq)[i] = monitorPosition;
2889  m_subMeanInter0Mean->setPositionValues(i,subMeanInter0Mean);
2890  m_subMeanInter0Clt95->setPositionValues(i,subMeanInter0Clt95);
2891  m_subMeanInter0Empirical90->setPositionValues(i,subMeanInter0Empirical90);
2892  m_subMeanInter0Min->setPositionValues(i,subMeanInter0Min);
2893  m_subMeanInter0Max->setPositionValues(i,subMeanInter0Max);
2894 
2895  return;
2896 }
2897 // --------------------------------------------------
2898 template <class V, class M>
2899 void
2900 SequenceOfVectors<V,M>::unifiedMeanMonitorStore(unsigned int i,
2901  unsigned int monitorPosition,
2902  V& unifiedMeanVec,
2903  V& unifiedMeanCltStd)
2904 {
2905  (*m_unifiedMeanMonitorPosSeq)[i] = monitorPosition;
2906  m_unifiedMeanVecSeq->setPositionValues(i,unifiedMeanVec);
2907  m_unifiedMeanCltStdSeq->setPositionValues(i,unifiedMeanCltStd);
2908 
2909  return;
2910 }
2911 // --------------------------------------------------
2912 template <class V, class M>
2913 void
2914 SequenceOfVectors<V,M>::subMeanMonitorWrite(std::ofstream& ofs)
2915 {
2916  m_subMeanMonitorPosSeq->subWriteContents(0,m_subMeanMonitorPosSeq->subSequenceSize(),ofs,UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT); // Yes, always ".m"
2917  m_subMeanVecSeq->subWriteContents (0,m_subMeanVecSeq->subSequenceSize(), ofs,UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT); // Yes, always ".m"
2918  m_subMeanCltStdSeq->subWriteContents (0,m_subMeanCltStdSeq->subSequenceSize(), ofs,UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT); // Yes, always ".m"
2919 
2920  return;
2921 }
2922 // --------------------------------------------------
2923 template <class V, class M>
2924 void
2925 SequenceOfVectors<V,M>::subMeanInter0MonitorWrite(std::ofstream& ofs)
2926 {
2927  m_subMeanInter0MonitorPosSeq->subWriteContents(0,m_subMeanInter0MonitorPosSeq->subSequenceSize(),ofs,UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT); // Yes, always ".m"
2928  m_subMeanInter0Mean->subWriteContents (0,m_subMeanInter0Mean->subSequenceSize(), ofs,UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT); // Yes, always ".m"
2929  m_subMeanInter0Clt95->subWriteContents (0,m_subMeanInter0Clt95->subSequenceSize(), ofs,UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT); // Yes, always ".m"
2930  m_subMeanInter0Empirical90->subWriteContents (0,m_subMeanInter0Empirical90->subSequenceSize(), ofs,UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT); // Yes, always ".m"
2931  m_subMeanInter0Min->subWriteContents (0,m_subMeanInter0Min->subSequenceSize(), ofs,UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT); // Yes, always ".m"
2932  m_subMeanInter0Max->subWriteContents (0,m_subMeanInter0Max->subSequenceSize(), ofs,UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT); // Yes, always ".m"
2933 
2934  return;
2935 }
2936 // --------------------------------------------------
2937 template <class V, class M>
2938 void
2939 SequenceOfVectors<V,M>::unifiedMeanMonitorWrite(std::ofstream& ofs)
2940 {
2941  // std::set<unsigned int> tmpSet;
2942  // tmpSet.insert(0);
2943  m_unifiedMeanMonitorPosSeq->subWriteContents(0,m_unifiedMeanMonitorPosSeq->subSequenceSize(),ofs,UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT); // Yes, always ".m" // Yes, 'subWriteContents()'
2944  m_unifiedMeanVecSeq->subWriteContents (0,m_unifiedMeanVecSeq->subSequenceSize(), ofs,UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT); // Yes, always ".m" // Yes, 'subWriteContents()'
2945  m_unifiedMeanCltStdSeq->subWriteContents (0,m_unifiedMeanCltStdSeq->subSequenceSize(), ofs,UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT); // Yes, always ".m" // Yes, 'subWriteContents()'
2946 
2947  return;
2948 }
2949 // --------------------------------------------------
2950 template <class V, class M>
2951 void
2952 SequenceOfVectors<V,M>::subMeanMonitorFree()
2953 {
2954  delete m_subMeanMonitorPosSeq;
2955  m_subMeanMonitorPosSeq = NULL;
2956  delete m_subMeanVecSeq;
2957  m_subMeanVecSeq = NULL;
2958  delete m_subMeanCltStdSeq;
2959  m_subMeanCltStdSeq = NULL;
2960 
2961  return;
2962 }
2963 // --------------------------------------------------
2964 template <class V, class M>
2965 void
2966 SequenceOfVectors<V,M>::subMeanInter0MonitorFree()
2967 {
2968  delete m_subMeanInter0MonitorPosSeq;
2969  m_subMeanInter0MonitorPosSeq = NULL;
2970  delete m_subMeanInter0Mean;
2971  m_subMeanInter0Mean = NULL;
2972  delete m_subMeanInter0Clt95;
2973  m_subMeanInter0Clt95 = NULL;
2974  delete m_subMeanInter0Empirical90;
2975  m_subMeanInter0Empirical90 = NULL;
2976  delete m_subMeanInter0Min;
2977  m_subMeanInter0Min = NULL;
2978  delete m_subMeanInter0Max;
2979  m_subMeanInter0Max = NULL;
2980 
2981  return;
2982 }
2983 // --------------------------------------------------
2984 template <class V, class M>
2985 void
2986 SequenceOfVectors<V,M>::unifiedMeanMonitorFree()
2987 {
2988  delete m_unifiedMeanMonitorPosSeq;
2989  m_unifiedMeanMonitorPosSeq = NULL;
2990  delete m_unifiedMeanVecSeq;
2991  m_unifiedMeanVecSeq = NULL;
2992  delete m_unifiedMeanCltStdSeq;
2993  m_unifiedMeanCltStdSeq = NULL;
2994 
2995  return;
2996 }
2997 #endif // #ifdef UQ_CODE_HAS_MONITORS
2998 
2999 } // End namespace QUESO
3000 
3001 template class QUESO::SequenceOfVectors<QUESO::GslVector, QUESO::GslMatrix>;
QUESO::SequenceOfVectors::meanStacc
void meanStacc(unsigned int initialPos, V &meanStaccVec) const
Definition: SequenceOfVectors.C:2558
QUESO::ScalarSequence::unifiedScaleForKde
T unifiedScaleForKde(bool useOnlyInter0Comm, unsigned int initialPos, const T &unifiedIqrValue, unsigned int kdeDimension) const
Selects the scales (bandwidth) for the kernel density estimation, considering the unified sequence...
Definition: ScalarSequence.C:2207
QUESO::ArrayOfOneDGrids
Class to accommodate arrays of one-dimensional grid.
Definition: ArrayOfOneDGrids.h:47
QUESO::SequenceOfVectors::unifiedMeanCltStd
void unifiedMeanCltStd(unsigned int initialPos, unsigned int numPos, const V &unifiedMeanVec, V &unifiedSamVec) const
Definition: SequenceOfVectors.C:2379
QUESO::BaseVectorSequence::copy
void copy(const BaseVectorSequence< V, M > &src)
Copies vector sequence src to this.
Definition: VectorSequence.C:652
QUESO::ScalarSequence::subInterQuantileRange
T subInterQuantileRange(unsigned int initialPos) const
Returns the interquartile range of the values in the sub-sequence.
Definition: ScalarSequence.C:1963
QUESO::ScalarSequence::resizeSequence
void resizeSequence(unsigned int newSequenceSize)
Resizes the size of the sequence of scalars.
Definition: ScalarSequence.C:173
QUESO::SequenceOfVectors::writeSubMatlabHeader
void writeSubMatlabHeader(std::ofstream &ofs, double sequenceSize, double vectorSizeLocal) const
Helper function to write matlab-specific header info for vectors.
Definition: SequenceOfVectors.C:1459
QUESO::SequenceOfVectors::unifiedPopulationVariance
void unifiedPopulationVariance(unsigned int initialPos, unsigned int numPos, const V &unifiedMeanVec, V &unifiedPopVec) const
Finds the population variance of the unified sequence, considering numPos positions starting at posit...
Definition: SequenceOfVectors.C:698
QUESO::ScalarSequence::unifiedMeanExtra
T unifiedMeanExtra(bool useOnlyInter0Comm, unsigned int initialPos, unsigned int localNumPos) const
Finds the mean value of the unified sequence of numPos positions starting at position initialPos...
Definition: ScalarSequence.C:844
QUESO::SequenceOfVectors::subPopulationVariance
void subPopulationVariance(unsigned int initialPos, unsigned int numPos, const V &meanVec, V &popVec) const
Finds the population variance of the sub-sequence, considering numPos positions starting at position ...
Definition: SequenceOfVectors.C:666
QUESO::SequenceOfVectors::unifiedWriteContents
void unifiedWriteContents(const std::string &fileName, const std::string &fileType) const
Definition: SequenceOfVectors.C:1493
QUESO::VectorSpace
A class representing a vector space.
Definition: VectorSet.h:49
QUESO::SequenceOfVectors::subMeanCltStd
void subMeanCltStd(unsigned int initialPos, unsigned int numPos, const V &meanVec, V &samVec) const
Definition: SequenceOfVectors.C:2347
QUESO::SequenceOfVectors::unifiedUniformlySampledCdf
void unifiedUniformlySampledCdf(const V &numEvaluationPointsVec, ArrayOfOneDGrids< V, M > &unifiedCdfGrids, ArrayOfOneDTables< V, M > &unifiedCdfValues) const
Uniformly samples from the CDF from the sub-sequence.
Definition: SequenceOfVectors.C:262
QUESO::ScalarSequence::subMeanExtra
T subMeanExtra(unsigned int initialPos, unsigned int numPos) const
Finds the mean value of the sub-sequence, considering numPos positions starting at position initialPo...
Definition: ScalarSequence.C:806
QUESO::ScalarSequence::unifiedSampleVarianceExtra
T unifiedSampleVarianceExtra(bool useOnlyInter0Comm, unsigned int initialPos, unsigned int localNumPos, const T &unifiedMeanValue) const
Finds the sample variance of the unified sequence, considering numPos positions starting at position ...
Definition: ScalarSequence.C:1039
QUESO::SequenceOfVectors::geweke
void geweke(unsigned int initialPos, double ratioNa, double ratioNb, V &gewVec) const
Definition: SequenceOfVectors.C:2528
QUESO::ScalarSequence::unifiedGaussian1dKde
void unifiedGaussian1dKde(bool useOnlyInter0Comm, unsigned int initialPos, double unifiedScaleValue, const std::vector< T > &unifiedEvaluationPositions, std::vector< double > &unifiedDensityValues) const
Gaussian kernel for the KDE estimate of the unified sequence.
Definition: ScalarSequence.C:2308
QUESO::SequenceOfVectors::subGaussian1dKde
void subGaussian1dKde(unsigned int initialPos, const V &scaleVec, const std::vector< V * > &evalParamVecs, std::vector< V * > &densityVecs) const
Gaussian kernel for the KDE estimate of the sub-sequence.
Definition: SequenceOfVectors.C:1165
QUESO::SequenceOfVectors::resetValues
void resetValues(unsigned int initialPos, unsigned int numPos)
Resets a total of numPos values of the sequence starting at position initialPos.
Definition: SequenceOfVectors.C:126
QUESO::SequenceOfVectors::copy
void copy(const SequenceOfVectors< V, M > &src)
Copies vector sequence src to this.
Definition: SequenceOfVectors.C:2234
QUESO::ArrayOfOneDTables
Class to accommodate arrays of one-dimensional tables.
Definition: ArrayOfOneDTables.h:49
QUESO::UQ_OK_RC
const int UQ_OK_RC
Definition: Defines.h:92
QUESO::SequenceOfVectors::subMinMaxExtra
void subMinMaxExtra(unsigned int initialPos, unsigned int numPos, V &minVec, V &maxVec) const
Finds the minimum and the maximum values of the sub-sequence, considering numPos positions starting a...
Definition: SequenceOfVectors.C:892
QUESO::ScalarSequence::unifiedCdfPercentageRange
void unifiedCdfPercentageRange(bool useOnlyInter0Comm, unsigned int initialPos, unsigned int numPos, double range, T &lowerValue, T &upperValue) const
Definition: ScalarSequence.C:3764
QUESO::SequenceOfVectors::unifiedScalesForKde
void unifiedScalesForKde(unsigned int initialPos, const V &unifiedIqrVec, unsigned int kdeDimension, V &unifiedScaleVec) const
Selects the scales (bandwidth) for the kernel density estimation, considering the unified sequence...
Definition: SequenceOfVectors.C:1133
QUESO::ScalarSequence::unifiedSampleStd
T unifiedSampleStd(bool useOnlyInter0Comm, unsigned int initialPos, unsigned int localNumPos, const T &unifiedMeanValue) const
Finds the sample standard deviation of the unified sequence, considering localnumPos positions starti...
Definition: ScalarSequence.C:1126
QUESO::ScalarSequence::subScaleForKde
T subScaleForKde(unsigned int initialPos, const T &iqrValue, unsigned int kdeDimension) const
Selects the scales (output value) for the kernel density estimation, considering only the sub-sequenc...
Definition: ScalarSequence.C:2167
QUESO::ScalarSequence::unifiedPopulationVariance
T unifiedPopulationVariance(bool useOnlyInter0Comm, unsigned int initialPos, unsigned int numPos, const T &unifiedMeanValue) const
Finds the population variance of the unified sequence, considering numPos positions starting at posit...
Definition: ScalarSequence.C:1213
QUESO::ScalarSequence::unifiedMeanCltStd
T unifiedMeanCltStd(bool useOnlyInter0Comm, unsigned int initialPos, unsigned int localNumPos, const T &unifiedMeanValue) const
Definition: ScalarSequence.C:3434
QUESO::SequenceOfVectors::subInterQuantileRange
void subInterQuantileRange(unsigned int initialPos, V &iqrVec) const
Returns the interquartile range of the values in the sub-sequence.
Definition: SequenceOfVectors.C:1049
QUESO::SequenceOfVectors::unifiedReadContents
void unifiedReadContents(const std::string &fileName, const std::string &fileType, const unsigned int subSequenceSize)
Reads the unified sequence from a file.
Definition: SequenceOfVectors.C:1763
QUESO::ArrayOfOneDTables::setOneDTable
void setOneDTable(unsigned int rowId, const std::vector< double > &values)
Sets the one-dimensional table.
Definition: ArrayOfOneDTables.C:59
QUESO::ScalarSequence::unifiedMedianExtra
T unifiedMedianExtra(bool useOnlyInter0Comm, unsigned int initialPos, unsigned int localNumPos) const
Finds the median value of the unified sequence, considering numPos positions starting at position ini...
Definition: ScalarSequence.C:965
QUESO::SequenceOfVectors::subUniformlySampledMdf
void subUniformlySampledMdf(const V &numEvaluationPointsVec, ArrayOfOneDGrids< V, M > &mdfGrids, ArrayOfOneDTables< V, M > &mdfValues) const
Definition: SequenceOfVectors.C:2306
QUESO::ScalarSequence::subHistogram
void subHistogram(unsigned int initialPos, const T &minHorizontalValue, const T &maxHorizontalValue, std::vector< T > &centers, std::vector< unsigned int > &bins) const
Calculates the histogram of the sub-sequence.
Definition: ScalarSequence.C:1569
QUESO::SequenceOfVectors::unifiedCdfPercentageRange
void unifiedCdfPercentageRange(unsigned int initialPos, unsigned int numPos, double range, V &lowerVec, V &upperVec) const
Definition: SequenceOfVectors.C:2618
QUESO::ScalarSequence::autoCorrViaFft
void autoCorrViaFft(unsigned int initialPos, unsigned int numPos, unsigned int maxLag, std::vector< T > &autoCorrs) const
Calculates the autocorrelation via Fast Fourier transforms (FFT).
Definition: ScalarSequence.C:1331
QUESO::SequenceOfVectors::SequenceOfVectors
SequenceOfVectors(const VectorSpace< V, M > &vectorSpace, unsigned int subSequenceSize, const std::string &name)
Default constructor.
Definition: SequenceOfVectors.C:33
QUESO::SequenceOfVectors::unifiedSampleVarianceExtra
void unifiedSampleVarianceExtra(unsigned int initialPos, unsigned int numPos, const V &unifiedMeanVec, V &unifiedSamVec) const
Finds the sample variance of the unified sequence, considering numPos positions starting at position ...
Definition: SequenceOfVectors.C:568
QUESO::SequenceOfVectors::subMeanExtra
void subMeanExtra(unsigned int initialPos, unsigned int numPos, V &meanVec) const
Finds the mean value of the sub-sequence, considering numPos positions starting at position initialPo...
Definition: SequenceOfVectors.C:308
QUESO::SequenceOfVectors::subSequenceSize
unsigned int subSequenceSize() const
Size of the sub-sequence of vectors.
Definition: SequenceOfVectors.C:103
QUESO::ScalarSequence::psd
void psd(unsigned int initialPos, unsigned int numBlocks, double hopSizeRatio, std::vector< double > &psdSequence) const
Definition: ScalarSequence.C:3530
QUESO::SequenceOfVectors::subMeanInter0MonitorAlloc
void subMeanInter0MonitorAlloc(unsigned int numberOfMonitorPositions)
Definition: SequenceOfVectors.C:2768
QUESO::SequenceOfVectors::subMeanInter0MonitorStore
void subMeanInter0MonitorStore(unsigned int i, unsigned int monitorPosition, const V &subMeanInter0Mean, const V &subMeanInter0Clt95, const V &subMeanInter0Empirical90, const V &subMeanInter0Min, const V &subMeanInter0Max)
Definition: SequenceOfVectors.C:2880
QUESO::BaseVectorSequence::deleteStoredVectors
void deleteStoredVectors()
Deletes all the stored vectors.
Definition: VectorSequence.C:298
QUESO::SequenceOfVectors::unifiedInterQuantileRange
void unifiedInterQuantileRange(unsigned int initialPos, V &unifiedIqrVec) const
Returns the interquartile range of the values in the unified sequence.
Definition: SequenceOfVectors.C:1075
QUESO::ScalarSequence::subCdfStacc
void subCdfStacc(unsigned int initialPos, std::vector< double > &cdfStaccValues, std::vector< double > &cdfStaccValuesup, std::vector< double > &cdfStaccValueslow, const ScalarSequence< T > &sortedDataValues) const
Definition: ScalarSequence.C:3808
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::SequenceOfVectors::subMeanMonitorStore
void subMeanMonitorStore(unsigned int i, unsigned int monitorPosition, const V &subMeanVec, const V &subMeanCltStd)
Definition: SequenceOfVectors.C:2866
QUESO::SequenceOfVectors::unifiedMeanMonitorAlloc
void unifiedMeanMonitorAlloc(unsigned int numberOfMonitorPositions)
Definition: SequenceOfVectors.C:2782
QUESO::ScalarSequence::subPopulationVariance
T subPopulationVariance(unsigned int initialPos, unsigned int numPos, const T &meanValue) const
Finds the population variance of the sub-sequence, considering numPos positions starting at position ...
Definition: ScalarSequence.C:1186
QUESO::SequenceOfVectors::operator[]
const V * operator[](unsigned int posId) const
Definition: SequenceOfVectors.C:2282
QUESO::SequenceOfVectors::subCdfStacc
void subCdfStacc(unsigned int initialPos, std::vector< V * > &cdfStaccVecs, std::vector< V * > &cdfStaccVecsUp, std::vector< V * > &cdfStaccVecsLow, std::vector< V * > &sortedDataVecs) const
Definition: SequenceOfVectors.C:2653
QUESO::SequenceOfVectors::subScalesForKde
void subScalesForKde(unsigned int initialPos, const V &iqrVec, unsigned int kdeDimension, V &scaleVec) const
Selects the scales (bandwidth, scaleVec) for the kernel density estimation, considering only the sub-...
Definition: SequenceOfVectors.C:1102
QUESO::SequenceOfVectors::unifiedMeanMonitorWrite
void unifiedMeanMonitorWrite(std::ofstream &ofs)
Definition: SequenceOfVectors.C:2939
QUESO::SequenceOfVectors::subMeanMonitorAlloc
void subMeanMonitorAlloc(unsigned int numberOfMonitorPositions)
Definition: SequenceOfVectors.C:2757
QUESO::SequenceOfVectors::writeUnifiedMatlabHeader
void writeUnifiedMatlabHeader(std::ofstream &ofs, double sequenceSize, double vectorSizeLocal) const
Definition: SequenceOfVectors.C:1471
QUESO::ScalarSequence::unifiedHistogram
void unifiedHistogram(bool useOnlyInter0Comm, unsigned int initialPos, const T &unifiedMinHorizontalValue, const T &unifiedMaxHorizontalValue, std::vector< T > &unifiedCenters, std::vector< unsigned int > &unifiedBins) const
Calculates the histogram of the unified sequence.
Definition: ScalarSequence.C:1616
QUESO::SequenceOfVectors
Class for handling vector samples (sequence of vectors).
Definition: SequenceOfVectors.h:49
QUESO::ScalarSequence::subMinMaxExtra
void subMinMaxExtra(unsigned int initialPos, unsigned int numPos, T &minValue, T &maxValue) const
Finds the minimum and the maximum values of the sub-sequence, considering numPos positions starting a...
Definition: ScalarSequence.C:1469
QUESO::SequenceOfVectors::extractScalarSeq
void extractScalarSeq(unsigned int initialPos, unsigned int spacing, unsigned int numPos, unsigned int paramId, ScalarSequence< double > &scalarSeq) const
Extracts a sequence of scalars.
Definition: SequenceOfVectors.C:2210
QUESO::ArrayOfOneDGrids::setUniformGrids
void setUniformGrids(const V &sizesVec, const V &minPositionsVec, const V &maxPositionsVec)
Sets an array of uniform grids.
Definition: ArrayOfOneDGrids.C:102
QUESO::SequenceOfVectors::unifiedGaussian1dKde
void unifiedGaussian1dKde(unsigned int initialPos, const V &unifiedScaleVec, const std::vector< V * > &unifiedEvalParamVecs, std::vector< V * > &unifiedDensityVecs) const
Gaussian kernel for the KDE estimate of the unified sequence.
Definition: SequenceOfVectors.C:1214
QUESO::BaseVectorSequence
Base class for handling vector and array samples (sequence of vectors or arrays). ...
Definition: VectorSequence.h:56
QUESO::SequenceOfVectors::~SequenceOfVectors
~SequenceOfVectors()
Destructor.
QUESO::SequenceOfVectors::unifiedMinMaxExtra
void unifiedMinMaxExtra(unsigned int initialPos, unsigned int numPos, V &unifiedMinVec, V &unifiedMaxVec) const
Finds the minimum and the maximum values of the unified sequence, considering numPos positions starti...
Definition: SequenceOfVectors.C:922
QUESO::SequenceOfVectors::m_seq
std::vector< const V * > m_seq
Sequence of vectors.
Definition: SequenceOfVectors.h:435
QUESO::SequenceOfVectors::filter
void filter(unsigned int initialPos, unsigned int spacing)
Filters positions in the sequence of vectors.
Definition: SequenceOfVectors.C:2042
QUESO::SequenceOfVectors::unifiedMeanMonitorStore
void unifiedMeanMonitorStore(unsigned int i, unsigned int monitorPosition, V &unifiedMeanVec, V &unifiedMeanCltStd)
Definition: SequenceOfVectors.C:2900
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::ScalarSequence::autoCovariance
T autoCovariance(unsigned int initialPos, unsigned int numPos, const T &meanValue, unsigned int lag) const
Calculates the autocovariance.
Definition: ScalarSequence.C:1272
QUESO::SequenceOfVectors::writeTxtHeader
void writeTxtHeader(std::ofstream &ofs, double sequenceSize, double vectorSizeLocal) const
Helper function to write plain txt info for vectors.
Definition: SequenceOfVectors.C:1483
QUESO::SequenceOfVectors::seqVectorPositionIteratorTypedef
std::vector< const V * >::iterator seqVectorPositionIteratorTypedef
Definition: SequenceOfVectors.h:56
QUESO::ScalarSequence::unifiedMinMaxExtra
void unifiedMinMaxExtra(bool useOnlyInter0Comm, unsigned int initialPos, unsigned int numPos, T &unifiedMinValue, T &unifiedMaxValue) const
Finds the minimum and the maximum values of the unified sequence, considering numPos positions starti...
Definition: ScalarSequence.C:1498
QUESO::SequenceOfVectors::subCdfPercentageRange
void subCdfPercentageRange(unsigned int initialPos, unsigned int numPos, double range, V &lowerVec, V &upperVec) const
Definition: SequenceOfVectors.C:2584
QUESO::FilePtrSetStruct::ifsVar
std::ifstream * ifsVar
Provides a stream interface to read data from files.
Definition: Environment.h:89
QUESO::SequenceOfVectors::autoCorrViaFft
void autoCorrViaFft(unsigned int initialPos, unsigned int numPos, const std::vector< unsigned int > &lags, std::vector< V * > &corrVecs) const
Calculates the autocorrelation via Fast Fourier transforms (FFT).
Definition: SequenceOfVectors.C:798
QUESO::SequenceOfVectors::subMeanMonitorRun
void subMeanMonitorRun(unsigned int monitorPosition, V &subMeanVec, V &subMeanCltStd)
Definition: SequenceOfVectors.C:2793
QUESO::ScalarSequence::unifiedInterQuantileRange
T unifiedInterQuantileRange(bool useOnlyInter0Comm, unsigned int initialPos) const
Returns the interquartile range of the values in the unified sequence.
Definition: ScalarSequence.C:2076
QUESO::SequenceOfVectors::extractRawData
void extractRawData(unsigned int initialPos, unsigned int spacing, unsigned int numPos, unsigned int paramId, std::vector< double > &rawData) const
Extracts the raw data.
Definition: SequenceOfVectors.C:2253
QUESO::SequenceOfVectors::select
void select(const std::vector< unsigned int > &idsOfUniquePositions)
TODO: It shall select positions in the sequence of vectors.
Definition: SequenceOfVectors.C:2033
QUESO::ScalarSequence::bmm
T bmm(unsigned int initialPos, unsigned int batchLength) const
Definition: ScalarSequence.C:3494
QUESO::SequenceOfVectors::subMeanInter0MonitorRun
void subMeanInter0MonitorRun(unsigned int currentPosition, V &subMeanInter0Mean, V &subMeanInter0Clt95, V &subMeanInter0Empirical90, V &subMeanInter0Min, V &subMeanInter0Max)
Definition: SequenceOfVectors.C:2811
QUESO::SequenceOfVectors::unifiedMedianExtra
void unifiedMedianExtra(unsigned int initialPos, unsigned int numPos, V &unifiedMedianVec) const
Finds the median value of the unfed sequence, considering numPos positions starting at position initi...
Definition: SequenceOfVectors.C:477
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::SequenceOfVectors::operator=
SequenceOfVectors< V, M > & operator=(const SequenceOfVectors< V, M > &rhs)
Copies values from rhs to this.
Definition: SequenceOfVectors.C:94
QUESO::matrixProduct
GslMatrix matrixProduct(const GslVector &v1, const GslVector &v2)
Definition: GslMatrix.C:2036
QUESO::ScalarSequence::subGaussian1dKde
void subGaussian1dKde(unsigned int initialPos, double scaleValue, const std::vector< T > &evaluationPositions, std::vector< double > &densityValues) const
Gaussian kernel for the KDE estimate of the sub-sequence.
Definition: ScalarSequence.C:2278
QUESO::queso_require_equal_to_msg
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"))
Definition: MonteCarloSGOptions.C:460
QUESO::SequenceOfVectors::subWriteContents
void subWriteContents(unsigned int initialPos, unsigned int numPos, const std::string &fileName, const std::string &fileType, const std::set< unsigned int > &allowedSubEnvIds) const
Writes the sub-sequence to a file.
Definition: SequenceOfVectors.C:1264
QUESO::ScalarSequence::subMeanCltStd
T subMeanCltStd(unsigned int initialPos, unsigned int numPos, const T &meanValue) const
Definition: ScalarSequence.C:3405
QUESO::ScalarSequence::subSort
void subSort(unsigned int initialPos, ScalarSequence< T > &sortedSequence) const
Sorts the sub-sequence of scalars.
Definition: ScalarSequence.C:1842
QUESO::ScalarSequence::subMedianExtra
T subMedianExtra(unsigned int initialPos, unsigned int numPos) const
Finds the median value of the sub-sequence, considering numPos positions starting at position initial...
Definition: ScalarSequence.C:922
QUESO::SequenceOfVectors::subUniformlySampledCdf
void subUniformlySampledCdf(const V &numEvaluationPointsVec, ArrayOfOneDGrids< V, M > &cdfGrids, ArrayOfOneDTables< V, M > &cdfValues) const
Uniformly samples from the CDF from the sub-sequence.
Definition: SequenceOfVectors.C:227
QUESO::SequenceOfVectors::subHistogram
void subHistogram(unsigned int initialPos, const V &minVec, const V &maxVec, std::vector< V * > &centersForAllBins, std::vector< V * > &quanttsForAllBins) const
Calculates the histogram of the sub-sequence.
Definition: SequenceOfVectors.C:956
QUESO::ScalarSequence::meanStacc
T meanStacc(unsigned int initialPos) const
Definition: ScalarSequence.C:3718
QUESO::SequenceOfVectors::unifiedMeanMonitorRun
void unifiedMeanMonitorRun(unsigned int currentPosition, V &unifiedMeanVec, V &unifiedMeanCltStd)
Definition: SequenceOfVectors.C:2849
QUESO::ScalarSequence::autoCorrViaDef
T autoCorrViaDef(unsigned int initialPos, unsigned int numPos, unsigned int lag) const
Calculates the autocorrelation via definition.
Definition: ScalarSequence.C:1302
QUESO::ScalarSequence::subSampleStd
T subSampleStd(unsigned int initialPos, unsigned int numPos, const T &meanValue) const
Finds the sample standard deviation of the unified sequence, considering numPos positions starting at...
Definition: ScalarSequence.C:1098
QUESO::SequenceOfVectors::subMeanInter0MonitorWrite
void subMeanInter0MonitorWrite(std::ofstream &ofs)
Definition: SequenceOfVectors.C:2925
QUESO::SequenceOfVectors::estimateConvBrooksGelman
double estimateConvBrooksGelman(unsigned int initialPos, unsigned int numPos) const
Estimates convergence rate using Brooks &amp; Gelman method.
Definition: SequenceOfVectors.C:2083
QUESO::ScalarSequence::subCdfPercentageRange
void subCdfPercentageRange(unsigned int initialPos, unsigned int numPos, double range, T &lowerValue, T &upperValue) const
Definition: ScalarSequence.C:3730
QUESO::SequenceOfVectors::resizeSequence
void resizeSequence(unsigned int newSubSequenceSize)
Resizes the sequence.
Definition: SequenceOfVectors.C:110
QUESO::SequenceOfVectors::unifiedHistogram
void unifiedHistogram(unsigned int initialPos, const V &unifiedMinVec, const V &unifiedMaxVec, std::vector< V * > &unifiedCentersForAllBins, std::vector< V * > &unifiedQuanttsForAllBins) const
Calculates the histogram of the unified sequence.
Definition: SequenceOfVectors.C:1002
QUESO::FilePtrSetStruct::ofsVar
std::ofstream * ofsVar
Provides a stream interface to write data to files.
Definition: Environment.h:86
QUESO::SequenceOfVectors::unifiedMeanExtra
void unifiedMeanExtra(unsigned int initialPos, unsigned int numPos, V &unifiedMeanVec) const
Finds the mean value of the unified sequence, considering numPos positions starting at position initi...
Definition: SequenceOfVectors.C:363
QUESO::SequenceOfVectors::subSampleVarianceExtra
void subSampleVarianceExtra(unsigned int initialPos, unsigned int numPos, const V &meanVec, V &samVec) const
Finds the sample variance of the sub-sequence, considering numPos positions starting at position init...
Definition: SequenceOfVectors.C:536
QUESO::SequenceOfVectors::subMedianExtra
void subMedianExtra(unsigned int initialPos, unsigned int numPos, V &medianVec) const
Finds the median value of the sub-sequence, considering numPos positions starting at position initial...
Definition: SequenceOfVectors.C:422
QUESO::FilePtrSetStruct
Struct for handling data input and output from files.
Definition: Environment.h:77
QUESO::SequenceOfVectors::bmm
void bmm(unsigned int initialPos, unsigned int batchLength, V &bmmVec) const
Definition: SequenceOfVectors.C:2412
QUESO::SequenceOfVectors::subMeanMonitorWrite
void subMeanMonitorWrite(std::ofstream &ofs)
Definition: SequenceOfVectors.C:2914
QUESO::ScalarSequence::geweke
T geweke(unsigned int initialPos, double ratioNa, double ratioNb) const
Definition: ScalarSequence.C:3653
QUESO::MiscGetEllapsedSeconds
double MiscGetEllapsedSeconds(struct timeval *timeval0)
Definition: Miscellaneous.C:282
QUESO::SequenceOfVectors::psdAtZero
void psdAtZero(unsigned int initialPos, unsigned int numBlocks, double hopSizeRatio, V &psdVec) const
Definition: SequenceOfVectors.C:2495
QUESO::SequenceOfVectors::autoCovariance
void autoCovariance(unsigned int initialPos, unsigned int numPos, const V &meanVec, unsigned int lag, V &covVec) const
Calculates the autocovariance.
Definition: SequenceOfVectors.C:731
QUESO::SequenceOfVectors::psd
void psd(unsigned int initialPos, unsigned int numBlocks, double hopSizeRatio, unsigned int paramId, std::vector< double > &psdResult) const
Definition: SequenceOfVectors.C:2467
QUESO::SequenceOfVectors::autoCorrViaDef
void autoCorrViaDef(unsigned int initialPos, unsigned int numPos, unsigned int lag, V &corrVec) const
Calculates the autocorrelation via definition.
Definition: SequenceOfVectors.C:766
QUESO::SequenceOfVectors::fftForward
void fftForward(unsigned int initialPos, unsigned int fftSize, unsigned int paramId, std::vector< std::complex< double > > &fftResult) const
Definition: SequenceOfVectors.C:2440
QUESO::ScalarSequence::subSampleVarianceExtra
T subSampleVarianceExtra(unsigned int initialPos, unsigned int numPos, const T &meanValue) const
Finds the sample variance of the sub-sequence, considering numPos positions starting at position init...
Definition: ScalarSequence.C:1012
QUESO::SequenceOfVectors::erasePositions
void erasePositions(unsigned int initialPos, unsigned int numPos)
Erases numPos elements of the sequence starting at position initialPos.
Definition: SequenceOfVectors.C:154
QUESO::SequenceOfVectors::unifiedSampleStd
void unifiedSampleStd(unsigned int initialPos, unsigned int numPos, const V &unifiedMeanVec, V &unifiedStdVec) const
Finds the sample standard deviation of the unified sequence, considering numPos positions starting at...
Definition: SequenceOfVectors.C:633
Generated on Tue Jun 5 2018 19:48:54 for queso-0.57.1 by  doxygen 1.8.5