queso-0.51.1
TeuchosVector.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,2009,2010,2011,2012,2013 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/Defines.h>
26 
27 #ifdef QUESO_HAS_TRILINOS
28 
29 #include <queso/TeuchosVector.h>
30 
31 namespace QUESO {
32 
33 using std:: cout;
34 using std:: endl;
35 
36 // ---------------------------------------------------
37 // default constructor -------------------------------
38 TeuchosVector::TeuchosVector() :
39  Vector()
40 {
41  UQ_FATAL_TEST_MACRO(true,
42  m_env.worldRank(),
43  "TeuchosVector::constructor(), default",
44  "should not be used by user");
45 };
46 
47 
48 // constructor with dimension -----------------------
49 TeuchosVector::TeuchosVector(const BaseEnvironment& env, const Map& map)
50  :
51  Vector(env,map)
52 {
53  m_vec.size(map.NumGlobalElements());
54 
55  UQ_FATAL_TEST_MACRO((m_vec == NULL),
56  m_env.worldRank(),
57  "TeuchosVector::constructor(1)",
58  "null vector generated");
59 
60  UQ_FATAL_TEST_MACRO(m_vec.length() != map.NumMyElements(),
61  m_env.worldRank(),
62  "TeuchosVector::constructor(1)",
63  "incompatible local vec size");
64 
65  UQ_FATAL_TEST_MACRO(m_vec.length() != map.NumGlobalElements(),
66  m_env.worldRank(),
67  "TeuchosVector::constructor(1)",
68  "incompatible global vec size");
69 
70  UQ_FATAL_TEST_MACRO(m_vec.length() != m_map.NumMyElements(),
71  m_env.worldRank(),
72  "TeuchosVector::constructor(1)",
73  "incompatible own vec size");
74 
75  //std::cout << "In TeuchosVector::constructor(env,map)"
76  // << "\n m_vec.length() = " << m_vec.length()
77  // << "\n map.NumGlobalElements() = " << map.NumGlobalElements()
78  // << "\n map.NumMyElements() = " << map.NumMyElements()
79  // << std::endl;
80 }
81 
82 // ---------------------------------------------------
83 TeuchosVector::TeuchosVector(const BaseEnvironment& env, const Map& map, double value)
84  :
85  Vector(env,map)
86 {
87  m_vec.size(map.NumGlobalElements());
88  m_vec = value;
89 
90  UQ_FATAL_TEST_MACRO((m_vec == NULL),
91  m_env.worldRank(),
92  "TeuchosVector::constructor(2)",
93  "null vector generated");
94 
95  UQ_FATAL_TEST_MACRO(m_vec.length() != map.NumMyElements(),
96  m_env.worldRank(),
97  "TeuchosVector::constructor(2)",
98  "incompatible local vec size");
99 
100  UQ_FATAL_TEST_MACRO(m_vec.length() != map.NumGlobalElements(),
101  m_env.worldRank(),
102  "TeuchosVector::constructor(2)",
103  "incompatible global vec size");
104 
105  UQ_FATAL_TEST_MACRO(m_vec.length() != m_map.NumMyElements(),
106  m_env.worldRank(),
107  "TeuchosVector::constructor(2)",
108  "incompatible own vec size");
109 }
110 
111 // ---------------------------------------------------
112 TeuchosVector::TeuchosVector(const BaseEnvironment& env, double d1, double d2, const Map& map)
113  :
114  Vector(env,map)
115  {
116  m_vec.size(map.NumGlobalElements());
117 
118  UQ_FATAL_TEST_MACRO((m_vec == NULL),
119  m_env.worldRank(),
120  "TeuchosVector::constructor(3), linspace",
121  "null vector generated");
122 
123  UQ_FATAL_TEST_MACRO(m_vec.length() != map.NumMyElements(),
124  m_env.worldRank(),
125  "TeuchosVector::constructor(3)",
126  "incompatible local vec size");
127 
128  UQ_FATAL_TEST_MACRO(m_vec.length() != map.NumGlobalElements(),
129  m_env.worldRank(),
130  "TeuchosVector::constructor(3)",
131  "incompatible global vec size");
132 
133  for (int i = 0; i < m_vec.length(); ++i) {
134  double alpha = (double) i / ((double) m_vec.length() - 1.);
135  (*this)[i] = (1.-alpha)*d1 + alpha*d2;
136  }
137 
138  UQ_FATAL_TEST_MACRO(m_vec.length() != m_map.NumMyElements(),
139  m_env.worldRank(),
140  "TeuchosVector::constructor(3)",
141  "incompatible own vec size");
142 }
143 
144 // ---------------------------------------------------
145 TeuchosVector::TeuchosVector(const TeuchosVector& v, double d1, double d2)
146  :
147  Vector(v.env(),v.map())
148 {
149  m_vec.size(v.sizeLocal());
150 
151  UQ_FATAL_TEST_MACRO((m_vec == NULL),
152  m_env.worldRank(),
153  "TeuchosVector::constructor(4), linspace",
154  "null vector generated");
155 
156  UQ_FATAL_TEST_MACRO(m_vec.length() != v.map().NumMyElements(),
157  m_env.worldRank(),
158  "TeuchosVector::constructor(4)",
159  "incompatible local vec size");
160 
161  UQ_FATAL_TEST_MACRO(m_vec.length() != v.map().NumGlobalElements(),
162  m_env.worldRank(),
163  "TeuchosVector::constructor(4)",
164  "incompatible global vec size");
165 
166  for ( int i = 0; i < m_vec.length(); ++i) {
167  double alpha = (double) i / ((double) m_vec.length() - 1.);
168  (*this)[i] = (1.-alpha)*d1 + alpha*d2;
169  }
170 
171  UQ_FATAL_TEST_MACRO(m_vec.length() != m_map.NumMyElements(),
172  m_env.worldRank(),
173  "TeuchosVector::constructor(4)",
174  "incompatible own vec size");
175 }
176 
177 // ---------------------------------------------------
178 TeuchosVector::TeuchosVector(const TeuchosVector& v) // mox
179  :
180  Vector(v.env(),v.map())
181  {
182  m_vec.size(v.sizeLocal());
183 
184  UQ_FATAL_TEST_MACRO((m_vec == NULL),
185  m_env.worldRank(),
186  "TeuchosVector::constructor(5), copy",
187  "null vector generated");
188 
189  UQ_FATAL_TEST_MACRO(m_vec.length() != v.map().NumMyElements(),
190  m_env.worldRank(),
191  "TeuchosVector::constructor(5)",
192  "incompatible local vec size");
193 
194  UQ_FATAL_TEST_MACRO(m_vec.length() != v.map().NumGlobalElements(),
195  m_env.worldRank(),
196  "TeuchosVector::constructor(5)",
197  "incompatible global vec size");
198  this->copy(v);
199 
200  UQ_FATAL_TEST_MACRO(m_vec.length() != m_map.NumMyElements(),
201  m_env.worldRank(),
202  "TeuchosVector::constructor(5)",
203  "incompatible own vec size");
204 }
205 
206 
207 // destructor --------------------------------------
208 TeuchosVector::~TeuchosVector()
209 {
210 };
211 
212 // Set methods ------------------------------------
213 //-------------------------------------------------
214 // assigns values a to all entrances in the vector
215 TeuchosVector& TeuchosVector::operator=(double a)
216 {
217  m_vec.putScalar(a);
218  return *this;
219 }
220 
221 //-------------------------------------------------
222 TeuchosVector& TeuchosVector::operator=(const TeuchosVector& rhs)
223 {
224  unsigned int size1 = m_vec.length();
225  unsigned int size2 = rhs.sizeLocal();
226 
227  UQ_FATAL_TEST_MACRO(size1 != size2, // mox
228  m_env.worldRank(),
229  "TeuchosVector::operator=()",
230  "the vectors do NOT have the same size.\n");
231 
232  if (size1==size2){
233  for (unsigned int i=0;i<size1;i++){
234  m_vec[i]=rhs[i];
235  }
236  }
237 
238  return *this;
239 }
240 
241 //-------------------------------------------------
242 TeuchosVector& TeuchosVector::operator*=(double a)
243 {
244  m_vec.scale(a); //Scale this vector by a; *this = a*this.
245  return *this;
246 }
247 
248  //-------------------------------------------------
249 TeuchosVector& TeuchosVector::operator/=(double a)
250 {
251  (*this) *= (1.0/a);
252  return *this;
253 }
254 
255 //-------------------------------------------------
256 TeuchosVector& TeuchosVector::operator*=(const TeuchosVector& rhs)
257 {
258  unsigned int size1 = this->sizeLocal();
259  unsigned int size2 = rhs.sizeLocal();
260  UQ_FATAL_TEST_MACRO((size1 != size2),
261  m_env.worldRank(),
262  "TeuchosVector::operator*=()",
263  "the vectors do NOT have the same size.\n");
264  if (size1==size2){
265  for (unsigned int i = 0; i < size1; ++i) {
266  (*this)[i] *= rhs[i];
267  }
268  }
269  return *this;
270 }
271 
272 //-------------------------------------------------
273 TeuchosVector& TeuchosVector::operator/=(const TeuchosVector& rhs)
274 {
275  unsigned int size1 = this->sizeLocal();
276  unsigned int size2 = rhs.sizeLocal();
277 
278  UQ_FATAL_TEST_MACRO((size1 != size2),
279  m_env.worldRank(),
280  "TeuchosVector::operator/=()",
281  "the vectors do NOT have the same size.\n");
282  if (size1==size2){
283  for (unsigned int i = 0; i < size1; ++i) {
284  (*this)[i] /= rhs[i];
285  }
286  }
287  return *this;
288 }
289 
290 //-------------------------------------------------
291 TeuchosVector& TeuchosVector::operator+=(const TeuchosVector& rhs)
292 {
293  unsigned int size1 = this->sizeLocal();
294  unsigned int size2 = rhs.sizeLocal();
295 
296  UQ_FATAL_TEST_MACRO((size1 != size2),
297  m_env.worldRank(),
298  "TeuchosVector::operator+=()",
299  "the vectors do NOT have the same size.\n");
300 
301  if (size1==size2){
302  for (unsigned int i = 0; i < size1; ++i) {
303  (*this)[i] += rhs[i];
304  }
305  }
306  return *this;
307 }
308 
309 //-------------------------------------------------
310 TeuchosVector& TeuchosVector::operator-=(const TeuchosVector& rhs)
311 {
312  unsigned int size1 = this->sizeLocal();
313  unsigned int size2 = rhs.sizeLocal();
314 
315  UQ_FATAL_TEST_MACRO((size1 != size2),
316  m_env.worldRank(),
317  "TeuchosVector::operator-=()",
318  "the vectors do NOT have the same size.\n");
319 
320  if (size1==size2){
321  for (unsigned int i = 0; i < size1; ++i) {
322  (*this)[i] -= rhs[i];
323  }
324  }
325 
326  return *this;
327 }
328 
329 
330 // Accessor methods --------------------------------
331 //-------------------------------------------------
332 double& TeuchosVector::operator[](unsigned int i)
333 {
334  return m_vec[i];
335 }
336 
337 //-------------------------------------------------
338 const double& TeuchosVector::operator[](unsigned int i) const
339 {
340  return m_vec[i];
341 }
342 
343 // Attribute methods ------------------------------
344 //-------------------------------------------------
345 unsigned int TeuchosVector::sizeLocal() const
346 {
347  UQ_FATAL_TEST_MACRO(m_vec.length() != (int) m_map.NumMyElements(),
348  m_env.worldRank(),
349  "TeuchosVector::sizeLocal()",
350  "incompatible vec size");
351 
352  return m_vec.length();
353 }
354 
355 //-------------------------------------------------
356 unsigned int TeuchosVector::sizeGlobal() const
357 {
358  UQ_FATAL_TEST_MACRO(m_vec.length() != (int) m_map.NumGlobalElements(),
359  m_env.worldRank(),
360  "TeuchosVector::sizeGlobal()",
361  "incompatible vec size");
362  return m_vec.length();
363 }
364 
365 //-------------------------------------------------
366 // TODO: needs to be checked. It may not be used at all. Kemelli 4/30/13.
367 double*
368 TeuchosVector::values()
369 {
370  return m_vec.values();
371 };
372 
373 // Getting Max and Min values -------------------
374 // Max ------------------------------------------
375 double TeuchosVector::getMaxValue( ) const //dummy
376 {
377  const unsigned int size = this->sizeLocal();
378  std::vector<double> aux;
379 
380  for (unsigned int i=0; i<size; i++ ) {
381  aux.push_back((*this)[i]) ;
382  }
383 
384  return *max_element (aux.begin(),aux.end());
385 }
386 
387 // Min ------------------------------------------
388 double TeuchosVector::getMinValue( ) const //dummy
389 {
390  const unsigned int size = this->sizeLocal();
391  std::vector<double> aux;
392 
393  for (unsigned int i=0; i<size; i++ ) {
394  aux.push_back((*this)[i]) ;
395  }
396 
397  return *min_element (aux.begin(),aux.end());
398 
399 }
400 
401 // Max index -----------------------------------
402 int TeuchosVector::getMaxValueIndex( ) const //dummy
403 {
404  const unsigned int size = this->sizeLocal();
405  std::vector<double> vect;
406 
407  for (unsigned int i=0; i<size; i++ ) {
408  vect.push_back((*this)[i]) ;
409  }
410  std::vector<double>::iterator iter_max = max_element(vect.begin(), vect.end());
411  return distance(vect.begin(), iter_max);
412 }
413 
414 // Min index -----------------------------------
415 int TeuchosVector::getMinValueIndex( ) const //dummy
416 {
417  const unsigned int size = this->sizeLocal();
418  std::vector<double> vect;
419 
420  for (unsigned int i=0; i<size; i++ ) {
421  vect.push_back((*this)[i]) ;
422  }
423  std::vector<double>::iterator iter_min = min_element(vect.begin(), vect.end());
424  return distance(vect.begin(), iter_min);
425 }
426 
427 // Max and index -------------------------------
428 void TeuchosVector::getMaxValueAndIndex( double& max_value, int& max_value_index ) //dummy
429 {
430  const unsigned int size = this->sizeLocal();
431  std::vector<double> vect;
432 
433  for (unsigned int i=0; i<size; i++ ) {
434  vect.push_back((*this)[i]) ;
435  }
436  std::vector<double>::iterator iter_max = max_element(vect.begin(), vect.end());
437 
438  max_value = *iter_max;
439  max_value_index = distance(vect.begin(), iter_max);
440 
441  return;
442 }
443 
444 // Min and index -------------------------------
445 void TeuchosVector::getMinValueAndIndex( double& min_value, int& min_value_index ) //dummy
446 {
447  const unsigned int size = this->sizeLocal();
448  std::vector<double> vect;
449 
450  for (unsigned int i=0; i<size; i++ ) {
451  vect.push_back((*this)[i]) ;
452  }
453  std::vector<double>::iterator iter_min = min_element(vect.begin(), vect.end());
454 
455  min_value = *iter_min;
456  min_value_index = distance(vect.begin(), iter_min);
457 
458  return;
459 }
460 
461 // Norm methods ------------------------------------
462 // -------------------------------------------------
463 double TeuchosVector::norm2Sq() const
464 {
465  return (m_vec).dot(m_vec );
466 }
467 
468 //-------------------------------------------------
469 double TeuchosVector::norm2() const
470 {
471  return std::sqrt(this->norm2Sq());
472 }
473 
474 //-------------------------------------------------
475 double TeuchosVector::norm1() const
476 {
477  double result = 0.;
478 
479  unsigned int size = this->sizeLocal();
480  for (unsigned int i = 0; i < size; ++i) {
481  result += fabs((*this)[i]);
482  }
483 
484  return result;
485 }
486 
487 //-------------------------------------------------
488 double TeuchosVector::normInf() const
489 {
490  double result = 0.;
491 
492  unsigned int size = this->sizeLocal();
493  double aux = 0.;
494  for (unsigned int i = 0; i < size; ++i) {
495  aux = fabs((*this)[i]);
496  if (aux > result) result = aux;
497  }
498 
499  return result;
500 }
501 
502 // Comparison methods -----------------------------
503 //-------------------------------------------------
504 bool
505 TeuchosVector::atLeastOneComponentSmallerThan(const TeuchosVector& rhs) const
506 {
507  UQ_FATAL_TEST_MACRO((this->sizeLocal() != rhs.sizeLocal()),
508  m_env.worldRank(),
509  "TeuchosVector::atLeastOneComponentSmallerThan()",
510  "vectors have different sizes");
511 
512  bool result = false;
513  unsigned int i = 0;
514  unsigned int size = this->sizeLocal();
515  while ((i < size) && (result == false)) {
516  result = ( (*this)[i] < rhs[i] );
517  i++;
518  };
519 
520  return result;
521 }
522 
523 //-------------------------------------------------
524 bool
525 TeuchosVector::atLeastOneComponentBiggerThan(const TeuchosVector& rhs) const
526 {
527  UQ_FATAL_TEST_MACRO((this->sizeLocal() != rhs.sizeLocal()),
528  m_env.worldRank(),
529  "TeuchosVector::atLeastOneComponentBiggerThan()",
530  "vectors have different sizes");
531 
532  bool result = false;
533  unsigned int i = 0;
534  unsigned int size = this->sizeLocal();
535  while ((i < size) && (result == false)) {
536  result = ( (*this)[i] > rhs[i] );
537  i++;
538  };
539 
540  return result;
541 }
542 
543 //-------------------------------------------------
544 bool
545 TeuchosVector::atLeastOneComponentSmallerOrEqualThan(const TeuchosVector& rhs) const
546 {
547  UQ_FATAL_TEST_MACRO((this->sizeLocal() != rhs.sizeLocal()),
548  m_env.worldRank(),
549  "TeuchosVector::atLeastOneComponentSmallerOrEqualThan()",
550  "vectors have different sizes");
551 
552  bool result = false;
553  unsigned int i = 0;
554  unsigned int size = this->sizeLocal();
555  while ((i < size) && (result == false)) {
556  result = ( (*this)[i] <= rhs[i] ); // prudencio 2012-02-06
557  i++;
558  };
559 
560  return result;
561 }
562 
563 //-------------------------------------------------
564 bool
565 TeuchosVector::atLeastOneComponentBiggerOrEqualThan(const TeuchosVector& rhs) const
566 {
567  UQ_FATAL_TEST_MACRO((this->sizeLocal() != rhs.sizeLocal()),
568  m_env.worldRank(),
569  "TeuchosVector::atLeastOneComponentBiggerOrEqualThan()",
570  "vectors have different sizes");
571 
572  bool result = false;
573  unsigned int i = 0;
574  unsigned int size = this->sizeLocal();
575  while ((i < size) && (result == false)) {
576  result = ( (*this)[i] >= rhs[i] ); // prudencio 2012-02-06
577  i++;
578  };
579 
580  return result;
581 }
582 
583 
584 // Get/Set methods -----------------------------------
585 //----------------------------------------------------
586 void TeuchosVector::cwSet(double value)
587 {
588  (*this)=value;
589 
590  return;
591 }
592 
593 //----------------------------------------------------
594 void TeuchosVector::cwSet(unsigned int initialPos, const TeuchosVector& vec)
595 {
596  UQ_FATAL_TEST_MACRO(initialPos >= this->sizeLocal(),
597  m_env.worldRank(),
598  "TeuchosVector::cwSet()",
599  "invalid initialPos");
600 
601  UQ_FATAL_TEST_MACRO((initialPos +vec.sizeLocal()) > this->sizeLocal(),
602  m_env.worldRank(),
603  "TeuchosVector::cwSet()",
604  "invalid vec.sizeLocal()");
605 
606  for (unsigned int i = 0; i < vec.sizeLocal(); ++i) {
607  (*this)[initialPos+i] = vec[i];
608  }
609 
610  return;
611 }
612 
613 
614 //----------------------------------------------------
615 /* extracts elements from vector (*this) starting at position initialPos and save in vec */
616 void TeuchosVector::cwExtract(unsigned int initialPos, TeuchosVector& vec) const
617 {
618  UQ_FATAL_TEST_MACRO(initialPos >= this->sizeLocal(),
619  m_env.worldRank(),
620  "TeuchosVector::cwExtract()",
621  "invalid initialPos");
622 
623  UQ_FATAL_TEST_MACRO((initialPos +vec.sizeLocal()) > this->sizeLocal(),
624  m_env.worldRank(),
625  "TeuchosVector::cwExtract()",
626  "invalid vec.sizeLocal()");
627 
628  for (unsigned int i = 0; i < vec.sizeLocal(); ++i) {
629  vec[i] = (*this)[initialPos+i];
630  }
631 
632  return;
633 }
634 
635 //----------------------------------------------------
636 void TeuchosVector::cwInvert()
637 {
638  unsigned int size = this->sizeLocal();
639  for (unsigned int i = 0; i < size; ++i) {
640  (*this)[i] = 1./(*this)[i];
641  }
642 
643  return;
644 }
645 
646 //----------------------------------------------------
647 void TeuchosVector::cwSqrt()
648 {
649  unsigned int size = this->sizeLocal();
650  for (unsigned int i = 0; i < size; ++i) {
651  (*this)[i] = sqrt((*this)[i]);
652  }
653 
654  return;
655 }
656 
657 //----------------------------------------------------
658 void
659 TeuchosVector::cwSetConcatenated(const TeuchosVector& v1, const TeuchosVector& v2)
660 {
661  UQ_FATAL_TEST_MACRO(this->sizeLocal() != v1.sizeLocal() + v2.sizeLocal(),
662  m_env.worldRank(),
663  "TeuchosVector::cwSetConcatenated(1)",
664  "incompatible vector sizes");
665 
666 // if ( this->sizeLocal() != v1.sizeLocal() + v2.sizeLocal() ) {
667 // std::cout << "ERROR in TeuchosVector:: cwSetConcatenated ---> the vectors' sizes are not compatible.\n";
668 // cout << " in TeuchosVector:: cwSetConcatenated ---> resizing resulting vector... new size = "
669 // << v1.sizeLocal()+v1.sizeLocal() <<endl;
670 //
671 // m_vec.resize(v1.sizeLocal()+v1.sizeLocal());
672 // }
673 
674  for (unsigned int i = 0; i < v1.sizeLocal(); ++i) {
675  (*this)[i] = v1[i];
676  }
677 
678  for (unsigned int i = 0; i < v2.sizeLocal(); ++i) {
679  (*this)[v1.sizeLocal()+i] = v2[i];
680  }
681 
682  return;
683 }
684 
685 // -------------------------------------------------
686 //updated on 3/18, to use the RngBase+Boost
687 void TeuchosVector::cwSetGaussian(double mean, double stdDev)
688 {
689  for (unsigned int i = 0; i < this->sizeLocal(); ++i) {
690  (*this)[i] = mean + m_env.rngObject()->gaussianSample(stdDev);
691  }
692  return;
693 };
694 
695 // -------------------------------------------------
696 //updated on 3/18, to use the RngBase+Boost
697 void TeuchosVector::cwSetGaussian(const TeuchosVector& meanVec, const TeuchosVector& stdDevVec)
698 {
699  for (unsigned int i = 0; i < this->sizeLocal(); ++i) {
700  (*this)[i] = meanVec[i] + m_env.rngObject()->gaussianSample(stdDevVec[i]);
701  }
702  return;
703 };
704 
705 
706 //----------------------------------------------------
707 //implemented/checked 2/26/13
708 //updated on 3/18, to use the RngBase+Boost
709  void TeuchosVector::cwSetUniform(const TeuchosVector& aVec, const TeuchosVector& bVec)
710 {
711  for (unsigned int i = 0; i < this->sizeLocal(); ++i) {
712  (*this)[i] = aVec[i] + (bVec[i]-aVec[i])*m_env.rngObject()->uniformSample();
713  }
714  return;
715 }
716 
717 
718 // -------------------------------------------------
719 //updated on 3/18, to use the RngBase+Boost
720 void TeuchosVector::cwSetBeta(const TeuchosVector& alpha, const TeuchosVector& beta)
721 {
722  UQ_FATAL_TEST_MACRO(this->sizeLocal() != alpha.sizeLocal(),
723  m_env.worldRank(),
724  "TeuchosVector::cwSetBeta()",
725  "incompatible alpha size");
726 
727  UQ_FATAL_TEST_MACRO(this->sizeLocal() != beta.sizeLocal(),
728  m_env.worldRank(),
729  "TeuchosVector::cwSetBeta()",
730  "incompatible beta size");
731 
732  for (unsigned int i = 0; i < this->sizeLocal(); ++i)
733  {
734  (*this)[i] = m_env.rngObject()->betaSample(alpha[i],beta[i]);
735 
736  if ((m_env.subDisplayFile()) && (m_env.displayVerbosity() >= 99))
737  {
738  *m_env.subDisplayFile() << "In TeuchosVector::cwSetBeta()"
739  << ": fullRank " << m_env.fullRank()
740  << ", i = " << i
741  << ", alpha[i] = " << alpha[i]
742  << ", beta[i] = " << beta[i]
743  << ", sample = " << (*this)[i]
744  << std::endl;
745  }
746  }
747  return;
748 };
749 
750 // -------------------------------------------------
751 //updated on 3/18, to use the RngBase+Boost
752 void TeuchosVector::cwSetGamma(const TeuchosVector& aVec, const TeuchosVector& bVec)
753 {
754  UQ_FATAL_TEST_MACRO(this->sizeLocal() != aVec.sizeLocal(),
755  m_env.worldRank(),
756  "TeuchosVector::cwSetGamma()",
757  "incompatible a size");
758 
759  UQ_FATAL_TEST_MACRO(this->sizeLocal() != bVec.sizeLocal(),
760  m_env.worldRank(),
761  "TeuchosVector::cwSetGamma()",
762  "incompatible b size");
763 
764  for (unsigned int i = 0; i < this->sizeLocal(); ++i) {
765  (*this)[i] = m_env.rngObject()->gammaSample(aVec[i],bVec[i]);
766  }
767  return;
768 }
769 
770 // -------------------------------------------------
771 //updated on 3/18, to use the RngBase+Boost
772 // Using Gamma Distribution to calculate InverseGamma.
773 // Note the divisions: 1.0/b and the 1.0/generator; they are crucial
774 void TeuchosVector::cwSetInverseGamma(const TeuchosVector& alpha, const TeuchosVector& beta)
775 {
776  UQ_FATAL_TEST_MACRO(this->sizeLocal() != alpha.sizeLocal(),
777  m_env.worldRank(),
778  "TeuchosVector::cwSetInverseGamma()",
779  "incompatible alpha size");
780 
781  UQ_FATAL_TEST_MACRO(this->sizeLocal() != beta.sizeLocal(),
782  m_env.worldRank(),
783  "TeuchosVector::cwSetInverseGamma()",
784  "incompatible beta size");
785 
786  for (unsigned int i = 0; i < this->sizeLocal(); ++i) {
787  (*this)[i] = 1./m_env.rngObject()->gammaSample(alpha[i],1./beta[i]);
788  }
789  return;
790 }
791 
792 
793 // Miscellaneous methods -------------------------
794 // absolute values -------------------------------
795 // tested 1/8/13
796 TeuchosVector
797 TeuchosVector::abs() const
798 {
799  TeuchosVector abs_of_this_vec( *this );
800 
801  unsigned int size = abs_of_this_vec.sizeLocal();
802 
803  for( unsigned int i = 0; i < size; ++i )
804  {
805  abs_of_this_vec[i] = std::fabs( (*this)[i] );
806  }
807 
808  return abs_of_this_vec;
809 
810 }
811 
812 // -------------------------------------------------
813 void
814 TeuchosVector::copy_to_std_vector(std::vector<double>& vec)
815 {
816  unsigned int size = this->sizeLocal();
817  vec.resize(size);
818 
819  for (unsigned int i = 0; i < size; ++i)
820  vec[i] = m_vec[i];
821 
822  return;
823 }
824 
825 // -------------------------------------------------
826 void
827 TeuchosVector::copy_from_std_vector(const std::vector<double> vec)
828 {
829  unsigned int size1 = vec.size(), size2= this->sizeLocal();
830 
831  UQ_FATAL_TEST_MACRO((size1 != size2),
832  m_env.worldRank(),
833  "In TeuchosVector::copy_from_std_vector()",
834  "vectors have different sizes");
835 
836  for (unsigned int i = 0; i < size1; ++i)
837  m_vec[i] = vec[i];
838 
839  return;
840 }
841 
842 // -------------------------------------------------
843 void TeuchosVector::sort()
844 {
845  std::vector<double> vec;
846 
847  (*this).copy_to_std_vector(vec);
848 
849  // using default comparison (operator <):
850  std::sort (vec.begin(), vec.end());
851 
852  (*this).copy_from_std_vector(vec);
853 };
854 
855 // -------------------------------------------------
856 double TeuchosVector::sumOfComponents() const
857 {
858  double result = 0.;
859  unsigned int size = this->sizeLocal();
860  for (unsigned int i = 0; i < size; ++i) {
861  result += (*this)[i];
862  }
863 
864  return result;
865 }
866 
867 // -------------------------------------------------
868 // added 2/28/13
869 void
870 TeuchosVector::mpiBcast(int srcRank, const MpiComm& bcastComm)
871 {
872  // Filter out those nodes that should not participate
873  if (bcastComm.MyPID() < 0) return;
874 
875  // Check 'srcRank'
876  UQ_FATAL_TEST_MACRO((srcRank < 0) || (srcRank >= bcastComm.NumProc()),
877  m_env.worldRank(),
878  "TeuchosVector::mpiBcast()",
879  "invalud srcRank");
880 
881  // Check number of participant nodes
882  double localNumNodes = 1.;
883  double totalNumNodes = 0.;
884  bcastComm.Allreduce((void *) &localNumNodes, (void *) &totalNumNodes, (int) 1, RawValue_MPI_DOUBLE, RawValue_MPI_SUM,
885  "TeuchosVector::mpiBcast()",
886  "failed MPI.Allreduce() for numNodes");
887  UQ_FATAL_TEST_MACRO(((int) totalNumNodes) != bcastComm.NumProc(),
888  m_env.worldRank(),
889  "TeuchosVector::mpiBcast()",
890  "inconsistent numNodes");
891 
892  // Check that all participant nodes have the same vector size
893  double localVectorSize = this->sizeLocal();
894  double sumOfVectorSizes = 0.;
895  bcastComm.Allreduce((void *) &localVectorSize, (void *) &sumOfVectorSizes, (int) 1, RawValue_MPI_DOUBLE, RawValue_MPI_SUM,
896  "TeuchosVector::mpiBcast()",
897  "failed MPI.Allreduce() for vectorSize");
898 
899  if ( ((unsigned int) sumOfVectorSizes) != ((unsigned int)(totalNumNodes*localVectorSize)) ) {
900  std::cerr << "rank " << bcastComm.MyPID()
901  << ": sumOfVectorSizes = " << sumOfVectorSizes
902  << ", totalNumNodes = " << totalNumNodes
903  << ", localVectorSize = " << localVectorSize
904  << std::endl;
905  }
906  bcastComm.Barrier();
907  UQ_FATAL_TEST_MACRO(((unsigned int) sumOfVectorSizes) != ((unsigned int)(totalNumNodes*localVectorSize)),
908  m_env.worldRank(),
909  "TeuchosVector::mpiBcast()",
910  "inconsistent vectorSize");
911 
912  // Ok, bcast data
913  std::vector<double> dataBuffer((unsigned int) localVectorSize, 0.);
914  if (bcastComm.MyPID() == srcRank) {
915  for (unsigned int i = 0; i < dataBuffer.size(); ++i) {
916  dataBuffer[i] = (*this)[i];
917  }
918  }
919 
920  bcastComm.Bcast((void *) &dataBuffer[0], (int) localVectorSize, RawValue_MPI_DOUBLE, srcRank,
921  "TeuchosVector::mpiBcast()",
922  "failed MPI.Bcast()");
923 
924  if (bcastComm.MyPID() != srcRank) {
925  for (unsigned int i = 0; i < dataBuffer.size(); ++i) {
926  (*this)[i] = dataBuffer[i];
927  }
928  }
929 
930  return;
931 }
932 
933 // -------------------------------------------------
934 // added 2/28/13
935 void
936 TeuchosVector::mpiAllReduce(RawType_MPI_Op mpiOperation, const MpiComm& opComm, TeuchosVector& resultVec) const
937 {
938  // Filter out those nodes that should not participate
939  if (opComm.MyPID() < 0) return;
940 
941  unsigned int size = this->sizeLocal();
942  UQ_FATAL_TEST_MACRO(size != resultVec.sizeLocal(),
943  m_env.worldRank(),
944  "TeuchosVector::mpiAllReduce()",
945  "different vector sizes");
946 
947  for (unsigned int i = 0; i < size; ++i) {
948  double srcValue = (*this)[i];
949  double resultValue = 0.;
950  opComm.Allreduce((void *) &srcValue, (void *) &resultValue, (int) 1, RawValue_MPI_DOUBLE, mpiOperation,
951  "TeuchosVector::mpiAllReduce()",
952  "failed MPI.Allreduce()");
953  resultVec[i] = resultValue;
954  }
955 
956  return;
957 }
958 
959 // -------------------------------------------------
960 // added 2/28/13
961 void
962 TeuchosVector::mpiAllQuantile(double probability, const MpiComm& opComm, TeuchosVector& resultVec) const
963 {
964  // Filter out those nodes that should not participate
965  if (opComm.MyPID() < 0) return;
966 
967  UQ_FATAL_TEST_MACRO((probability < 0.) || (1. < probability),
968  m_env.worldRank(),
969  "TeuchosVector::mpiAllQuantile()",
970  "invalid input");
971 
972  unsigned int size = this->sizeLocal();
973  UQ_FATAL_TEST_MACRO(size != resultVec.sizeLocal(),
974  m_env.worldRank(),
975  "TeuchosVector::mpiAllQuantile()",
976  "different vector sizes");
977 
978  for (unsigned int i = 0; i < size; ++i) {
979  double auxDouble = (int) (*this)[i];
980  std::vector<double> vecOfDoubles(opComm.NumProc(),0.);
981  opComm.Gather((void *) &auxDouble, 1, RawValue_MPI_DOUBLE, (void *) &vecOfDoubles[0], (int) 1, RawValue_MPI_DOUBLE, 0,
982  "TeuchosVector::mpiAllQuantile()",
983  "failed MPI.Gather()");
984 
985  std::sort(vecOfDoubles.begin(), vecOfDoubles.end());
986 
987  double result = vecOfDoubles[(unsigned int)( probability*((double)(vecOfDoubles.size()-1)) )];
988 
989  opComm.Bcast((void *) &result, (int) 1, RawValue_MPI_DOUBLE, 0,
990  "TeuchosVector::mpiAllQuantile()",
991  "failed MPI.Bcast()");
992 
993  resultVec[i] = result;
994  }
995 
996  return;
997 }
998 
999 // -------------------------------------------------
1000 // added/tested 2/28/13
1001 void
1002 TeuchosVector::matlabLinearInterpExtrap(
1003  const TeuchosVector& x1Vec,
1004  const TeuchosVector& y1Vec,
1005  const TeuchosVector& x2Vec)
1006 {
1007  UQ_FATAL_TEST_MACRO(x1Vec.sizeLocal() <= 1,
1008  m_env.worldRank(),
1009  "TeuchosVector::matlabLinearInterpExtrap()",
1010  "invalid 'x1' size");
1011 
1012  UQ_FATAL_TEST_MACRO(x1Vec.sizeLocal() != y1Vec.sizeLocal(),
1013  m_env.worldRank(),
1014  "TeuchosVector::matlabLinearInterpExtrap()",
1015  "invalid 'x1' and 'y1' sizes");
1016 
1017  UQ_FATAL_TEST_MACRO(x2Vec.sizeLocal() != this->sizeLocal(),
1018  m_env.worldRank(),
1019  "TeuchosVector::matlabLinearInterpExtrap()",
1020  "invalid 'x2' and 'this' sizes");
1021 
1022  for (unsigned int i = 1; i < x1Vec.sizeLocal(); ++i) { // Yes, '1'
1023  UQ_FATAL_TEST_MACRO(x1Vec[i] <= x1Vec[i-1],
1024  m_env.worldRank(),
1025  "TeuchosVector::matlabLinearInterpExtrap()",
1026  "invalid 'x1' values");
1027  }
1028 
1029  for (unsigned int id2 = 0; id2 < x2Vec.sizeLocal(); ++id2) {
1030  double x2 = x2Vec[id2];
1031  unsigned int id1 = 0;
1032  bool found1 = false;
1033  for (id1 = 0; id1 < x1Vec.sizeLocal(); ++id1) {
1034  if (x2 <= x1Vec[id1]) {
1035  found1 = true;
1036  break;
1037  }
1038  }
1039  bool makeLinearModel = false;
1040  double xa = 0.;
1041  double xb = 0.;
1042  double ya = 0.;
1043  double yb = 0.;
1044  if (x2 == x1Vec[id1]) {
1045  (*this)[id2] = y1Vec[id1];
1046  }
1047  else if (x2 < x1Vec[0]) {
1048  // Extrapolation case
1049  makeLinearModel = true;
1050  xa = x1Vec[0];
1051  xb = x1Vec[1];
1052  ya = y1Vec[0];
1053  yb = y1Vec[1];
1054  }
1055  else if (found1 == true) {
1056  // Interpolation case
1057  makeLinearModel = true;
1058  xa = x1Vec[id1-1];
1059  xb = x1Vec[id1];
1060  ya = y1Vec[id1-1];
1061  yb = y1Vec[id1];
1062  }
1063  else {
1064  // Extrapolation case
1065  makeLinearModel = true;
1066  xa = x1Vec[x1Vec.sizeLocal()-2];
1067  xb = x1Vec[x1Vec.sizeLocal()-1];
1068  ya = y1Vec[x1Vec.sizeLocal()-2];
1069  yb = y1Vec[x1Vec.sizeLocal()-1];
1070  }
1071 
1072  if (makeLinearModel) {
1073  double rate = (yb-ya)/(xb-xa);
1074  (*this)[id2] = ya + (x2-xa)*rate;
1075  }
1076  }
1077  return;
1078 }
1079 
1080 // -------------------------------------------------
1081 // added 2/28/13
1082 void
1083 TeuchosVector::matlabDiff(
1084  unsigned int firstPositionToStoreDiff,
1085  double valueForRemainderPosition,
1086  TeuchosVector& outputVec) const
1087 {
1088  unsigned int size = this->sizeLocal();
1089 
1090  UQ_FATAL_TEST_MACRO(firstPositionToStoreDiff > 1,
1091  m_env.worldRank(),
1092  "TeuchosVector::matlabDiff()",
1093  "invalid firstPositionToStoreDiff");
1094 
1095  UQ_FATAL_TEST_MACRO(size != outputVec.sizeLocal(),
1096  m_env.worldRank(),
1097  "TeuchosVector::matlabDiff()",
1098  "invalid size of outputVecs");
1099 
1100  for (unsigned int i = 0; i < (size-1); ++i) {
1101  outputVec[firstPositionToStoreDiff+i] = (*this)[i+1]-(*this)[i];
1102  }
1103  if (firstPositionToStoreDiff == 0) {
1104  outputVec[size-1] = valueForRemainderPosition;
1105  }
1106  else {
1107  outputVec[0] = valueForRemainderPosition;
1108  }
1109 
1110  return;
1111 }
1112 
1113 // I/O methods -------------------------------------
1114 // -------------------------------------------------
1115 void
1116 TeuchosVector::print(std::ostream& os) const
1117 {
1118  unsigned int size = this->sizeLocal();
1119 
1120  std::ostream::fmtflags curr_fmt = os.flags();
1121 
1122  if (m_printScientific) {
1123  unsigned int savedPrecision = os.precision();
1124  os.precision(16);
1125 
1126  if (m_printHorizontally) {
1127  for (unsigned int i = 0; i < size; ++i) {
1128  os << std::scientific << (*this)[i]
1129  << " ";
1130  }
1131  }
1132  else {
1133  for (unsigned int i = 0; i < size; ++i) {
1134  os << std::scientific << (*this)[i]
1135  << std::endl;
1136  }
1137  }
1138 
1139  os.precision(savedPrecision);
1140  }
1141  else {
1142  if (m_printHorizontally) {
1143  for (unsigned int i = 0; i < size; ++i) {
1144  os << std::dec << (*this)[i]
1145  << " ";
1146  }
1147  }
1148  else {
1149  for (unsigned int i = 0; i < size; ++i) {
1150  os << std::dec << (*this)[i]
1151  << std::endl;
1152  }
1153  }
1154  }
1155 
1156  os.flags(curr_fmt);
1157  return;
1158 }
1159 
1160 // -------------------------------------------------
1161 void
1162 TeuchosVector::subReadContents(
1163  const std::string& fileName,
1164  const std::string& fileType,
1165  const std::set<unsigned int>& allowedSubEnvIds)
1166 {
1167  UQ_FATAL_TEST_MACRO(m_env.subRank() < 0,
1168  m_env.worldRank(),
1169  "TeuchosVector::subReadContents()",
1170  "unexpected subRank");
1171 
1172  UQ_FATAL_TEST_MACRO(this->numOfProcsForStorage() > 1,
1173  m_env.worldRank(),
1174  "TeuchosVector::subReadContents()",
1175  "implemented just for sequential vectors for now");
1176 
1177  FilePtrSetStruct filePtrSet;
1178  if (m_env.openInputFile(fileName,
1179  fileType, // "m or hdf"
1180  allowedSubEnvIds,
1181  filePtrSet)) {
1182  double subReadSize = this->sizeLocal();
1183 
1184  // In the logic below, the id of a line' begins with value 0 (zero)
1185  unsigned int idOfMyFirstLine = 1;
1186  unsigned int idOfMyLastLine = this->sizeLocal();
1187  unsigned int numParams = 1; // Yes, just '1'
1188 
1189  // Read number of chain positions in the file by taking care of the first line,
1190  // which resembles something like 'variable_name = zeros(n_positions,m_params);'
1191  std::string tmpString;
1192 
1193  // Read 'variable name' string
1194  *filePtrSet.ifsVar >> tmpString;
1195  //std::cout << "Just read '" << tmpString << "'" << std::endl;
1196 
1197  // Read '=' sign
1198  *filePtrSet.ifsVar >> tmpString;
1199  //std::cout << "Just read '" << tmpString << "'" << std::endl;
1200  UQ_FATAL_TEST_MACRO(tmpString != "=",
1201  m_env.worldRank(),
1202  "TeuchosVector::subReadContents()",
1203  "string should be the '=' sign");
1204 
1205  // Read 'zeros(n_positions,n_params)' string
1206  *filePtrSet.ifsVar >> tmpString;
1207  //std::cout << "Just read '" << tmpString << "'" << std::endl;
1208  unsigned int posInTmpString = 6;
1209 
1210  // Isolate 'n_positions' in a string
1211  char nPositionsString[tmpString.size()-posInTmpString+1];
1212  unsigned int posInPositionsString = 0;
1213  do {
1214  UQ_FATAL_TEST_MACRO(posInTmpString >= tmpString.size(),
1215  m_env.worldRank(),
1216  "TeuchosVector::subReadContents()",
1217  "symbol ',' not found in first line of file");
1218  nPositionsString[posInPositionsString++] = tmpString[posInTmpString++];
1219  } while (tmpString[posInTmpString] != ',');
1220  nPositionsString[posInPositionsString] = '\0';
1221 
1222  // Isolate 'n_params' in a string
1223  posInTmpString++; // Avoid reading ',' char
1224  char nParamsString[tmpString.size()-posInTmpString+1];
1225  unsigned int posInParamsString = 0;
1226  do {
1227  UQ_FATAL_TEST_MACRO(posInTmpString >= tmpString.size(),
1228  m_env.worldRank(),
1229  "TeuchosVector::subReadContents()",
1230  "symbol ')' not found in first line of file");
1231  nParamsString[posInParamsString++] = tmpString[posInTmpString++];
1232  } while (tmpString[posInTmpString] != ')');
1233  nParamsString[posInParamsString] = '\0';
1234 
1235  // Convert 'n_positions' and 'n_params' strings to numbers
1236  unsigned int sizeOfVecInFile = (unsigned int) strtod(nPositionsString,NULL);
1237  unsigned int numParamsInFile = (unsigned int) strtod(nParamsString, NULL);
1238  if (m_env.subDisplayFile()) {
1239  *m_env.subDisplayFile() << "In TeuchosVector::subReadContents()"
1240  << ": fullRank " << m_env.fullRank()
1241  << ", sizeOfVecInFile = " << sizeOfVecInFile
1242  << ", numParamsInFile = " << numParamsInFile
1243  << ", this->sizeLocal() = " << this->sizeLocal()
1244  << std::endl;
1245  }
1246 
1247  // Check if [size of vec in file] >= [requested sub vec size]
1248  UQ_FATAL_TEST_MACRO(sizeOfVecInFile < subReadSize,
1249  m_env.worldRank(),
1250  "TeuchosVector::subReadContents()",
1251  "size of vec in file is not big enough");
1252 
1253  // Check if [num params in file] == [num params in current vec]
1254  UQ_FATAL_TEST_MACRO(numParamsInFile != numParams,
1255  m_env.worldRank(),
1256  "TeuchosVector::subReadContents()",
1257  "number of parameters of vec in file is different than number of parameters in this vec object");
1258 
1259  // Code common to any core in a communicator
1260  unsigned int maxCharsPerLine = 64*numParams; // Up to about 60 characters to represent each parameter value
1261 
1262  unsigned int lineId = 0;
1263  while (lineId < idOfMyFirstLine) {
1264  filePtrSet.ifsVar->ignore(maxCharsPerLine,'\n');
1265  lineId++;
1266  };
1267 
1268  if (m_env.subDisplayFile()) {
1269  *m_env.subDisplayFile() << "In TeuchosVector::subReadContents()"
1270  << ": beginning to read input actual data"
1271  << std::endl;
1272  }
1273 
1274  // Take care of initial part of the first data line,
1275  // which resembles something like 'variable_name = [value1 value2 ...'
1276  // Read 'variable name' string
1277  *filePtrSet.ifsVar >> tmpString;
1278  //std::cout << "Core 0 just read '" << tmpString << "'" << std::endl;
1279 
1280  // Read '=' sign
1281  *filePtrSet.ifsVar >> tmpString;
1282  //std::cout << "Core 0 just read '" << tmpString << "'" << std::endl;
1283  UQ_FATAL_TEST_MACRO(tmpString != "=",
1284  m_env.worldRank(),
1285  "TeuchosVector::subReadContents()",
1286  "in core 0, string should be the '=' sign");
1287 
1288  // Take into account the ' [' portion
1289  std::streampos tmpPos = filePtrSet.ifsVar->tellg();
1290  filePtrSet.ifsVar->seekg(tmpPos+(std::streampos)2);
1291 
1292  if (m_env.subDisplayFile()) {
1293  *m_env.subDisplayFile() << "In TeuchosVector::subReadContents()"
1294  << ": beginning to read lines with numbers only"
1295  << ", lineId = " << lineId
1296  << ", idOfMyFirstLine = " << idOfMyFirstLine
1297  << ", idOfMyLastLine = " << idOfMyLastLine
1298  << std::endl;
1299  }
1300 
1301  while (lineId <= idOfMyLastLine) {
1302  *filePtrSet.ifsVar >> (*this)[lineId - idOfMyFirstLine];
1303  lineId++;
1304  };
1305 
1306  m_env.closeFile(filePtrSet,fileType);
1307  }
1308 
1309  return;
1310 }
1311 
1312 // -------------------------------------------------
1313 void
1314 TeuchosVector::subWriteContents(
1315  const std::string& varNamePrefix,
1316  const std::string& fileName,
1317  const std::string& fileType,
1318  const std::set<unsigned int>& allowedSubEnvIds) const
1319 {
1320  UQ_FATAL_TEST_MACRO(m_env.subRank() < 0,
1321  m_env.worldRank(),
1322  "TeuchosVector::subWriteContents()",
1323  "unexpected subRank");
1324 
1325  UQ_FATAL_TEST_MACRO(this->numOfProcsForStorage() > 1,
1326  m_env.worldRank(),
1327  "TeuchosVector::subWriteContents()",
1328  "implemented just for sequential vectors for now");
1329 
1330  FilePtrSetStruct filePtrSet;
1331  if (m_env.openOutputFile(fileName,
1332  fileType, // "m or hdf"
1333  allowedSubEnvIds,
1334  false,
1335  filePtrSet)) {
1336  *filePtrSet.ofsVar << varNamePrefix << "_sub" << m_env.subIdString() << " = zeros(" << this->sizeLocal()
1337  << "," << 1
1338  << ");"
1339  << std::endl;
1340  *filePtrSet.ofsVar << varNamePrefix << "_sub" << m_env.subIdString() << " = [";
1341 
1342  bool savedVectorPrintScientific = this->getPrintScientific();
1343  bool savedVectorPrintHorizontally = this->getPrintHorizontally();
1344  this->setPrintScientific (true);
1345  this->setPrintHorizontally(false);
1346  *filePtrSet.ofsVar << *this;
1347  //<< std::endl; // No need for 'endl' because horizontally = 'false'
1348  this->setPrintHorizontally(savedVectorPrintHorizontally);
1349  this->setPrintScientific (savedVectorPrintScientific);
1350 
1351  *filePtrSet.ofsVar << "];\n";
1352 
1353  m_env.closeFile(filePtrSet,fileType);
1354  }
1355 
1356  return;
1357 }
1358 
1359 
1360 // -------------------------------------------------
1361 // Private metfods ---------------------------------
1362 // -------------------------------------------------
1363 
1364 void
1365 TeuchosVector::copy(const TeuchosVector& rhs)
1366 {
1367  this->Vector::copy(rhs);
1368 
1369  unsigned int size1 = m_vec.length();
1370  unsigned int size2 = rhs.sizeLocal();
1371 
1372  if (size1==size2){
1373  for (unsigned int i=0;i<size1;i++){
1374  m_vec[i]=rhs[i];
1375  }
1376  }
1377  return;
1378 }
1379 
1380 
1381 // -------------------------------------------------
1382 // Operators outside class definition --------------
1383 // -------------------------------------------------
1384 
1385 // -------------------------------------------------
1386 TeuchosVector operator/(double a, const TeuchosVector& x)
1387 {
1388  TeuchosVector answer(x); //copy x to answer
1389  answer.cwInvert();
1390  answer *= a;
1391  return answer;
1392 }
1393 
1394 // -------------------------------------------------
1395 TeuchosVector operator/(const TeuchosVector& x, const TeuchosVector& y)
1396 {
1397  TeuchosVector answer(x);
1398  answer /= y;
1399  return answer;
1400 }
1401 
1402 // -------------------------------------------------
1403 TeuchosVector operator*(double a, const TeuchosVector& x)
1404 {
1405  TeuchosVector answer(x);
1406  answer *= a;
1407  return answer;
1408 }
1409 
1410 // -------------------------------------------------
1411 TeuchosVector operator*(const TeuchosVector& x, const TeuchosVector& y)
1412 {
1413  TeuchosVector answer(x);
1414  answer *= y;
1415  return answer;
1416 }
1417 
1418 // -------------------------------------------------
1419 double scalarProduct(const TeuchosVector& x, const TeuchosVector& y)
1420 {
1421  unsigned int size1 = x.sizeLocal();
1422  unsigned int size2 = y.sizeLocal();
1423 
1424  UQ_FATAL_TEST_MACRO((size1 != size2),
1425  x.env().worldRank(),
1426  "scalarProduct()",
1427  "different sizes of x and y");
1428  double result = 0.;
1429  for (unsigned int i = 0; i < size1; ++i) {
1430  result += x[i]*y[i];
1431  }
1432 
1433  return result;
1434 }
1435 
1436 // -------------------------------------------------
1437 TeuchosVector operator+(const TeuchosVector& x, const TeuchosVector& y)
1438 {
1439  TeuchosVector answer(x);
1440  answer += y;
1441  return answer;
1442 }
1443 
1444 // -------------------------------------------------
1445 TeuchosVector operator-(const TeuchosVector& x, const TeuchosVector& y)
1446 {
1447  TeuchosVector answer(x);
1448  answer -= y;
1449  return answer;
1450 }
1451 
1452 // -------------------------------------------------
1453 bool
1454 operator== (const TeuchosVector& lhs, const TeuchosVector& rhs)
1455 {
1456  bool answer = true;
1457  unsigned int size1 = lhs.sizeLocal();
1458  unsigned int size2 = rhs.sizeLocal();
1459 
1460  UQ_FATAL_TEST_MACRO((size1 != size2),
1461  lhs.env().worldRank(),
1462  "operator==()",
1463  "different sizes of lhs and rhs");
1464 
1465  for (unsigned int i = 0; i < size1; ++i) {
1466  if (lhs[i] != rhs[i]) {
1467  answer = false;
1468  break;
1469  }
1470  }
1471 
1472  return answer;
1473 }
1474 
1475 // -------------------------------------------------
1476 std::ostream&
1477 operator<<(std::ostream& os, const TeuchosVector& obj)
1478 {
1479  obj.print(os);
1480 
1481  return os;
1482 }
1483 
1484 // -------------------------------------------------
1485 
1486 } // End namespace QUESO
1487 
1488 #endif // ifdef QUESO_HAS_TRILINOS
RawValue_MPI_DOUBLE
#define RawValue_MPI_DOUBLE
Definition: MpiComm.h:48
QUESO::operator+
GslMatrix operator+(const GslMatrix &m1, const GslMatrix &m2)
Definition: GslMatrix.C:2395
QUESO::operator-
GslMatrix operator-(const GslMatrix &m1, const GslMatrix &m2)
Definition: GslMatrix.C:2402
QUESO::scalarProduct
double scalarProduct(const GslVector &x, const GslVector &y)
Definition: GslVector.C:1360
QUESO::operator<<
std::ostream & operator<<(std::ostream &os, const BaseEnvironment &obj)
RawValue_MPI_SUM
#define RawValue_MPI_SUM
Definition: MpiComm.h:52
QUESO::operator==
bool operator==(const GslVector &lhs, const GslVector &rhs)
Definition: GslVector.C:1394
QUESO::operator/
GslVector operator/(double a, const GslVector &x)
Definition: GslVector.C:1327
QUESO::operator*
GslMatrix operator*(double a, const GslMatrix &mat)
Definition: GslMatrix.C:2352
QUESO::Vector::copy
virtual void copy(const Vector &src)
Copies vector src to this matrix.
Definition: Vector.C:183
UQ_FATAL_TEST_MACRO
#define UQ_FATAL_TEST_MACRO(test, givenRank, where, what)
Definition: Defines.h:223
QUESO::RawType_MPI_Op
MPI_Op RawType_MPI_Op
Definition: MpiComm.h:41
Generated on Thu Apr 23 2015 19:26:15 for queso-0.51.1 by  doxygen 1.8.5