25 #include <queso/SequenceOfVectors.h>
26 #include <queso/GslVector.h>
27 #include <queso/GslMatrix.h>
32 template <
class V,
class M>
35 unsigned int subSequenceSize,
36 const std::string& name)
39 m_seq (subSequenceSize,NULL)
40 #ifdef UQ_CODE_HAS_MONITORS
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)
67 template <
class V,
class M>
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;
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;
82 if (m_unifiedMeanMonitorPosSeq)
delete m_unifiedMeanMonitorPosSeq;
83 if (m_unifiedMeanVecSeq )
delete m_unifiedMeanVecSeq;
84 if (m_unifiedMeanCltStdSeq )
delete m_unifiedMeanCltStdSeq;
87 for (
unsigned int i = 0; i < (
unsigned int) m_seq.size(); ++i) {
88 if (m_seq[i])
delete m_seq[i];
92 template <
class V,
class M>
101 template <
class V,
class M>
108 template <
class V,
class M>
112 if (newSubSequenceSize != this->subSequenceSize()) {
113 if (newSubSequenceSize < this->subSequenceSize()) {
114 this->resetValues(newSubSequenceSize,this->subSequenceSize()-newSubSequenceSize);
116 m_seq.resize(newSubSequenceSize,NULL);
117 std::vector<const V*>(m_seq).swap(m_seq);
124 template <
class V,
class M>
128 bool bRC = ((initialPos < this->subSequenceSize()) &&
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
140 "SequenceOfVectors<V,M>::resetValues()",
141 "invalid input data");
143 for (
unsigned int j = 0; j < numPos; ++j) {
144 if (m_seq[initialPos+j] != NULL) {
145 delete m_seq[initialPos+j];
146 m_seq[initialPos+j] = NULL;
155 template <
class V,
class M>
159 bool bRC = ((initialPos < this->subSequenceSize()) &&
161 ((initialPos+numPos) <= this->subSequenceSize()));
164 "SequenceOfVectors<V,M>::erasePositions()",
165 "invalid input data");
167 for (
unsigned int j = 0; j < numPos; ++j) {
168 if (m_seq[initialPos+j] != NULL) {
169 delete m_seq[initialPos+j];
170 m_seq[initialPos+j] = NULL;
175 if (initialPos < this->subSequenceSize()) std::advance(posIteratorBegin,initialPos);
176 else posIteratorBegin = m_seq.end();
178 unsigned int posEnd = initialPos + numPos - 1;
180 if (posEnd < this->subSequenceSize()) std::advance(posIteratorEnd,posEnd);
181 else posIteratorEnd = m_seq.end();
183 unsigned int oldSubSequenceSize = this->subSequenceSize();
184 m_seq.erase(posIteratorBegin,posIteratorEnd);
187 "SequenceOfVectors::erasePositions()",
188 "(oldSubSequenceSize - numPos) != this->subSequenceSize()");
195 template <
class V,
class M>
201 "SequenceOfVectorss<V,M>::getPositionValues()",
202 "posId > subSequenceSize()");
206 "SequenceOfVectorss<V,M>::getPositionValues()",
214 vec = *(m_seq[posId]);
219 template <
class V,
class M>
225 "SequenceOfVectorss<V,M>::setPositionValues()",
226 "posId > subSequenceSize()");
230 "SequenceOfVectorss<V,M>::setPositionValues()",
233 if (m_seq[posId] != NULL)
delete m_seq[posId];
234 m_seq[posId] =
new V(vec);
246 template <
class V,
class M>
249 const V& numEvaluationPointsVec,
253 V minDomainValues(m_vectorSpace.zeroVector());
254 V maxDomainValues(m_vectorSpace.zeroVector());
258 unsigned int numParams = this->vectorSizeLocal();
259 for (
unsigned int i = 0; i < numParams; ++i) {
260 this->extractScalarSeq(0,
266 std::vector<double> aCdf(0);
267 data.subUniformlySampledCdf((
unsigned int) numEvaluationPointsVec[i],
281 template <
class V,
class M>
284 const V& numEvaluationPointsVec,
288 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 10)) {
289 *m_env.subDisplayFile() <<
"Entering SequenceOfVectors<V,M>::unifiedUniformlySampledCdf()"
293 V unifiedMinDomainValues(m_vectorSpace.zeroVector());
294 V unifiedMaxDomainValues(m_vectorSpace.zeroVector());
298 unsigned int numParams = this->vectorSizeLocal();
299 for (
unsigned int i = 0; i < numParams; ++i) {
300 this->extractScalarSeq(0,
306 std::vector<double> aCdf(0);
307 data.unifiedUniformlySampledCdf(m_vectorSpace.numOfProcsForStorage() == 1,
308 (
unsigned int) numEvaluationPointsVec[i],
309 unifiedMinDomainValues[i],
310 unifiedMaxDomainValues[i],
316 unifiedMinDomainValues,
317 unifiedMaxDomainValues);
319 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 10)) {
320 *m_env.subDisplayFile() <<
"Leaving SequenceOfVectors<V,M>::unifiedUniformlySampledCdf()"
327 template <
class V,
class M>
330 unsigned int initialPos,
334 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 5)) {
335 *m_env.subDisplayFile() <<
"Entering SequenceOfVectors<V,M>::subMeanExtra()"
336 <<
": initialPos = " << initialPos
337 <<
", numPos = " << numPos
338 <<
", sub sequence size = " << this->subSequenceSize()
342 bool bRC = ((initialPos < this->subSequenceSize()) &&
344 ((initialPos+numPos) <= this->subSequenceSize()) &&
345 (this->vectorSizeLocal() == meanVec.sizeLocal() ));
346 if ((bRC ==
false) && (m_env.subDisplayFile())) {
347 *m_env.subDisplayFile() <<
"In SequenceOfVectors<V,M>::subMeanExtra()"
348 <<
", initialPos = " << initialPos
349 <<
", this->subSequenceSize() = " << this->subSequenceSize()
350 <<
", numPos = " << numPos
351 <<
", this->vectorSizeLocal() = " << this->vectorSizeLocal()
352 <<
", meanVec.sizeLocal() = " << meanVec.sizeLocal()
357 "SequenceOfVectors<V,M>::subMeanExtra()",
358 "invalid input data");
362 unsigned int numParams = this->vectorSizeLocal();
363 for (
unsigned int i = 0; i < numParams; ++i) {
364 this->extractScalarSeq(initialPos,
369 meanVec[i] = data.subMeanExtra(0,
373 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 5)) {
374 *m_env.subDisplayFile() <<
"Leaving SequenceOfVectors<V,M>::subMeanExtra()"
375 <<
": initialPos = " << initialPos
376 <<
", numPos = " << numPos
377 <<
", sub sequence size = " << this->subSequenceSize()
378 <<
", meanVec = " << meanVec
385 template <
class V,
class M>
388 unsigned int initialPos,
390 V& unifiedMeanVec)
const
392 unsigned int tmpUnif = this->unifiedSequenceSize();
393 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 5)) {
394 *m_env.subDisplayFile() <<
"Entering SequenceOfVectors<V,M>::unifiedMeanExtra()"
395 <<
": initialPos = " << initialPos
396 <<
", numPos = " << numPos
397 <<
", sub sequence size = " << this->subSequenceSize()
398 <<
", unified sequence size = " << tmpUnif
402 bool bRC = ((initialPos < this->subSequenceSize() ) &&
404 ((initialPos+numPos) <= this->subSequenceSize() ) &&
405 (this->vectorSizeLocal() == unifiedMeanVec.sizeLocal()));
406 if ((bRC ==
false) && (m_env.subDisplayFile())) {
407 *m_env.subDisplayFile() <<
"In SequenceOfVectors<V,M>::unifiedMeanExtra()"
408 <<
", initialPos = " << initialPos
409 <<
", this->subSequenceSize() = " << this->subSequenceSize()
410 <<
", numPos = " << numPos
411 <<
", this->vectorSizeLocal() = " << this->vectorSizeLocal()
412 <<
", unifiedMeanVec.sizeLocal() = " << unifiedMeanVec.sizeLocal()
417 "SequenceOfVectors<V,M>::unifiedMeanExtra()",
418 "invalid input data");
422 unsigned int numParams = this->vectorSizeLocal();
423 for (
unsigned int i = 0; i < numParams; ++i) {
424 this->extractScalarSeq(initialPos,
429 unifiedMeanVec[i] = data.unifiedMeanExtra(m_vectorSpace.numOfProcsForStorage() == 1,
434 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 5)) {
435 *m_env.subDisplayFile() <<
"Leaving SequenceOfVectors<V,M>::unifiedMeanExtra()"
436 <<
": initialPos = " << initialPos
437 <<
", numPos = " << numPos
438 <<
", sub sequence size = " << this->subSequenceSize()
439 <<
", unified sequence size = " << tmpUnif
440 <<
", unifiedMeanVec = " << unifiedMeanVec
447 template <
class V,
class M>
450 unsigned int initialPos,
454 if (this->subSequenceSize() == 0)
return;
456 bool bRC = ((initialPos < this->subSequenceSize()) &&
458 ((initialPos+numPos) <= this->subSequenceSize()));
460 std::cerr <<
"In SequenceOfVectors<V,M>::subMedianExtra()"
461 <<
": ERROR at fullRank " << m_env.fullRank()
462 <<
", initialPos = " << initialPos
463 <<
", numPos = " << numPos
464 <<
", this->subSequenceSize() = " << this->subSequenceSize()
466 if (m_env.subDisplayFile()) {
467 *m_env.subDisplayFile() <<
"In SequenceOfVectors<V,M>::subMedianExtra()"
468 <<
": ERROR at fullRank " << m_env.fullRank()
469 <<
", initialPos = " << initialPos
470 <<
", numPos = " << numPos
471 <<
", this->subSequenceSize() = " << this->subSequenceSize()
477 "SequenceOfVectors<V,M>::subMedianExtra()",
478 "invalid input data");
482 unsigned int numParams = this->vectorSizeLocal();
483 for (
unsigned int i = 0; i < numParams; ++i) {
484 this->extractScalarSeq(initialPos,
489 medianVec[i] = data.subMedianExtra(0,
493 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 5)) {
494 *m_env.subDisplayFile() <<
"Leaving SequenceOfVectors<V,M>::subMedianExtra()"
495 <<
": initialPos = " << initialPos
496 <<
", numPos = " << numPos
497 <<
", sub sequence size = " << this->subSequenceSize()
498 <<
", medianVec = " << medianVec
505 template <
class V,
class M>
508 unsigned int initialPos,
510 V& unifiedMedianVec)
const
512 unsigned int tmpUnif = this->unifiedSequenceSize();
513 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 5)) {
514 *m_env.subDisplayFile() <<
"Entering SequenceOfVectors<V,M>::unifiedMedianExtra()"
515 <<
": initialPos = " << initialPos
516 <<
", numPos = " << numPos
517 <<
", sub sequence size = " << this->subSequenceSize()
518 <<
", unified sequence size = " << tmpUnif
522 bool bRC = ((initialPos < this->subSequenceSize() ) &&
524 ((initialPos+numPos) <= this->subSequenceSize() ) &&
525 (this->vectorSizeLocal() == unifiedMedianVec.sizeLocal()));
526 if ((bRC ==
false) && (m_env.subDisplayFile())) {
527 *m_env.subDisplayFile() <<
"In SequenceOfVectors<V,M>::unifiedMedianExtra()"
528 <<
", initialPos = " << initialPos
529 <<
", this->subSequenceSize() = " << this->subSequenceSize()
530 <<
", numPos = " << numPos
531 <<
", this->vectorSizeLocal() = " << this->vectorSizeLocal()
532 <<
", unifiedMedianVec.sizeLocal() = " << unifiedMedianVec.sizeLocal()
537 "SequenceOfVectors<V,M>::unifiedMedianExtra()",
538 "invalid input data");
542 unsigned int numParams = this->vectorSizeLocal();
543 for (
unsigned int i = 0; i < numParams; ++i) {
544 this->extractScalarSeq(initialPos,
549 unifiedMedianVec[i] = data.unifiedMedianExtra(m_vectorSpace.numOfProcsForStorage() == 1,
554 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 5)) {
555 *m_env.subDisplayFile() <<
"Leaving SequenceOfVectors<V,M>::unifiedMedianExtra()"
556 <<
": initialPos = " << initialPos
557 <<
", numPos = " << numPos
558 <<
", sub sequence size = " << this->subSequenceSize()
559 <<
", unified sequence size = " << tmpUnif
560 <<
", unifiedMedianVec = " << unifiedMedianVec
567 template <
class V,
class M>
570 unsigned int initialPos,
575 bool bRC = ((initialPos < this->subSequenceSize()) &&
577 ((initialPos+numPos) <= this->subSequenceSize()) &&
578 (this->vectorSizeLocal() == meanVec.sizeLocal() ) &&
579 (this->vectorSizeLocal() == samVec.sizeLocal() ));
582 "SequenceOfVectors<V,M>::subSampleVarianceExtra()",
583 "invalid input data");
587 unsigned int numParams = this->vectorSizeLocal();
588 for (
unsigned int i = 0; i < numParams; ++i) {
589 this->extractScalarSeq(initialPos,
594 samVec[i] = data.subSampleVarianceExtra(0,
602 template <
class V,
class M>
605 unsigned int initialPos,
607 const V& unifiedMeanVec,
608 V& unifiedSamVec)
const
610 bool bRC = ((initialPos < this->subSequenceSize() ) &&
612 ((initialPos+numPos) <= this->subSequenceSize() ) &&
613 (this->vectorSizeLocal() == unifiedMeanVec.sizeLocal()) &&
614 (this->vectorSizeLocal() == unifiedSamVec.sizeLocal() ));
617 "SequenceOfVectors<V,M>::unifiedSampleVarianceExtra()",
618 "invalid input data");
622 unsigned int numParams = this->vectorSizeLocal();
623 for (
unsigned int i = 0; i < numParams; ++i) {
624 this->extractScalarSeq(initialPos,
629 unifiedSamVec[i] = data.unifiedSampleVarianceExtra(m_vectorSpace.numOfProcsForStorage() == 1,
638 template <
class V,
class M>
641 unsigned int initialPos,
646 bool bRC = ((initialPos < this->subSequenceSize()) &&
648 ((initialPos+numPos) <= this->subSequenceSize()) &&
649 (this->vectorSizeLocal() == meanVec.sizeLocal() ) &&
650 (this->vectorSizeLocal() == stdvec.sizeLocal() ));
653 "SequenceOfVectors<V,M>::subSampleStd()",
654 "invalid input data");
658 unsigned int numParams = this->vectorSizeLocal();
659 for (
unsigned int i = 0; i < numParams; ++i) {
660 this->extractScalarSeq(initialPos,
665 stdvec[i] = data.subSampleStd(0,
673 template <
class V,
class M>
676 unsigned int initialPos,
678 const V& unifiedMeanVec,
679 V& unifiedStdVec)
const
681 bool bRC = ((initialPos < this->subSequenceSize() ) &&
683 ((initialPos+numPos) <= this->subSequenceSize() ) &&
684 (this->vectorSizeLocal() == unifiedMeanVec.sizeLocal()) &&
685 (this->vectorSizeLocal() == unifiedStdVec.sizeLocal() ));
688 "SequenceOfVectors<V,M>::unifiedSampleStd()",
689 "invalid input data");
693 unsigned int numParams = this->vectorSizeLocal();
694 for (
unsigned int i = 0; i < numParams; ++i) {
695 this->extractScalarSeq(initialPos,
700 unifiedStdVec[i] = data.unifiedSampleStd(m_vectorSpace.numOfProcsForStorage() == 1,
709 template <
class V,
class M>
712 unsigned int initialPos,
717 bool bRC = ((initialPos < this->subSequenceSize()) &&
719 ((initialPos+numPos) <= this->subSequenceSize()) &&
720 (this->vectorSizeLocal() == meanVec.sizeLocal() ) &&
721 (this->vectorSizeLocal() == popVec.sizeLocal() ));
724 "SequenceOfVectors<V,M>::subPopulationVariance()",
725 "invalid input data");
729 unsigned int numParams = this->vectorSizeLocal();
730 for (
unsigned int i = 0; i < numParams; ++i) {
731 this->extractScalarSeq(initialPos,
736 popVec[i] = data.subPopulationVariance(0,
744 template <
class V,
class M>
747 unsigned int initialPos,
749 const V& unifiedMeanVec,
750 V& unifiedPopVec)
const
752 bool bRC = ((initialPos < this->subSequenceSize() ) &&
754 ((initialPos+numPos) <= this->subSequenceSize() ) &&
755 (this->vectorSizeLocal() == unifiedMeanVec.sizeLocal()) &&
756 (this->vectorSizeLocal() == unifiedPopVec.sizeLocal() ));
759 "SequenceOfVectors<V,M>::unifiedPopulationVariance()",
760 "invalid input data");
764 unsigned int numParams = this->vectorSizeLocal();
765 for (
unsigned int i = 0; i < numParams; ++i) {
766 this->extractScalarSeq(initialPos,
771 unifiedPopVec[i] = data.unifiedPopulationVariance(m_vectorSpace.numOfProcsForStorage() == 1,
780 template <
class V,
class M>
783 unsigned int initialPos,
789 bool bRC = ((initialPos < this->subSequenceSize()) &&
791 ((initialPos+numPos) <= this->subSequenceSize()) &&
792 (this->vectorSizeLocal() == meanVec.sizeLocal() ) &&
794 (this->vectorSizeLocal() == covVec.sizeLocal() ));
797 "SequenceOfVectors<V,M>::autoCovariance()",
798 "invalid input data");
802 unsigned int numParams = this->vectorSizeLocal();
803 for (
unsigned int i = 0; i < numParams; ++i) {
804 this->extractScalarSeq(initialPos,
809 covVec[i] = data.autoCovariance(0,
818 template <
class V,
class M>
821 unsigned int initialPos,
826 bool bRC = ((initialPos < this->subSequenceSize()) &&
828 ((initialPos+numPos) <= this->subSequenceSize()) &&
830 (this->vectorSizeLocal() == corrVec.sizeLocal() ));
833 "SequenceOfVectors<V,M>::autoCorrViaDef()",
834 "invalid input data");
838 unsigned int numParams = this->vectorSizeLocal();
839 for (
unsigned int i = 0; i < numParams; ++i) {
840 this->extractScalarSeq(initialPos,
845 corrVec[i] = data.autoCorrViaDef(0,
853 template <
class V,
class M>
856 unsigned int initialPos,
858 const std::vector<unsigned int>& lags,
859 std::vector<V*>& corrVecs)
const
861 bool bRC = ((initialPos < this->subSequenceSize()) &&
863 ((initialPos+numPos) <= this->subSequenceSize()) &&
864 (0 < lags.size() ) &&
865 (lags[lags.size()-1] < numPos ));
868 "SequenceOfVectors<V,M>::autoCorrViaFft()",
869 "invalid input data");
871 for (
unsigned int j = lags.size(); j < corrVecs.size(); ++j) {
872 if (corrVecs[j] != NULL) {
877 corrVecs.resize(lags.size(),NULL);
878 for (
unsigned int j = 0; j < corrVecs.size(); ++j) {
879 if (corrVecs[j] == NULL) corrVecs[j] =
new V(m_vectorSpace.zeroVector());
883 unsigned int maxLag = lags[lags.size()-1];
884 std::vector<double> autoCorrs(maxLag+1,0.);
886 unsigned int numParams = this->vectorSizeLocal();
887 for (
unsigned int i = 0; i < numParams; ++i) {
888 this->extractScalarSeq(initialPos,
907 for (
unsigned int j = 0; j < lags.size(); ++j) {
908 (*(corrVecs[j]))[i] = autoCorrs[lags[j]];
915 template <
class V,
class M>
918 unsigned int initialPos,
921 V& autoCorrsSumVec)
const
923 bool bRC = ((initialPos < this->subSequenceSize()) &&
925 ((initialPos+numPos) <= this->subSequenceSize()) &&
927 (numSum <= numPos ) &&
928 (autoCorrsSumVec.sizeLocal() == this->vectorSizeLocal()));
931 "SequenceOfVectors<V,M>::autoCorrViaFft(), for sum",
932 "invalid input data");
936 unsigned int numParams = this->vectorSizeLocal();
937 for (
unsigned int i = 0; i < numParams; ++i) {
938 this->extractScalarSeq(initialPos,
944 data.autoCorrViaFft(0,
953 template <
class V,
class M>
956 unsigned int initialPos,
961 bool bRC = ((0 < numPos ) &&
962 ((initialPos+numPos) <= this->subSequenceSize()) &&
963 (this->vectorSizeLocal() == minVec.sizeLocal() ) &&
964 (this->vectorSizeLocal() == maxVec.sizeLocal() ));
967 "SequenceOfVectors<V,M>::subMinMaxExtra()",
968 "invalid input data");
971 unsigned int numParams = this->vectorSizeLocal();
974 for (
unsigned int i = 0; i < numParams; ++i) {
975 this->extractScalarSeq(initialPos,
986 template <
class V,
class M>
989 unsigned int initialPos,
992 V& unifiedMaxVec)
const
994 bool bRC = ((0 < numPos ) &&
995 ((initialPos+numPos) <= this->subSequenceSize() ) &&
996 (this->vectorSizeLocal() == unifiedMinVec.sizeLocal()) &&
997 (this->vectorSizeLocal() == unifiedMaxVec.sizeLocal()));
1000 "SequenceOfVectors<V,M>::unifiedMinMaxExtra()",
1001 "invalid input data");
1004 unsigned int numParams = this->vectorSizeLocal();
1007 for (
unsigned int i = 0; i < numParams; ++i) {
1008 this->extractScalarSeq(initialPos,
1023 template <
class V,
class M>
1026 unsigned int initialPos,
1029 std::vector<V*>& centersForAllBins,
1030 std::vector<V*>& quanttsForAllBins)
const
1032 bool bRC = ((initialPos < this->subSequenceSize() ) &&
1033 (this->vectorSizeLocal() == minVec.sizeLocal() ) &&
1034 (this->vectorSizeLocal() == maxVec.sizeLocal() ) &&
1035 (0 < centersForAllBins.size()) &&
1036 (centersForAllBins.size() == quanttsForAllBins.size()));
1039 "SequenceOfVectors<V,M>::subHistogram()",
1040 "invalid input data");
1042 for (
unsigned int j = 0; j < quanttsForAllBins.size(); ++j) {
1043 centersForAllBins[j] =
new V(m_vectorSpace.zeroVector());
1044 quanttsForAllBins [j] =
new V(m_vectorSpace.zeroVector());
1047 unsigned int dataSize = this->subSequenceSize() - initialPos;
1048 unsigned int numParams = this->vectorSizeLocal();
1049 for (
unsigned int i = 0; i < numParams; ++i) {
1051 for (
unsigned int j = 0; j < dataSize; ++j) {
1052 data[j] = (*(m_seq[initialPos+j]))[i];
1055 std::vector<double > centers(centersForAllBins.size(),0.);
1056 std::vector<unsigned int> quantts(quanttsForAllBins.size(), 0 );
1063 for (
unsigned int j = 0; j < quantts.size(); ++j) {
1064 (*(centersForAllBins[j]))[i] = centers[j];
1065 (*(quanttsForAllBins[j]))[i] = (double) quantts[j];
1072 template <
class V,
class M>
1075 unsigned int initialPos,
1076 const V& unifiedMinVec,
1077 const V& unifiedMaxVec,
1078 std::vector<V*>& unifiedCentersForAllBins,
1079 std::vector<V*>& unifiedQuanttsForAllBins)
const
1081 bool bRC = ((initialPos < this->subSequenceSize() ) &&
1082 (this->vectorSizeLocal() == unifiedMinVec.sizeLocal() ) &&
1083 (this->vectorSizeLocal() == unifiedMaxVec.sizeLocal() ) &&
1084 (0 < unifiedCentersForAllBins.size()) &&
1085 (unifiedCentersForAllBins.size() == unifiedQuanttsForAllBins.size()));
1088 "SequenceOfVectors<V,M>::unifiedHistogram()",
1089 "invalid input data");
1091 for (
unsigned int j = 0; j < unifiedQuanttsForAllBins.size(); ++j) {
1092 unifiedCentersForAllBins[j] =
new V(m_vectorSpace.zeroVector());
1093 unifiedQuanttsForAllBins [j] =
new V(m_vectorSpace.zeroVector());
1096 unsigned int dataSize = this->subSequenceSize() - initialPos;
1097 unsigned int numParams = this->vectorSizeLocal();
1098 for (
unsigned int i = 0; i < numParams; ++i) {
1100 for (
unsigned int j = 0; j < dataSize; ++j) {
1101 data[j] = (*(m_seq[initialPos+j]))[i];
1104 std::vector<double > unifiedCenters(unifiedCentersForAllBins.size(),0.);
1105 std::vector<unsigned int> unifiedQuantts(unifiedQuanttsForAllBins.size(), 0 );
1113 for (
unsigned int j = 0; j < unifiedQuantts.size(); ++j) {
1114 (*(unifiedCentersForAllBins[j]))[i] = unifiedCenters[j];
1115 (*(unifiedQuanttsForAllBins[j]))[i] = (double) unifiedQuantts[j];
1122 template <
class V,
class M>
1125 unsigned int initialPos,
1128 bool bRC = ((initialPos < this->subSequenceSize()) &&
1129 (this->vectorSizeLocal() == iqrVec.sizeLocal() ));
1132 "SequenceOfVectors<V,M>::subInterQuantileRange()",
1133 "invalid input data");
1135 unsigned int numPos = this->subSequenceSize() - initialPos;
1138 unsigned int numParams = this->vectorSizeLocal();
1139 for (
unsigned int i = 0; i < numParams; ++i) {
1140 this->extractScalarSeq(initialPos,
1151 template <
class V,
class M>
1154 unsigned int initialPos,
1155 V& unifiedIqrVec)
const
1157 bool bRC = ((initialPos < this->subSequenceSize() ) &&
1158 (this->vectorSizeLocal() == unifiedIqrVec.sizeLocal()));
1161 "SequenceOfVectors<V,M>::unifiedInterQuantileRange()",
1162 "invalid input data");
1164 unsigned int numPos = this->subSequenceSize() - initialPos;
1167 unsigned int numParams = this->vectorSizeLocal();
1168 for (
unsigned int i = 0; i < numParams; ++i) {
1169 this->extractScalarSeq(initialPos,
1181 template <
class V,
class M>
1184 unsigned int initialPos,
1186 unsigned int kdeDimension,
1189 bool bRC = ((initialPos < this->subSequenceSize()) &&
1190 (this->vectorSizeLocal() == iqrVec.sizeLocal() ) &&
1191 (this->vectorSizeLocal() == scaleVec.sizeLocal() ));
1194 "SequenceOfVectors<V,M>::subScalesForKde()",
1195 "invalid input data");
1197 unsigned int numPos = this->subSequenceSize() - initialPos;
1200 unsigned int numParams = this->vectorSizeLocal();
1201 for (
unsigned int i = 0; i < numParams; ++i) {
1202 this->extractScalarSeq(initialPos,
1215 template <
class V,
class M>
1218 unsigned int initialPos,
1219 const V& unifiedIqrVec,
1220 unsigned int kdeDimension,
1221 V& unifiedScaleVec)
const
1223 bool bRC = ((initialPos < this->subSequenceSize() ) &&
1224 (this->vectorSizeLocal() == unifiedIqrVec.sizeLocal() ) &&
1225 (this->vectorSizeLocal() == unifiedScaleVec.sizeLocal()));
1228 "SequenceOfVectors<V,M>::unifiedScalesForKde()",
1229 "invalid input data");
1231 unsigned int numPos = this->subSequenceSize() - initialPos;
1234 unsigned int numParams = this->vectorSizeLocal();
1235 for (
unsigned int i = 0; i < numParams; ++i) {
1236 this->extractScalarSeq(initialPos,
1241 unifiedScaleVec[i] = data.
unifiedScaleForKde(m_vectorSpace.numOfProcsForStorage() == 1,
1250 template <
class V,
class M>
1253 unsigned int initialPos,
1255 const std::vector<V*>& evalParamVecs,
1256 std::vector<V*>& densityVecs)
const
1258 bool bRC = ((initialPos < this->subSequenceSize()) &&
1259 (this->vectorSizeLocal() == scaleVec.sizeLocal() ) &&
1260 (0 < evalParamVecs.size() ) &&
1261 (evalParamVecs.size() == densityVecs.size() ));
1264 "SequenceOfVectors<V,M>::subGaussian1dKde()",
1265 "invalid input data");
1267 unsigned int numPos = this->subSequenceSize() - initialPos;
1270 unsigned int numEvals = evalParamVecs.size();
1271 for (
unsigned int j = 0; j < numEvals; ++j) {
1272 densityVecs[j] =
new V(m_vectorSpace.zeroVector());
1274 std::vector<double> evalParams(numEvals,0.);
1275 std::vector<double> densities (numEvals,0.);
1277 unsigned int numParams = this->vectorSizeLocal();
1278 for (
unsigned int i = 0; i < numParams; ++i) {
1279 this->extractScalarSeq(initialPos,
1285 for (
unsigned int j = 0; j < numEvals; ++j) {
1286 evalParams[j] = (*evalParamVecs[j])[i];
1294 for (
unsigned int j = 0; j < numEvals; ++j) {
1295 (*densityVecs[j])[i] = densities[j];
1302 template <
class V,
class M>
1305 unsigned int initialPos,
1306 const V& unifiedScaleVec,
1307 const std::vector<V*>& unifiedEvalParamVecs,
1308 std::vector<V*>& unifiedDensityVecs)
const
1310 bool bRC = ((initialPos < this->subSequenceSize() ) &&
1311 (this->vectorSizeLocal() == unifiedScaleVec.sizeLocal()) &&
1312 (0 < unifiedEvalParamVecs.size()) &&
1313 (unifiedEvalParamVecs.size() == unifiedDensityVecs.size() ));
1316 "SequenceOfVectors<V,M>::unifiedGaussian1dKde()",
1317 "invalid input data");
1319 unsigned int numPos = this->subSequenceSize() - initialPos;
1322 unsigned int numEvals = unifiedEvalParamVecs.size();
1323 for (
unsigned int j = 0; j < numEvals; ++j) {
1324 unifiedDensityVecs[j] =
new V(m_vectorSpace.zeroVector());
1326 std::vector<double> unifiedEvalParams(numEvals,0.);
1327 std::vector<double> unifiedDensities (numEvals,0.);
1329 unsigned int numParams = this->vectorSizeLocal();
1330 for (
unsigned int i = 0; i < numParams; ++i) {
1331 this->extractScalarSeq(initialPos,
1337 for (
unsigned int j = 0; j < numEvals; ++j) {
1338 unifiedEvalParams[j] = (*unifiedEvalParamVecs[j])[i];
1347 for (
unsigned int j = 0; j < numEvals; ++j) {
1348 (*unifiedDensityVecs[j])[i] = unifiedDensities[j];
1355 template <
class V,
class M>
1358 unsigned int initialPos,
1359 unsigned int numPos,
1360 const std::string& fileName,
1361 const std::string& fileType,
1362 const std::set<unsigned int>& allowedSubEnvIds)
const
1366 "SequenceOfVectors<V,M>::subWriteContents(1)",
1367 "unexpected subRank");
1370 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 99)) {
1371 *m_env.subDisplayFile() <<
"In SequenceOfVectors<V,M>::subWriteContents()"
1372 <<
": about to try to open file '" << fileName <<
"." << fileType
1374 <<
", initialPos = " << initialPos
1375 <<
", numPos = " << numPos
1378 if (m_env.openOutputFile(fileName,
1384 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 99)) {
1385 *m_env.subDisplayFile() <<
"In SequenceOfVectors<V,M>::subWriteContents()"
1386 <<
": successfully opened file '" << fileName <<
"." << fileType
1390 this->subWriteContents(initialPos,
1394 m_env.closeFile(filePtrSet,fileType);
1396 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 99)) {
1397 *m_env.subDisplayFile() <<
"In SequenceOfVectors<V,M>::subWriteContents()"
1398 <<
": before Barrier()"
1401 m_env.subComm().Barrier();
1406 template <
class V,
class M>
1409 unsigned int initialPos,
1410 unsigned int numPos,
1412 const std::string& fileType)
const
1416 "SequenceOfVectors<V,M>::subWriteContents(2)",
1417 "filePtrSet.ofsVar should not be NULL");
1420 this->subWriteContents(initialPos,
1428 "SequenceOfVectors<V,M>::subWriteContents(2)",
1429 "hdf file type not supported yet");
1434 "SequenceOfVectors<V,M>::subWriteContents(2)",
1435 "invalid file type");
1441 template <
class V,
class M>
1444 unsigned int initialPos,
1445 unsigned int numPos,
1447 const std::string& fileType)
const
1451 "SequenceOfVectors<V,M>::subWriteContents(3)",
1452 "invalid routine input parameters");
1454 if (fileType.c_str()) {};
1456 if (initialPos == 0) {
1457 ofs << m_name <<
"_sub" << m_env.subIdString() <<
" = zeros(" << this->subSequenceSize()
1458 <<
"," << this->vectorSizeLocal()
1461 ofs << m_name <<
"_sub" << m_env.subIdString() <<
" = [";
1464 for (
unsigned int j = initialPos; j < initialPos+numPos; ++j) {
1465 bool savedVectorPrintScientific = m_seq[j]->getPrintScientific();
1466 bool savedVectorPrintState = m_seq[j]->getPrintHorizontally();
1467 m_seq[j]->setPrintScientific (
true);
1468 m_seq[j]->setPrintHorizontally(
true);
1473 m_seq[j]->setPrintHorizontally(savedVectorPrintState);
1474 m_seq[j]->setPrintScientific (savedVectorPrintScientific);
1476 if ((initialPos+numPos) == this->subSequenceSize()) {
1483 template <
class V,
class M>
1486 const std::string& fileName,
1487 const std::string& inputFileType)
const
1489 std::string fileType(inputFileType);
1490 #ifdef QUESO_HAS_HDF5
1494 if (m_env.subDisplayFile()) {
1495 *m_env.subDisplayFile() <<
"WARNING in SequenceOfVectors<V,M>::unifiedWriteContents()"
1497 <<
"' has been requested, but this QUESO library has not been built with 'hdf5'"
1502 if (m_env.subRank() == 0) {
1503 std::cerr <<
"WARNING in SequenceOfVectors<V,M>::unifiedWriteContents()"
1505 <<
"' has been requested, but this QUESO library has not been built with 'hdf5'"
1517 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 10)) {
1518 *m_env.subDisplayFile() <<
"Entering SequenceOfVectors<V,M>::unifiedWriteContents()"
1519 <<
": worldRank " << m_env.worldRank()
1520 <<
", fullRank " << m_env.fullRank()
1521 <<
", subEnvironment " << m_env.subId()
1522 <<
", subRank " << m_env.subRank()
1523 <<
", inter0Rank " << m_env.inter0Rank()
1525 <<
", fileName = " << fileName
1529 if (m_env.inter0Rank() >= 0) {
1530 for (
unsigned int r = 0; r < (
unsigned int) m_env.inter0Comm().NumProc(); ++r) {
1531 if (m_env.inter0Rank() == (int) r) {
1533 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 10)) {
1534 *m_env.subDisplayFile() <<
"In SequenceOfVectors<V,M>::unifiedWriteContents()"
1535 <<
": worldRank " << m_env.worldRank()
1536 <<
", fullRank " << m_env.fullRank()
1537 <<
", subEnvironment " << m_env.subId()
1538 <<
", subRank " << m_env.subRank()
1539 <<
", inter0Rank " << m_env.inter0Rank()
1541 <<
", fileName = " << fileName
1542 <<
", about to open file for r = " << r
1547 bool writeOver =
false;
1550 if (m_env.openUnifiedOutputFile(fileName,
1553 unifiedFilePtrSet)) {
1554 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 10)) {
1555 *m_env.subDisplayFile() <<
"In SequenceOfVectors<V,M>::unifiedWriteContents()"
1556 <<
": worldRank " << m_env.worldRank()
1557 <<
", fullRank " << m_env.fullRank()
1558 <<
", subEnvironment " << m_env.subId()
1559 <<
", subRank " << m_env.subRank()
1560 <<
", inter0Rank " << m_env.inter0Rank()
1562 <<
", fileName = " << fileName
1563 <<
", just opened file for r = " << r
1567 unsigned int chainSize = this->subSequenceSize();
1570 *unifiedFilePtrSet.
ofsVar << m_name <<
"_unified" <<
" = zeros(" << this->subSequenceSize()*m_env.inter0Comm().NumProc()
1571 <<
"," << this->vectorSizeLocal()
1574 *unifiedFilePtrSet.
ofsVar << m_name <<
"_unified" <<
" = [";
1577 for (
unsigned int j = 0; j < chainSize; ++j) {
1590 bool savedVectorPrintScientific = m_seq[j]->getPrintScientific();
1591 bool savedVectorPrintState = m_seq[j]->getPrintHorizontally();
1592 m_seq[j]->setPrintScientific (
true);
1593 m_seq[j]->setPrintHorizontally(
true);
1595 *unifiedFilePtrSet.
ofsVar << *(m_seq[j])
1598 m_seq[j]->setPrintHorizontally(savedVectorPrintState);
1599 m_seq[j]->setPrintScientific (savedVectorPrintScientific);
1602 #ifdef QUESO_HAS_HDF5
1604 unsigned int numParams = m_vectorSpace.dimLocal();
1606 hid_t datatype = H5Tcopy(H5T_NATIVE_DOUBLE);
1609 dimsf[0] = numParams;
1610 dimsf[1] = chainSize;
1611 hid_t dataspace = H5Screate_simple(2, dimsf, NULL);
1613 hid_t dataset = H5Dcreate2(unifiedFilePtrSet.h5Var,
1622 struct timeval timevalBegin;
1624 iRC = gettimeofday(&timevalBegin,NULL);
1628 std::vector<double*> dataOut((
size_t) numParams,NULL);
1629 dataOut[0] = (
double*) malloc(numParams*chainSize*
sizeof(
double));
1630 for (
unsigned int i = 1; i < numParams; ++i) {
1631 dataOut[i] = dataOut[i-1] + chainSize;
1634 for (
unsigned int j = 0; j < chainSize; ++j) {
1635 V tmpVec(*(m_seq[j]));
1636 for (
unsigned int i = 0; i < numParams; ++i) {
1637 dataOut[i][j] = tmpVec[i];
1643 status = H5Dwrite(dataset,
1648 (
void*) dataOut[0]);
1653 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 2)) {
1654 *m_env.subDisplayFile() <<
"In SequenceOfVectors<V,M>::unifiedWriteContents()"
1655 <<
": worldRank " << m_env.worldRank()
1656 <<
", fullRank " << m_env.fullRank()
1657 <<
", subEnvironment " << m_env.subId()
1658 <<
", subRank " << m_env.subRank()
1659 <<
", inter0Rank " << m_env.inter0Rank()
1660 <<
", fileName = " << fileName
1661 <<
", numParams = " << numParams
1662 <<
", chainSize = " << chainSize
1663 <<
", writeTime = " << writeTime <<
" seconds"
1669 H5Sclose(dataspace);
1674 for (
unsigned int tmpIndex = 0; tmpIndex < dataOut.size(); tmpIndex++) {
1675 free (dataOut[tmpIndex]);
1681 "SequenceOfVectors<V,M>::unifiedWriteContents()",
1682 "hdf file type not supported for multiple sub-environments yet");
1689 "SequenceOfVectors<V,M>::unifiedWriteContents()",
1690 "invalid file type");
1693 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 10)) {
1694 *m_env.subDisplayFile() <<
"In SequenceOfVectors<V,M>::unifiedWriteContents()"
1695 <<
": worldRank " << m_env.worldRank()
1696 <<
", fullRank " << m_env.fullRank()
1697 <<
", subEnvironment " << m_env.subId()
1698 <<
", subRank " << m_env.subRank()
1699 <<
", inter0Rank " << m_env.inter0Rank()
1701 <<
", fileName = " << fileName
1702 <<
", about to close file for r = " << r
1706 m_env.closeFile(unifiedFilePtrSet,fileType);
1708 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 10)) {
1709 *m_env.subDisplayFile() <<
"In SequenceOfVectors<V,M>::unifiedWriteContents()"
1710 <<
": worldRank " << m_env.worldRank()
1711 <<
", fullRank " << m_env.fullRank()
1712 <<
", subEnvironment " << m_env.subId()
1713 <<
", subRank " << m_env.subRank()
1714 <<
", inter0Rank " << m_env.inter0Rank()
1716 <<
", fileName = " << fileName
1717 <<
", just closed file for r = " << r
1722 m_env.inter0Comm().Barrier();
1725 if (m_env.inter0Rank() == 0) {
1728 if (m_env.openUnifiedOutputFile(fileName,
1731 unifiedFilePtrSet)) {
1732 *unifiedFilePtrSet.
ofsVar <<
"];\n";
1733 m_env.closeFile(unifiedFilePtrSet,fileType);
1742 "SequenceOfVectors<V,M>::unifiedWriteContents(), final",
1743 "invalid file type");
1748 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 10)) {
1749 *m_env.subDisplayFile() <<
"Leaving SequenceOfVectors<V,M>::unifiedWriteContents()"
1750 <<
", fileName = " << fileName
1758 template <
class V,
class M>
1761 const std::string& fileName,
1762 const std::string& inputFileType,
1763 const unsigned int subReadSize)
1765 std::string fileType(inputFileType);
1766 #ifdef QUESO_HAS_HDF5
1770 if (m_env.subDisplayFile()) {
1771 *m_env.subDisplayFile() <<
"WARNING in SequenceOfVectors<V,M>::unifiedReadContents()"
1773 <<
"' has been requested, but this QUESO library has not been built with 'hdf5'"
1778 if (m_env.subRank() == 0) {
1779 std::cerr <<
"WARNING in SequenceOfVectors<V,M>::unifiedReadContents()"
1781 <<
"' has been requested, but this QUESO library has not been built with 'hdf5'"
1791 if (m_env.subDisplayFile()) {
1792 *m_env.subDisplayFile() <<
"Entering SequenceOfVectors<V,M>::unifiedReadContents()"
1793 <<
": worldRank " << m_env.worldRank()
1794 <<
", fullRank " << m_env.fullRank()
1795 <<
", subEnvironment " << m_env.subId()
1796 <<
", subRank " << m_env.subRank()
1797 <<
", inter0Rank " << m_env.inter0Rank()
1799 <<
", fileName = " << fileName
1800 <<
", subReadSize = " << subReadSize
1805 this->resizeSequence(subReadSize);
1807 if (m_env.inter0Rank() >= 0) {
1808 double unifiedReadSize = subReadSize*m_env.inter0Comm().NumProc();
1811 unsigned int idOfMyFirstLine = 1 + m_env.inter0Rank()*subReadSize;
1812 unsigned int idOfMyLastLine = (1 + m_env.inter0Rank())*subReadSize;
1813 unsigned int numParams = this->vectorSizeLocal();
1815 for (
unsigned int r = 0; r < (
unsigned int) m_env.inter0Comm().NumProc(); ++r) {
1816 if (m_env.inter0Rank() == (int) r) {
1819 if (m_env.openUnifiedInputFile(fileName,
1821 unifiedFilePtrSet)) {
1826 std::string tmpString;
1829 *unifiedFilePtrSet.
ifsVar >> tmpString;
1833 *unifiedFilePtrSet.
ifsVar >> tmpString;
1837 "SequenceOfVectors<V,M>::unifiedReadContents()",
1838 "string should be the '=' sign");
1841 *unifiedFilePtrSet.
ifsVar >> tmpString;
1843 unsigned int posInTmpString = 6;
1847 std::string nPositionsString((
size_t) (tmpString.size()-posInTmpString+1),
' ');
1848 unsigned int posInPositionsString = 0;
1852 "SequenceOfVectors<V,M>::unifiedReadContents()",
1853 "symbol ',' not found in first line of file");
1854 nPositionsString[posInPositionsString++] = tmpString[posInTmpString++];
1855 }
while (tmpString[posInTmpString] !=
',');
1856 nPositionsString[posInPositionsString] =
'\0';
1861 std::string nParamsString((
size_t) (tmpString.size()-posInTmpString+1),
' ');
1862 unsigned int posInParamsString = 0;
1866 "SequenceOfVectors<V,M>::unifiedReadContents()",
1867 "symbol ')' not found in first line of file");
1868 nParamsString[posInParamsString++] = tmpString[posInTmpString++];
1869 }
while (tmpString[posInTmpString] !=
')');
1870 nParamsString[posInParamsString] =
'\0';
1873 unsigned int sizeOfChainInFile = (
unsigned int) strtod(nPositionsString.c_str(),NULL);
1874 unsigned int numParamsInFile = (
unsigned int) strtod(nParamsString.c_str(), NULL);
1875 if (m_env.subDisplayFile()) {
1876 *m_env.subDisplayFile() <<
"In SequenceOfVectors<V,M>::unifiedReadContents()"
1877 <<
": worldRank " << m_env.worldRank()
1878 <<
", fullRank " << m_env.fullRank()
1879 <<
", sizeOfChainInFile = " << sizeOfChainInFile
1880 <<
", numParamsInFile = " << numParamsInFile
1887 "SequenceOfVectors<V,M>::unifiedReadContents()",
1888 "size of chain in file is not big enough");
1893 "SequenceOfVectors<V,M>::unifiedReadContents()",
1894 "number of parameters of chain in file is different than number of parameters in this chain object");
1898 unsigned int maxCharsPerLine = 64*numParams;
1900 unsigned int lineId = 0;
1901 while (lineId < idOfMyFirstLine) {
1902 unifiedFilePtrSet.
ifsVar->ignore(maxCharsPerLine,
'\n');
1909 std::string tmpString;
1912 *unifiedFilePtrSet.
ifsVar >> tmpString;
1916 *unifiedFilePtrSet.
ifsVar >> tmpString;
1920 "SequenceOfVectors<V,M>::unifiedReadContents()",
1921 "in core 0, string should be the '=' sign");
1924 std::streampos tmpPos = unifiedFilePtrSet.
ifsVar->tellg();
1925 unifiedFilePtrSet.
ifsVar->seekg(tmpPos+(std::streampos)2);
1928 V tmpVec(m_vectorSpace.zeroVector());
1929 while (lineId <= idOfMyLastLine) {
1930 for (
unsigned int i = 0; i < numParams; ++i) {
1931 *unifiedFilePtrSet.
ifsVar >> tmpVec[i];
1933 this->setPositionValues(lineId - idOfMyFirstLine, tmpVec);
1937 #ifdef QUESO_HAS_HDF5
1940 hid_t dataset = H5Dopen2(unifiedFilePtrSet.h5Var,
1943 hid_t datatype = H5Dget_type(dataset);
1944 H5T_class_t t_class = H5Tget_class(datatype);
1947 "SequenceOfVectors<V,M>::unifiedReadContents()",
1948 "t_class is not H5T_DOUBLE");
1949 hid_t dataspace = H5Dget_space(dataset);
1950 int rank = H5Sget_simple_extent_ndims(dataspace);
1953 "SequenceOfVectors<V,M>::unifiedReadContents()",
1954 "hdf rank is not 2");
1957 status_n = H5Sget_simple_extent_dims(dataspace, dims_in, NULL);
1965 "SequenceOfVectors<V,M>::unifiedReadContents()",
1966 "dims_in[0] is not equal to 'numParams'");
1969 "SequenceOfVectors<V,M>::unifiedReadContents()",
1970 "dims_in[1] is smaller that requested 'subReadSize'");
1972 struct timeval timevalBegin;
1974 iRC = gettimeofday(&timevalBegin,NULL);
1977 unsigned int chainSizeIn = (
unsigned int) dims_in[1];
1979 std::vector<double*> dataIn((
size_t) numParams,NULL);
1980 dataIn[0] = (
double*) malloc(numParams*chainSizeIn*
sizeof(
double));
1981 for (
unsigned int i = 1; i < numParams; ++i) {
1982 dataIn[i] = dataIn[i-1] + chainSizeIn;
1986 status = H5Dread(dataset,
1994 V tmpVec(m_vectorSpace.zeroVector());
1995 for (
unsigned int j = 0; j < subReadSize; ++j) {
1996 for (
unsigned int i = 0; i < numParams; ++i) {
1997 tmpVec[i] = dataIn[i][j];
1999 this->setPositionValues(j, tmpVec);
2003 if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 2)) {
2004 *m_env.subDisplayFile() <<
"In SequenceOfVectors<V,M>::unifiedReadContents()"
2005 <<
": worldRank " << m_env.worldRank()
2006 <<
", fullRank " << m_env.fullRank()
2007 <<
", subEnvironment " << m_env.subId()
2008 <<
", subRank " << m_env.subRank()
2009 <<
", inter0Rank " << m_env.inter0Rank()
2010 <<
", fileName = " << fileName
2011 <<
", numParams = " << numParams
2012 <<
", chainSizeIn = " << chainSizeIn
2013 <<
", subReadSize = " << subReadSize
2014 <<
", readTime = " << readTime <<
" seconds"
2018 H5Sclose(dataspace);
2022 for (
unsigned int tmpIndex = 0; tmpIndex < dataIn.size(); tmpIndex++) {
2023 free (dataIn[tmpIndex]);
2029 "SequenceOfVectors<V,M>::unifiedReadContents()",
2030 "hdf file type not supported for multiple sub-environments yet");
2037 "SequenceOfVectors<V,M>::unifiedReadContents()",
2038 "invalid file type");
2040 m_env.closeFile(unifiedFilePtrSet,fileType);
2043 m_env.inter0Comm().Barrier();
2047 V tmpVec(m_vectorSpace.zeroVector());
2048 for (
unsigned int i = 1; i < subReadSize; ++i) {
2049 this->setPositionValues(i,tmpVec);
2053 if (m_env.subDisplayFile()) {
2054 *m_env.subDisplayFile() <<
"Leaving SequenceOfVectors<V,M>::unifiedReadContents()"
2055 <<
", fileName = " << fileName
2063 template <
class V,
class M>
2069 "SequenceOfVectors<V,M>::select()",
2070 "Code is not complete yet");
2072 if (&idsOfUniquePositions) {};
2077 template <
class V,
class M>
2080 unsigned int initialPos,
2081 unsigned int spacing)
2083 if (m_env.subDisplayFile()) {
2084 *m_env.subDisplayFile() <<
"Entering SequenceOfVectors<V,M>::filter()"
2085 <<
": initialPos = " << initialPos
2086 <<
", spacing = " << spacing
2087 <<
", subSequenceSize = " << this->subSequenceSize()
2092 unsigned int j = initialPos;
2093 unsigned int originalSubSequenceSize = this->subSequenceSize();
2094 while (j < originalSubSequenceSize) {
2098 m_seq[i] =
new V(*(m_seq[j]));
2104 this->resetValues(i,originalSubSequenceSize-i);
2105 this->resizeSequence(i);
2107 if (m_env.subDisplayFile()) {
2108 *m_env.subDisplayFile() <<
"Leaving SequenceOfVectors<V,M>::filter()"
2109 <<
": initialPos = " << initialPos
2110 <<
", spacing = " << spacing
2111 <<
", subSequenceSize = " << this->subSequenceSize()
2118 template <
class V,
class M>
2121 unsigned int initialPos,
2122 unsigned int numPos)
const
2127 "SequenceOfVectors<V,M>::estimateConvBrooksGelman()",
2128 "At least two sequences required for Brooks-Gelman convergence test.");
2133 double convMeasure = -1.0;
2137 if( m_env.inter0Rank() >= 0 )
2143 V psi_j_dot = m_vectorSpace.zeroVector();
2144 V psi_dot_dot = m_vectorSpace.zeroVector();
2145 V work = m_vectorSpace.zeroVector();
2149 int m = m_env.numSubEnvironments();
2152 this->subMeanExtra ( initialPos, numPos, psi_j_dot );
2153 this->unifiedMeanExtra( initialPos, numPos, psi_dot_dot );
2156 std::cout <<
"psi_j_dot = " << psi_j_dot << std::endl;
2157 std::cout <<
"psi_dot_dot = " << psi_dot_dot << std::endl;
2163 M* W_local = m_vectorSpace.newDiagMatrix( m_vectorSpace.zeroVector() );
2164 M* W = m_vectorSpace.newDiagMatrix( m_vectorSpace.zeroVector() );
2165 V psi_j_t = m_vectorSpace.zeroVector();
2168 for(
unsigned int t = initialPos; t < initialPos+numPos; ++t )
2170 psi_j_t = *(m_seq[t]);
2172 work = psi_j_t - psi_j_dot;
2179 W_local->mpiSum( m_env.inter0Comm(), (*W) );
2181 (*W) = 1.0/(double(m)*(double(n)-1.0)) * (*W);
2184 std::cout <<
"n, m = " << n <<
", " << m << std::endl;
2185 std::cout <<
"W_local = " << *W_local << std::endl;
2186 std::cout <<
"W = " << *W << std::endl;
2196 M* B_over_n_local = m_vectorSpace.newDiagMatrix( m_vectorSpace.zeroVector() );
2197 M* B_over_n = m_vectorSpace.newDiagMatrix( m_vectorSpace.zeroVector() );
2199 work = psi_j_dot - psi_dot_dot;
2202 B_over_n_local->mpiSum( m_env.inter0Comm(), (*B_over_n) );
2205 delete B_over_n_local;
2207 (*B_over_n) = 1.0/(double(m)-1.0) * (*B_over_n);
2210 std::cout <<
"B_over_n = " << *B_over_n << std::endl;
2216 M* A = m_vectorSpace.newDiagMatrix( m_vectorSpace.zeroVector() );
2218 W->invertMultiply( *B_over_n, *A );
2221 std::cout <<
"A = " << *A << std::endl;
2229 V eigenVector = m_vectorSpace.zeroVector();
2231 A->largestEigen( eigenValue, eigenVector );
2237 convMeasure = (double(n)-1.0)/double(n) + (double(m)+1.0)/double(m)*eigenValue;
2248 template <
class V,
class M>
2251 unsigned int initialPos,
2252 unsigned int spacing,
2253 unsigned int numPos,
2254 unsigned int paramId,
2259 for (
unsigned int j = 0; j < numPos; ++j) {
2260 scalarSeq[j] = (*(m_seq[initialPos+j ]))[paramId];
2264 for (
unsigned int j = 0; j < numPos; ++j) {
2265 scalarSeq[j] = (*(m_seq[initialPos+j*spacing]))[paramId];
2272 template <
class V,
class M>
2277 for (
unsigned int i = 0; i < (
unsigned int) m_seq.size(); ++i) {
2284 for (
unsigned int i = 0; i < m_seq.size(); ++i) {
2285 m_seq[i] =
new V(*(src.
m_seq[i]));
2291 template <
class V,
class M>
2294 unsigned int initialPos,
2295 unsigned int spacing,
2296 unsigned int numPos,
2297 unsigned int paramId,
2298 std::vector<double>& rawData)
const
2300 rawData.resize(numPos);
2302 for (
unsigned int j = 0; j < numPos; ++j) {
2303 rawData[j] = (*(m_seq[initialPos+j ]))[paramId];
2307 for (
unsigned int j = 0; j < numPos; ++j) {
2308 rawData[j] = (*(m_seq[initialPos+j*spacing]))[paramId];
2319 #ifdef UQ_SEQ_VEC_USES_OPERATOR
2320 template <
class V,
class M>
2326 "SequenceOfVectorss<V,M>::operator[] const",
2327 "posId > subSequenceSize()");
2329 return (
const V*) (m_seq[posId]);
2332 template <
class V,
class M>
2334 SequenceOfVectors<V,M>::operator[](
unsigned int posId)
2338 "SequenceOfVectorss<V,M>::operator[] const",
2339 "posId > subSequenceSize()");
2341 return m_seq[posId];
2349 #ifdef UQ_ALSO_COMPUTE_MDFS_WITHOUT_KDE
2350 template <
class V,
class M>
2352 SequenceOfVectors<V,M>::subUniformlySampledMdf(
2353 const V& numEvaluationPointsVec,
2354 ArrayOfOneDGrids <V,M>& mdfGrids,
2355 ArrayOfOneDTables<V,M>& mdfValues)
const
2357 V minDomainValues(m_vectorSpace.zeroVector());
2358 V maxDomainValues(m_vectorSpace.zeroVector());
2360 ScalarSequence<double> data(m_env,0,
"");
2362 unsigned int numParams = this->vectorSizeLocal();
2363 for (
unsigned int i = 0; i < numParams; ++i) {
2364 this->extractScalarSeq(0,
2370 std::vector<double> aMdf(0);
2371 data.subUniformlySampledMdf((
unsigned int) numEvaluationPointsVec[i],
2375 mdfValues.setOneDTable(i,aMdf);
2378 mdfGrids.setUniformGrids(numEvaluationPointsVec,
2390 #ifdef QUESO_COMPUTES_EXTRA_POST_PROCESSING_STATISTICS
2391 template <
class V,
class M>
2393 SequenceOfVectors<V,M>::subMeanCltStd(
2394 unsigned int initialPos,
2395 unsigned int numPos,
2399 bool bRC = ((initialPos < this->subSequenceSize()) &&
2401 ((initialPos+numPos) <= this->subSequenceSize()) &&
2402 (this->vectorSizeLocal() == meanVec.sizeLocal() ) &&
2403 (this->vectorSizeLocal() == stdVec.sizeLocal() ));
2406 "SequenceOfVectors<V,M>::subMeanCltStd()",
2407 "invalid input data");
2409 ScalarSequence<double> data(m_env,0,
"");
2411 unsigned int numParams = this->vectorSizeLocal();
2412 for (
unsigned int i = 0; i < numParams; ++i) {
2413 this->extractScalarSeq(initialPos,
2418 stdVec[i] = data.subMeanCltStd(0,
2426 template <
class V,
class M>
2428 SequenceOfVectors<V,M>::unifiedMeanCltStd(
2429 unsigned int initialPos,
2430 unsigned int numPos,
2431 const V& unifiedMeanVec,
2432 V& unifiedSamVec)
const
2434 bool bRC = ((initialPos < this->subSequenceSize() ) &&
2436 ((initialPos+numPos) <= this->subSequenceSize() ) &&
2437 (this->vectorSizeLocal() == unifiedMeanVec.sizeLocal()) &&
2438 (this->vectorSizeLocal() == unifiedSamVec.sizeLocal() ));
2441 "SequenceOfVectors<V,M>::unifiedMeanCltStd()",
2442 "invalid input data");
2444 ScalarSequence<double> data(m_env,0,
"");
2446 unsigned int numParams = this->vectorSizeLocal();
2447 for (
unsigned int i = 0; i < numParams; ++i) {
2448 this->extractScalarSeq(initialPos,
2453 unifiedSamVec[i] = data.unifiedMeanCltStd(m_vectorSpace.numOfProcsForStorage() == 1,
2462 template <
class V,
class M>
2464 SequenceOfVectors<V,M>::bmm(
2465 unsigned int initialPos,
2466 unsigned int batchLength,
2469 bool bRC = ((initialPos < this->subSequenceSize() ) &&
2470 (batchLength < (this->subSequenceSize()-initialPos)) &&
2471 (this->vectorSizeLocal() == bmmVec.sizeLocal() ));
2474 "SequenceOfVectors<V,M>::bmm()",
2475 "invalid input data");
2477 ScalarSequence<double> data(m_env,0,
"");
2479 unsigned int numParams = this->vectorSizeLocal();
2480 for (
unsigned int i = 0; i < numParams; ++i) {
2481 this->extractScalarSeq(initialPos,
2483 this->subSequenceSize()-initialPos,
2486 bmmVec[i] = data.bmm(0,
2493 template <
class V,
class M>
2495 SequenceOfVectors<V,M>::fftForward(
2496 unsigned int initialPos,
2497 unsigned int fftSize,
2498 unsigned int paramId,
2499 std::vector<std::complex<double> >& fftResult)
const
2501 bool bRC = ((initialPos < this->subSequenceSize()) &&
2502 (paramId < this->vectorSizeLocal()) &&
2504 ((initialPos+fftSize) <= this->subSequenceSize()) &&
2505 (fftSize < this->subSequenceSize()));
2508 "SequenceOfVectors<V,M>::fftForward()",
2509 "invalid input data");
2511 std::vector<double> rawData(fftSize,0.);
2512 this->extractRawData(initialPos,
2518 m_fftObj->forward(rawData,fftSize,fftResult);
2523 template <
class V,
class M>
2525 SequenceOfVectors<V,M>::psd(
2526 unsigned int initialPos,
2527 unsigned int numBlocks,
2528 double hopSizeRatio,
2529 unsigned int paramId,
2530 std::vector<double>& psdResult)
const
2532 bool bRC = ((initialPos < this->subSequenceSize()) &&
2533 (paramId < this->vectorSizeLocal()));
2536 "SequenceOfVectors<V,M>::psd()",
2537 "invalid input data");
2539 ScalarSequence<double> data(m_env,0,
"");
2541 this->extractScalarSeq(initialPos,
2543 this->subSequenceSize()-initialPos,
2554 template <
class V,
class M>
2556 SequenceOfVectors<V,M>::psdAtZero(
2557 unsigned int initialPos,
2558 unsigned int numBlocks,
2559 double hopSizeRatio,
2562 bool bRC = ((initialPos < this->subSequenceSize()) &&
2563 (this->vectorSizeLocal() == psdVec.sizeLocal()));
2566 "SequenceOfVectors<V,M>::psdAtZero()",
2567 "invalid input data");
2569 ScalarSequence<double> data(m_env,0,
"");
2570 std::vector<double> psdResult(0,0.);
2572 unsigned int numParams = this->vectorSizeLocal();
2573 for (
unsigned int i = 0; i < numParams; ++i) {
2574 this->extractScalarSeq(initialPos,
2576 this->subSequenceSize()-initialPos,
2583 psdVec[i] = psdResult[0];
2590 template <
class V,
class M>
2592 SequenceOfVectors<V,M>::geweke(
2593 unsigned int initialPos,
2598 bool bRC = ((initialPos < this->subSequenceSize()) &&
2599 (this->vectorSizeLocal() == gewVec.sizeLocal() ));
2602 "SequenceOfVectors<V,M>::geweke()",
2603 "invalid input data");
2605 unsigned int numPos = this->subSequenceSize() - initialPos;
2606 ScalarSequence<double> data(m_env,0,
"");
2608 unsigned int numParams = this->vectorSizeLocal();
2609 for (
unsigned int i = 0; i < numParams; ++i) {
2610 this->extractScalarSeq(initialPos,
2615 gewVec[i] = data.geweke(0,
2623 template <
class V,
class M>
2625 SequenceOfVectors<V,M>::meanStacc(
2626 unsigned int initialPos,
2627 V& meanStaccVec)
const
2629 bool bRC = ((initialPos < this->subSequenceSize() ) &&
2630 (this->vectorSizeLocal() == meanStaccVec.sizeLocal()));
2633 "SequenceOfVectors<V,M>::meanStacc()",
2634 "invalid input data");
2636 unsigned int numPos = this->subSequenceSize() - initialPos;
2637 ScalarSequence<double> data(m_env,0,
"");
2639 unsigned int numParams = this->vectorSizeLocal();
2640 for (
unsigned int i = 0; i < numParams; ++i) {
2641 this->extractScalarSeq(initialPos,
2646 meanStaccVec[i] = data.meanStacc(0);
2652 template <
class V,
class M>
2654 SequenceOfVectors<V,M>::subCdfPercentageRange(
2655 unsigned int initialPos,
2656 unsigned int numPos,
2661 bool bRC = ((0 < numPos ) &&
2662 ((initialPos+numPos) <= this->subSequenceSize()) &&
2663 (this->vectorSizeLocal() == lowerVec.sizeLocal() ) &&
2664 (this->vectorSizeLocal() == upperVec.sizeLocal() ));
2667 "SequenceOfVectors<V,M>::subCdfPercentageRange()",
2668 "invalid input data");
2670 unsigned int numParams = this->vectorSizeLocal();
2671 ScalarSequence<double> data(m_env,0,
"");
2673 for (
unsigned int i = 0; i < numParams; ++i) {
2674 this->extractScalarSeq(initialPos,
2679 data.subCdfPercentageRange(0,
2689 template <
class V,
class M>
2691 SequenceOfVectors<V,M>::unifiedCdfPercentageRange(
2692 unsigned int initialPos,
2693 unsigned int numPos,
2698 bool bRC = ((0 < numPos ) &&
2699 ((initialPos+numPos) <= this->subSequenceSize()) &&
2700 (this->vectorSizeLocal() == lowerVec.sizeLocal() ) &&
2701 (this->vectorSizeLocal() == upperVec.sizeLocal() ));
2704 "SequenceOfVectors<V,M>::unifiedCdfPercentageRange()",
2705 "invalid input data");
2707 unsigned int numParams = this->vectorSizeLocal();
2708 ScalarSequence<double> data(m_env,0,
"");
2710 for (
unsigned int i = 0; i < numParams; ++i) {
2711 this->extractScalarSeq(initialPos,
2716 data.unifiedCdfPercentageRange(m_vectorSpace.numOfProcsForStorage() == 1,
2727 template <
class V,
class M>
2729 SequenceOfVectors<V,M>::subCdfStacc(
2730 unsigned int initialPos,
2731 std::vector<V*>& cdfStaccVecs,
2732 std::vector<V*>& cdfStaccVecsUp,
2733 std::vector<V*>& cdfStaccVecsLow,
2734 std::vector<V*>& sortedDataVecs)
const
2736 bool bRC = (initialPos < this->subSequenceSize());
2739 "SequenceOfVectors<V,M>::subCdfStacc()",
2740 "invalid input data");
2742 unsigned int numPos = this->subSequenceSize() - initialPos;
2743 unsigned int numEvals = numPos;
2744 for (
unsigned int j = 0; j < numEvals; ++j) {
2745 cdfStaccVecs [j] =
new V(m_vectorSpace.zeroVector());
2746 cdfStaccVecsUp [j] =
new V(m_vectorSpace.zeroVector());
2747 cdfStaccVecsLow[j] =
new V(m_vectorSpace.zeroVector());
2748 sortedDataVecs [j] =
new V(m_vectorSpace.zeroVector());
2750 std::vector<double> cdfStaccs (numEvals,0.);
2751 std::vector<double> cdfStaccsup (numEvals,0.);
2752 std::vector<double> cdfStaccslow(numEvals,0.);
2754 ScalarSequence<double> data (m_env,0,
"");
2755 ScalarSequence<double> sortedData(m_env,0,
"");
2756 unsigned int numParams = this->vectorSizeLocal();
2757 for (
unsigned int i = 0; i < numParams; ++i) {
2758 this->extractScalarSeq(initialPos,
2764 data.subSort(initialPos,sortedData);
2765 data.subCdfStacc(initialPos,
2771 for (
unsigned int j = 0; j < numEvals; ++j) {
2772 (*sortedDataVecs [j])[i] = sortedData [j];
2773 (*cdfStaccVecs [j])[i] = cdfStaccs [j];
2774 (*cdfStaccVecsUp [j])[i] = cdfStaccsup [j];
2775 (*cdfStaccVecsLow[j])[i] = cdfStaccslow[j];
2782 template <
class V,
class M>
2784 SequenceOfVectors<V,M>::subCdfStacc(
2785 unsigned int initialPos,
2786 const std::vector<V*>& evalPositionsVecs,
2787 std::vector<V*>& cdfStaccVecs)
const
2789 bool bRC = ((initialPos < this->subSequenceSize() ) &&
2790 (0 < evalPositionsVecs.size()) &&
2791 (evalPositionsVecs.size() == cdfStaccVecs.size() ));
2794 "SequenceOfVectors<V,M>::subCdfStacc()",
2795 "invalid input data");
2797 unsigned int numPos = this->subSequenceSize() - initialPos;
2798 ScalarSequence<double> data(m_env,0,
"");
2800 unsigned int numEvals = evalPositionsVecs.size();
2801 for (
unsigned int j = 0; j < numEvals; ++j) {
2802 cdfStaccVecs[j] =
new V(m_vectorSpace.zeroVector());
2804 std::vector<double> evalPositions(numEvals,0.);
2805 std::vector<double> cdfStaccs (numEvals,0.);
2807 unsigned int numParams = this->vectorSizeLocal();
2808 for (
unsigned int i = 0; i < numParams; ++i) {
2809 this->extractScalarSeq(initialPos,
2815 for (
unsigned int j = 0; j < numEvals; ++j) {
2816 evalPositions[j] = (*evalPositionsVecs[j])[i];
2823 for (
unsigned int j = 0; j < numEvals; ++j) {
2824 (*cdfStaccVecs[j])[i] = cdfStaccs[j];
2830 #endif // #ifdef QUESO_COMPUTES_EXTRA_POST_PROCESSING_STATISTICS
2836 #ifdef UQ_CODE_HAS_MONITORS
2837 template <
class V,
class M>
2839 SequenceOfVectors<V,M>::subMeanMonitorAlloc(
unsigned int numberOfMonitorPositions)
2841 m_subMeanMonitorPosSeq =
new ScalarSequence<double>(m_env, numberOfMonitorPositions,(m_name+
"_subMeanMonitorPosSeq").c_str());
2842 m_subMeanVecSeq =
new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+
"_subMeanVecSeq").c_str() );
2843 m_subMeanCltStdSeq =
new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+
"_subMeanCltStdSeq").c_str() );
2848 template <
class V,
class M>
2850 SequenceOfVectors<V,M>::subMeanInter0MonitorAlloc(
unsigned int numberOfMonitorPositions)
2852 m_subMeanInter0MonitorPosSeq =
new ScalarSequence<double>(m_env, numberOfMonitorPositions,(m_name+
"_subMeanInter0MonitorPosSeq").c_str() );
2853 m_subMeanInter0Mean =
new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+
"_subMeanInter0MeanSeq").c_str() );
2854 m_subMeanInter0Clt95 =
new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+
"_subMeanInter0Clt95Seq").c_str() );
2855 m_subMeanInter0Empirical90 =
new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+
"_subMeanInter0Empirical90Seq").c_str());
2856 m_subMeanInter0Min =
new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+
"_subMeanInter0MinSeq").c_str() );
2857 m_subMeanInter0Max =
new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+
"_subMeanInter0MaxSeq").c_str() );
2862 template <
class V,
class M>
2864 SequenceOfVectors<V,M>::unifiedMeanMonitorAlloc(
unsigned int numberOfMonitorPositions)
2866 m_unifiedMeanMonitorPosSeq =
new ScalarSequence<double>(m_env, numberOfMonitorPositions,(m_name+
"_unifiedMeanMonitorPosSeq").c_str());
2867 m_unifiedMeanVecSeq =
new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+
"_unifiedMeanVecSeq").c_str() );
2868 m_unifiedMeanCltStdSeq =
new SequenceOfVectors<V,M>(m_vectorSpace,numberOfMonitorPositions,(m_name+
"_unifiedMeanCltStdSeq").c_str() );
2873 template <
class V,
class M>
2875 SequenceOfVectors<V,M>::subMeanMonitorRun(
unsigned int monitorPosition,
2879 this->subMeanExtra(0,
2883 this->subMeanCltStd(0,
2891 template <
class V,
class M>
2893 SequenceOfVectors<V,M>::subMeanInter0MonitorRun(
unsigned int monitorPosition,
2894 V& subMeanInter0Mean,
2895 V& subMeanInter0Clt95,
2896 V& subMeanInter0Empirical90,
2897 V& subMeanInter0Min,
2898 V& subMeanInter0Max)
2900 V subMeanVec(m_vectorSpace.zeroVector());
2901 this->subMeanExtra(0,
2905 subMeanVec.mpiAllReduce(
RawValue_MPI_SUM,m_env.inter0Comm(),subMeanInter0Mean);
2906 subMeanInter0Mean /= ((double) m_env.inter0Comm().NumProc());
2908 V subMeanInter0CltVariance = subMeanVec-subMeanInter0Mean;
2909 subMeanInter0CltVariance *= subMeanInter0CltVariance;
2910 subMeanInter0CltVariance.mpiAllReduce(
RawValue_MPI_SUM,m_env.inter0Comm(),subMeanInter0Clt95);
2911 subMeanInter0Clt95 /= ((double) (m_env.inter0Comm().NumProc()-1));
2912 subMeanInter0Clt95 /= ((double) (m_env.inter0Comm().NumProc()-1));
2913 subMeanInter0Clt95.cwSqrt();
2914 subMeanInter0Clt95 *= 3.;
2916 V subMeanInter0Quantile5(m_vectorSpace.zeroVector());
2917 subMeanVec.mpiAllQuantile(.05,m_env.inter0Comm(),subMeanInter0Quantile5);
2918 V subMeanInter0Quantile95(m_vectorSpace.zeroVector());
2919 subMeanVec.mpiAllQuantile(.95,m_env.inter0Comm(),subMeanInter0Quantile95);
2920 subMeanInter0Empirical90 = subMeanInter0Quantile95 - subMeanInter0Quantile5;
2922 subMeanVec.mpiAllReduce(
RawValue_MPI_MIN,m_env.inter0Comm(),subMeanInter0Min);
2924 subMeanVec.mpiAllReduce(
RawValue_MPI_MAX,m_env.inter0Comm(),subMeanInter0Max);
2929 template <
class V,
class M>
2931 SequenceOfVectors<V,M>::unifiedMeanMonitorRun(
unsigned int monitorPosition,
2933 V& unifiedMeanCltStd)
2935 this->unifiedMeanExtra(0,
2939 this->unifiedMeanCltStd(0,
2946 template <
class V,
class M>
2948 SequenceOfVectors<V,M>::subMeanMonitorStore(
unsigned int i,
2949 unsigned int monitorPosition,
2950 const V& subMeanVec,
2951 const V& subMeanCltStd)
2953 (*m_subMeanMonitorPosSeq)[i] = monitorPosition;
2954 m_subMeanVecSeq->setPositionValues(i,subMeanVec);
2955 m_subMeanCltStdSeq->setPositionValues(i,subMeanCltStd);
2960 template <
class V,
class M>
2962 SequenceOfVectors<V,M>::subMeanInter0MonitorStore(
unsigned int i,
2963 unsigned int monitorPosition,
2964 const V& subMeanInter0Mean,
2965 const V& subMeanInter0Clt95,
2966 const V& subMeanInter0Empirical90,
2967 const V& subMeanInter0Min,
2968 const V& subMeanInter0Max)
2970 (*m_subMeanInter0MonitorPosSeq)[i] = monitorPosition;
2971 m_subMeanInter0Mean->setPositionValues(i,subMeanInter0Mean);
2972 m_subMeanInter0Clt95->setPositionValues(i,subMeanInter0Clt95);
2973 m_subMeanInter0Empirical90->setPositionValues(i,subMeanInter0Empirical90);
2974 m_subMeanInter0Min->setPositionValues(i,subMeanInter0Min);
2975 m_subMeanInter0Max->setPositionValues(i,subMeanInter0Max);
2980 template <
class V,
class M>
2982 SequenceOfVectors<V,M>::unifiedMeanMonitorStore(
unsigned int i,
2983 unsigned int monitorPosition,
2985 V& unifiedMeanCltStd)
2987 (*m_unifiedMeanMonitorPosSeq)[i] = monitorPosition;
2988 m_unifiedMeanVecSeq->setPositionValues(i,unifiedMeanVec);
2989 m_unifiedMeanCltStdSeq->setPositionValues(i,unifiedMeanCltStd);
2994 template <
class V,
class M>
2996 SequenceOfVectors<V,M>::subMeanMonitorWrite(std::ofstream& ofs)
3005 template <
class V,
class M>
3007 SequenceOfVectors<V,M>::subMeanInter0MonitorWrite(std::ofstream& ofs)
3019 template <
class V,
class M>
3021 SequenceOfVectors<V,M>::unifiedMeanMonitorWrite(std::ofstream& ofs)
3032 template <
class V,
class M>
3034 SequenceOfVectors<V,M>::subMeanMonitorFree()
3036 delete m_subMeanMonitorPosSeq;
3037 m_subMeanMonitorPosSeq = NULL;
3038 delete m_subMeanVecSeq;
3039 m_subMeanVecSeq = NULL;
3040 delete m_subMeanCltStdSeq;
3041 m_subMeanCltStdSeq = NULL;
3046 template <
class V,
class M>
3048 SequenceOfVectors<V,M>::subMeanInter0MonitorFree()
3050 delete m_subMeanInter0MonitorPosSeq;
3051 m_subMeanInter0MonitorPosSeq = NULL;
3052 delete m_subMeanInter0Mean;
3053 m_subMeanInter0Mean = NULL;
3054 delete m_subMeanInter0Clt95;
3055 m_subMeanInter0Clt95 = NULL;
3056 delete m_subMeanInter0Empirical90;
3057 m_subMeanInter0Empirical90 = NULL;
3058 delete m_subMeanInter0Min;
3059 m_subMeanInter0Min = NULL;
3060 delete m_subMeanInter0Max;
3061 m_subMeanInter0Max = NULL;
3066 template <
class V,
class M>
3068 SequenceOfVectors<V,M>::unifiedMeanMonitorFree()
3070 delete m_unifiedMeanMonitorPosSeq;
3071 m_unifiedMeanMonitorPosSeq = NULL;
3072 delete m_unifiedMeanVecSeq;
3073 m_unifiedMeanVecSeq = NULL;
3074 delete m_unifiedMeanCltStdSeq;
3075 m_unifiedMeanCltStdSeq = NULL;
3079 #endif // #ifdef UQ_CODE_HAS_MONITORS
T subInterQuantileRange(unsigned int initialPos) const
Returns the interquartile range of the values in the sub-sequence.
double MiscGetEllapsedSeconds(struct timeval *timeval0)
void subUniformlySampledCdf(const V &numEvaluationPointsVec, ArrayOfOneDGrids< V, M > &cdfGrids, ArrayOfOneDTables< V, M > &cdfValues) const
Uniformly samples from the CDF from the sub-sequence.
void deleteStoredVectors()
Deletes all the stored vectors.
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.
void unifiedReadContents(const std::string &fileName, const std::string &fileType, const unsigned int subSequenceSize)
Reads the unified sequence from a file.
Class to accommodate arrays of one-dimensional grid.
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...
void subHistogram(unsigned int initialPos, const V &minVec, const V &maxVec, std::vector< V * > ¢ersForAllBins, std::vector< V * > &quanttsForAllBins) const
Calculates the histogram of the sub-sequence.
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...
Struct for handling data input and output from files.
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).
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.
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...
void unifiedInterQuantileRange(unsigned int initialPos, V &unifiedIqrVec) const
Returns the interquartile range of the values in the unified sequence.
GslMatrix matrixProduct(const GslVector &v1, const GslVector &v2)
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...
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...
void extractScalarSeq(unsigned int initialPos, unsigned int spacing, unsigned int numPos, unsigned int paramId, ScalarSequence< double > &scalarSeq) const
Extracts a sequence of scalars.
#define UQ_FATAL_RC_MACRO(macroIRc, givenRank, where, what)
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...
unsigned int subSequenceSize() const
Size of the sub-sequence of vectors.
SequenceOfVectors< V, M > & operator=(const SequenceOfVectors< V, M > &rhs)
Copies values from rhs to this.
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...
std::ofstream * ofsVar
Provides a stream interface to write data to files.
void setUniformGrids(const V &sizesVec, const V &minPositionsVec, const V &maxPositionsVec)
Sets an array of uniform grids.
#define UQ_FILE_EXTENSION_FOR_MATLAB_FORMAT
double estimateConvBrooksGelman(unsigned int initialPos, unsigned int numPos) const
Estimates convergence rate using Brooks & Gelman method.
std::vector< const V * >::iterator seqVectorPositionIteratorTypedef
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...
void resizeSequence(unsigned int newSequenceSize)
Resizes the size of the sequence of scalars.
#define UQ_FATAL_TEST_MACRO(test, givenRank, where, what)
void unifiedWriteContents(const std::string &fileName, const std::string &fileType) const
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...
void autoCovariance(unsigned int initialPos, unsigned int numPos, const V &meanVec, unsigned int lag, V &covVec) const
Calculates the autocovariance.
void erasePositions(unsigned int initialPos, unsigned int numPos)
Erases numPos elements of the sequence starting at position initialPos.
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...
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.
void copy(const BaseVectorSequence< V, M > &src)
Copies vector sequence src to this.
void extractRawData(unsigned int initialPos, unsigned int spacing, unsigned int numPos, unsigned int paramId, std::vector< double > &rawData) const
Extracts the raw data.
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...
std::ifstream * ifsVar
Provides a stream interface to read data from files.
void setOneDTable(unsigned int rowId, const std::vector< double > &values)
Sets the one-dimensional table.
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...
T unifiedInterQuantileRange(bool useOnlyInter0Comm, unsigned int initialPos) const
Returns the interquartile range of the values in the unified sequence.
void unifiedUniformlySampledCdf(const V &numEvaluationPointsVec, ArrayOfOneDGrids< V, M > &unifiedCdfGrids, ArrayOfOneDTables< V, M > &unifiedCdfValues) const
Uniformly samples from the CDF from the sub-sequence.
void getPositionValues(unsigned int posId, V &vec) const
Gets the values of the sequence at position posId and stores them at vec.
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...
void autoCorrViaDef(unsigned int initialPos, unsigned int numPos, unsigned int lag, V &corrVec) const
Calculates the autocorrelation via definition.
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-...
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.
void resizeSequence(unsigned int newSubSequenceSize)
Resizes the sequence.
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.
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.
void filter(unsigned int initialPos, unsigned int spacing)
Filters positions in the sequence of vectors.
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.
Class to accommodate arrays of one-dimensional tables.
#define UQ_FILE_EXTENSION_FOR_HDF_FORMAT
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 ...
Base class for handling vector and array samples (sequence of vectors or arrays). ...
void subHistogram(unsigned int initialPos, const T &minHorizontalValue, const T &maxHorizontalValue, std::vector< T > ¢ers, std::vector< unsigned int > &bins) const
Calculates the histogram of the sub-sequence.
void subInterQuantileRange(unsigned int initialPos, V &iqrVec) const
Returns the interquartile range of the values in the sub-sequence.
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...
void copy(const SequenceOfVectors< V, M > &src)
Copies vector sequence src to this.
A class representing a vector space.
void select(const std::vector< unsigned int > &idsOfUniquePositions)
TODO: It shall select positions in the sequence of vectors.
Class for handling vector samples (sequence of vectors).
~SequenceOfVectors()
Destructor.
void setPositionValues(unsigned int posId, const V &vec)
Set the values in vec at position posId of the sequence.
void resetValues(unsigned int initialPos, unsigned int numPos)
Resets a total of numPos values of the sequence starting at position initialPos.
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...
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 ...
SequenceOfVectors(const VectorSpace< V, M > &vectorSpace, unsigned int subSequenceSize, const std::string &name)
Default constructor.
void autoCorrViaFft(unsigned int initialPos, unsigned int numPos, unsigned int maxLag, std::vector< T > &autoCorrs) const
Calculates the autocorrelation via Fast Fourier transforms (FFT).
std::vector< const V * > m_seq
Sequence of vectors.