Anasazi  Version of the Day
AnasaziSaddleContainer.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Anasazi: Block Eigensolvers Package
5 // Copyright (2004) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // This library is free software; you can redistribute it and/or modify
11 // it under the terms of the GNU Lesser General Public License as
12 // published by the Free Software Foundation; either version 2.1 of the
13 // License, or (at your option) any later version.
14 //
15 // This library is distributed in the hope that it will be useful, but
16 // WITHOUT ANY WARRANTY; without even the implied warranty of
17 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18 // Lesser General Public License for more details.
19 //
20 // You should have received a copy of the GNU Lesser General Public
21 // License along with this library; if not, write to the Free Software
22 // Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301
23 // USA
24 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
25 //
26 // ***********************************************************************
27 // @HEADER
28 
36 #ifndef ANASAZI_SADDLE_CONTAINER_HPP
37 #define ANASAZI_SADDLE_CONTAINER_HPP
38 
39 #include "AnasaziConfigDefs.hpp"
40 #include "Teuchos_VerboseObject.hpp"
41 
42 #ifdef HAVE_ANASAZI_BELOS
43 #include "BelosConfigDefs.hpp"
44 #include "BelosMultiVecTraits.hpp"
45 #endif
46 
47 using Teuchos::RCP;
48 using Teuchos::rcp;
49 
50 namespace Anasazi {
51 namespace Experimental {
52 
53 template <class ScalarType, class MV>
54 class SaddleContainer //: public Anasazi::SaddleContainer<ScalarType, MV>
55 {
56  template <class Scalar_, class MV_, class OP_> friend class SaddleOperator;
57 
58 private:
60  typedef Teuchos::SerialDenseMatrix<int,ScalarType> SerialDenseMatrix;
61  const ScalarType ONE, ZERO;
62  RCP<MV> upper_;
63  RCP<SerialDenseMatrix> lowerRaw_;
64  std::vector<int> indices_;
65 
66  RCP<SerialDenseMatrix> getLower() const;
67  void setLower(const RCP<SerialDenseMatrix> lower);
68 
69 public:
70  // Constructors
71  SaddleContainer( ) : ONE(Teuchos::ScalarTraits<ScalarType>::one()), ZERO(Teuchos::ScalarTraits<ScalarType>::zero()) { };
72  SaddleContainer( const RCP<MV> X, bool eye=false );
73 
74  // Things that are necessary for compilation
75  // Returns a clone of the current vector
76  SaddleContainer<ScalarType, MV> * Clone(const int nvecs) const;
77  // Returns a duplicate of the current vector
78  SaddleContainer<ScalarType, MV> * CloneCopy() const;
79  // Returns a duplicate of the specified vectors
80  SaddleContainer<ScalarType, MV> * CloneCopy(const std::vector<int> &index) const;
81  SaddleContainer<ScalarType, MV> * CloneViewNonConst(const std::vector<int> &index) const;
82  SaddleContainer<ScalarType, MV> * CloneViewNonConst(const Teuchos::Range1D& index) const;
83  const SaddleContainer<ScalarType, MV> * CloneView(const std::vector<int> &index) const;
84  const SaddleContainer<ScalarType, MV> * CloneView(const Teuchos::Range1D& index) const;
85  ptrdiff_t GetGlobalLength() const { return MVT::GetGlobalLength(*upper_) + lowerRaw_->numRows(); };
86  int GetNumberVecs() const { return MVT::GetNumberVecs(*upper_); };
87  // Update *this with alpha * A * B + beta * (*this)
88  void MvTimesMatAddMv(ScalarType alpha, const SaddleContainer<ScalarType,MV> &A,
90  ScalarType beta);
91  // Replace *this with alpha * A + beta * B
92  void MvAddMv(ScalarType alpha, const SaddleContainer<ScalarType,MV>& A,
93  ScalarType beta, const SaddleContainer<ScalarType,MV>& B);
94  // Scale the vectors by alpha
95  void MvScale( ScalarType alpha );
96  // Scale the i-th vector by alpha[i]
97  void MvScale( const std::vector<ScalarType>& alpha );
98  // Compute a dense matrix B through the matrix-matrix multiply alpha * A^H * (*this)
99  void MvTransMv (ScalarType alpha, const SaddleContainer<ScalarType, MV>& A,
101  // Compute a vector b where the components are the individual dot-products, i.e.b[i] = A[i]^H*this[i] where A[i] is the i-th column of A.
102  void MvDot (const SaddleContainer<ScalarType, MV>& A, std::vector<ScalarType> &b) const;
103  // Compute the 2-norm of each individual vector
104  void MvNorm ( std::vector<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType> &normvec) const;
105  // Copy the vectors in A to a set of vectors in *this. The numvecs vectors in
106  // A are copied to a subset of vectors in *this indicated by the indices given
107  // in index.
108  void SetBlock (const SaddleContainer<ScalarType, MV>& A, const std::vector<int> &index);
109  // Deep copy.
110  void Assign (const SaddleContainer<ScalarType, MV>&A);
111  // Fill the vectors in *this with random numbers.
112  void MvRandom ();
113  // Replace each element of the vectors in *this with alpha.
114  void MvInit (ScalarType alpha);
115  // Prints the multivector to an output stream
116  void MvPrint (std::ostream &os) const;
117 };
118 
119 
120 
121 // THIS IS NEW!
122 template <class ScalarType, class MV>
123 RCP< Teuchos::SerialDenseMatrix<int,ScalarType> > SaddleContainer<ScalarType, MV>::getLower() const
124 {
125  if(indices_.empty())
126  {
127  return lowerRaw_;
128  }
129 
130  int nrows = lowerRaw_->numRows();
131  int ncols = indices_.size();
132 
133  RCP<SerialDenseMatrix> lower = rcp(new SerialDenseMatrix(nrows,ncols,false));
134 
135  for(int r=0; r<nrows; r++)
136  {
137  for(int c=0; c<ncols; c++)
138  {
139  (*lower)(r,c) = (*lowerRaw_)(r,indices_[c]);
140  }
141  }
142 
143  return lower;
144 }
145 
146 
147 
148 // THIS IS NEW!
149 template <class ScalarType, class MV>
150 void SaddleContainer<ScalarType, MV>::setLower(const RCP<SerialDenseMatrix> lower)
151 {
152  // If the indices are empty, lower points to lowerRaw
153  if(indices_.empty())
154  {
155  return;
156  }
157 
158  int nrows = lowerRaw_->numRows();
159  int ncols = indices_.size();
160 
161  for(int r=0; r<nrows; r++)
162  {
163  for(int c=0; c<ncols; c++)
164  {
165  (*lowerRaw_)(r,indices_[c]) = (*lower)(r,c);
166  }
167  }
168 }
169 
170 
171 
172 // Constructor
173 template <class ScalarType, class MV>
174 SaddleContainer<ScalarType, MV>::SaddleContainer( const RCP<MV> X, bool eye )
175 : ONE(Teuchos::ScalarTraits<ScalarType>::one()), ZERO(Teuchos::ScalarTraits<ScalarType>::zero())
176 {
177  int nvecs = MVT::GetNumberVecs(*X);
178 
179  if(eye)
180  {
181  // Initialize upper_ as all 0s
182  upper_ = MVT::Clone(*X, nvecs);
183  MVT::MvInit(*upper_);
184 
185  // Initialize Y to be I
186  lowerRaw_ = rcp(new SerialDenseMatrix(nvecs,nvecs));
187  for(int i=0; i < nvecs; i++)
188  (*lowerRaw_)(i,i) = ONE;
189  }
190  else
191  {
192  // Point upper_ to X
193  upper_ = X;
194 
195  // Initialize Y to be 0
196  lowerRaw_ = rcp(new SerialDenseMatrix(nvecs,nvecs));
197  }
198 }
199 
200 
201 
202 // Returns a clone of the current vector
203 template <class ScalarType, class MV>
204 SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::Clone(const int nvecs) const
205 {
206  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
207 
208  newSC->upper_ = MVT::Clone(*upper_,nvecs);
209  newSC->lowerRaw_ = rcp(new SerialDenseMatrix(lowerRaw_->numRows(),nvecs));
210 
211  return newSC;
212 }
213 
214 
215 
216 // Returns a duplicate of the current vector
217 template <class ScalarType, class MV>
218 SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::CloneCopy() const
219 {
220  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
221 
222  newSC->upper_ = MVT::CloneCopy(*upper_);
223  newSC->lowerRaw_ = getLower();
224 
225  return newSC;
226 }
227 
228 
229 
230 // Returns a duplicate of the specified vectors
231 template <class ScalarType, class MV>
232 SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::CloneCopy(const std::vector< int > &index) const
233 {
234  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
235 
236  newSC->upper_ = MVT::CloneCopy(*upper_,index);
237 
238  int ncols = index.size();
239  int nrows = lowerRaw_->numRows();
240  RCP<SerialDenseMatrix> lower = getLower();
241  newSC->lowerRaw_ = rcp(new SerialDenseMatrix(nrows,ncols));
242  for(int c=0; c < ncols; c++)
243  {
244  for(int r=0; r < nrows; r++)
245  (*newSC->lowerRaw_)(r,c) = (*lower)(r,index[c]);
246  }
247 
248  return newSC;
249 }
250 
251 
252 
253 // THIS IS NEW!
254 template <class ScalarType, class MV>
255 SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::CloneViewNonConst(const std::vector<int> &index) const
256 {
257  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
258 
259  newSC->upper_ = MVT::CloneViewNonConst(*upper_,index);
260 
261  newSC->lowerRaw_ = lowerRaw_;
262 
263  if(!indices_.empty())
264  {
265  newSC->indices_.resize(index.size());
266  for(unsigned int i=0; i<index.size(); i++)
267  {
268  newSC->indices_[i] = indices_[index[i]];
269  }
270  }
271  else
272  {
273  newSC->indices_ = index;
274  }
275 
276  return newSC;
277 }
278 
279 
280 template <class ScalarType, class MV>
281 SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::CloneViewNonConst(const Teuchos::Range1D& index) const
282 {
283  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
284 
285  newSC->upper_ = MVT::CloneViewNonConst(*upper_,index);
286 
287  newSC->lowerRaw_ = lowerRaw_;
288 
289  newSC->indices_.resize(index.size());
290  for(unsigned int i=0; i<index.size(); i++)
291  {
292  newSC->indices_[i] = indices_[index.lbound()+i];
293  }
294 
295  return newSC;
296 }
297 
298 
299 
300 // THIS IS NEW!
301 template <class ScalarType, class MV>
302 const SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::CloneView(const std::vector<int> &index) const
303 {
304  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
305 
306  newSC->upper_ = MVT::CloneViewNonConst(*upper_,index);
307 
308  newSC->lowerRaw_ = lowerRaw_;
309 
310  if(!indices_.empty())
311  {
312  newSC->indices_.resize(index.size());
313  for(unsigned int i=0; i<index.size(); i++)
314  {
315  newSC->indices_[i] = indices_[index[i]];
316  }
317  }
318  else
319  {
320  newSC->indices_ = index;
321  }
322 
323  return newSC;
324 }
325 
326 
327 template <class ScalarType, class MV>
328 const SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::CloneView(const Teuchos::Range1D& index) const
329 {
330  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
331 
332  newSC->upper_ = MVT::CloneViewNonConst(*upper_,index);
333 
334  newSC->lowerRaw_ = lowerRaw_;
335 
336  newSC->indices_.resize(index.size());
337  for(unsigned int i=0; i<index.size(); i++)
338  {
339  newSC->indices_[i] = indices_[index.lbound()+i];
340  }
341 
342  return newSC;
343 }
344 
345 
346 
347 // Update *this with alpha * A * B + beta * (*this)
348 // THIS IS NEW!
349 template <class ScalarType, class MV>
350 void SaddleContainer<ScalarType, MV>::MvTimesMatAddMv(ScalarType alpha, const SaddleContainer<ScalarType,MV> &A,
352  ScalarType beta)
353 {
354  MVT::MvTimesMatAddMv(alpha,*(A.upper_),B,beta,*upper_);
355  RCP<SerialDenseMatrix> lower = getLower();
356  RCP<SerialDenseMatrix> Alower = A.getLower();
357  lower->multiply(Teuchos::NO_TRANS,Teuchos::NO_TRANS,alpha,*Alower,B,beta);
358  setLower(lower);
359 }
360 
361 
362 
363 // Replace *this with alpha * A + beta * B
364 template <class ScalarType, class MV>
365 void SaddleContainer<ScalarType, MV>::MvAddMv(ScalarType alpha, const SaddleContainer<ScalarType,MV>& A,
366  ScalarType beta, const SaddleContainer<ScalarType,MV>& B)
367 {
368  MVT::MvAddMv(alpha, *(A.upper_), beta, *(B.upper_), *upper_);
369 
370  RCP<SerialDenseMatrix> lower = getLower();
371  RCP<SerialDenseMatrix> Alower = A.getLower();
372  RCP<SerialDenseMatrix> Blower = B.getLower();
373 
374  //int ncolsA = Alower->numCols(); // unused
375  //int ncolsThis = lower->numCols(); // unused
376  //int nrows = lower->numRows(); // unused
377 
378  // Y = alpha A
379  lower->assign(*Alower);
380  if(alpha != ONE)
381  lower->scale(alpha);
382  // Y += beta B
383  if(beta == ONE)
384  *lower += *Blower;
385  else if(beta == -ONE)
386  *lower -= *Blower;
387  else if(beta != ZERO)
388  {
389  SerialDenseMatrix scaledB(*Blower);
390  scaledB.scale(beta);
391  *lower += *Blower;
392  }
393 
394  setLower(lower);
395 }
396 
397 
398 
399 // Scale the vectors by alpha
400 template <class ScalarType, class MV>
401 void SaddleContainer<ScalarType, MV>::MvScale( ScalarType alpha )
402 {
403  MVT::MvScale(*upper_, alpha);
404 
405 
406  RCP<SerialDenseMatrix> lower = getLower();
407  lower->scale(alpha);
408  setLower(lower);
409 }
410 
411 
412 
413 // Scale the i-th vector by alpha[i]
414 template <class ScalarType, class MV>
415 void SaddleContainer<ScalarType, MV>::MvScale( const std::vector<ScalarType>& alpha )
416 {
417  MVT::MvScale(*upper_, alpha);
418 
419  RCP<SerialDenseMatrix> lower = getLower();
420 
421  int nrows = lower->numRows();
422  int ncols = lower->numCols();
423 
424  for(int c=0; c<ncols; c++)
425  {
426  for(int r=0; r<nrows; r++)
427  (*lower)(r,c) *= alpha[c];
428  }
429 
430  setLower(lower);
431 }
432 
433 
434 
435 // Compute a dense matrix B through the matrix-matrix multiply alpha * A^H * (*this)
436 // THIS IS NEW!
437 template <class ScalarType, class MV>
438 void SaddleContainer<ScalarType, MV>::MvTransMv (ScalarType alpha, const SaddleContainer<ScalarType, MV>& A,
440 {
441  MVT::MvTransMv(alpha,*(A.upper_),*upper_,B);
442  RCP<SerialDenseMatrix> lower = getLower();
443  RCP<SerialDenseMatrix> Alower = A.getLower();
444  B.multiply(Teuchos::TRANS,Teuchos::NO_TRANS,alpha,*(Alower),*lower,ONE);
445 }
446 
447 
448 
449 // Compute a vector b where the components are the individual dot-products, i.e.b[i] = A[i]^H*this[i] where A[i] is the i-th column of A.
450 template <class ScalarType, class MV>
451 void SaddleContainer<ScalarType, MV>::MvDot (const SaddleContainer<ScalarType, MV>& A, std::vector<ScalarType> &b) const
452 {
453  MVT::MvDot(*upper_, *(A.upper_), b);
454 
455  RCP<SerialDenseMatrix> lower = getLower();
456  RCP<SerialDenseMatrix> Alower = A.getLower();
457 
458  int nrows = lower->numRows();
459  int ncols = lower->numCols();
460 
461  for(int c=0; c < ncols; c++)
462  {
463  for(int r=0; r < nrows; r++)
464  {
465  b[c] += ((*Alower)(r,c) * (*lower)(r,c));
466  }
467  }
468 }
469 
470 
471 
472 // Compute the 2-norm of each individual vector
473 // THIS IS NEW!
474 template <class ScalarType, class MV>
475 void SaddleContainer<ScalarType, MV>::MvNorm ( std::vector<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType> &normvec) const
476 {
477  // TODO: Make this better
478  MvDot(*this,normvec);
479  for(unsigned int i=0; i<normvec.size(); i++)
480  normvec[i] = sqrt(normvec[i]);
481 }
482 
483 
484 
485 // Copy the vectors in A to a set of vectors in *this. The numvecs vectors in
486 // A are copied to a subset of vectors in *this indicated by the indices given
487 // in index.
488 template <class ScalarType, class MV>
489 void SaddleContainer<ScalarType, MV>::SetBlock (const SaddleContainer<ScalarType, MV>& A, const std::vector<int> &index)
490 {
491  MVT::SetBlock(*(A.upper_), index, *upper_);
492 
493  RCP<SerialDenseMatrix> lower = getLower();
494  RCP<SerialDenseMatrix> Alower = A.getLower();
495 
496  int nrows = lower->numRows();
497 
498  int nvecs = index.size();
499  for(int c=0; c<nvecs; c++)
500  {
501  for(int r=0; r<nrows; r++)
502  (*lower)(r,index[c]) = (*Alower)(r,c);
503  }
504 
505  setLower(lower);
506 }
507 
508 
509 
510 // Deep copy.
511 template <class ScalarType, class MV>
512 void SaddleContainer<ScalarType, MV>::Assign (const SaddleContainer<ScalarType, MV>&A)
513 {
514  MVT::Assign(*(A.upper_),*(upper_));
515 
516  RCP<SerialDenseMatrix> lower = getLower();
517  RCP<SerialDenseMatrix> Alower = A.getLower();
518 
519  *lower = *Alower; // This is a well-defined operator for SerialDenseMatrix
520 
521  setLower(lower);
522 }
523 
524 
525 
526 // Fill the vectors in *this with random numbers.
527 // THIS IS NEW!
528 template <class ScalarType, class MV>
529 void SaddleContainer<ScalarType, MV>::MvRandom ()
530 {
531  MVT::MvRandom(*upper_);
532 
533  RCP<SerialDenseMatrix> lower = getLower();
534  lower->random();
535  setLower(lower);
536 }
537 
538 
539 
540 // Replace each element of the vectors in *this with alpha.
541 template <class ScalarType, class MV>
542 void SaddleContainer<ScalarType, MV>::MvInit (ScalarType alpha)
543 {
544  MVT::MvInit(*upper_,alpha);
545 
546  RCP<SerialDenseMatrix> lower = getLower();
547  lower->putScalar(alpha);
548  setLower(lower);
549 }
550 
551 
552 
553 // Prints the multivector to an output stream
554 template <class ScalarType, class MV>
555 void SaddleContainer<ScalarType, MV>::MvPrint (std::ostream &os) const
556 {
557  RCP<SerialDenseMatrix> lower = getLower();
558  //int nrows = lower->numRows(); // unused
559  //int ncols = lower->numCols(); // unused
560 
561  os << "Object SaddleContainer" << std::endl;
562  os << "X\n";
564 // os << "X\n" << *upper_ << std::endl;
565 
566  os << "Y\n" << *lower << std::endl;
567 }
568 
569 } // End namespace Experimental
570 
571 template<class ScalarType, class MV >
572 class MultiVecTraits<ScalarType,Experimental::SaddleContainer<ScalarType, MV> >
573 {
574 typedef Experimental::SaddleContainer<ScalarType,MV> SC;
575 
576 public:
577  static RCP<SC > Clone( const SC& mv, const int numvecs )
578  { return rcp( const_cast<SC&>(mv).Clone(numvecs) ); }
579 
580  static RCP<SC > CloneCopy( const SC& mv )
581  { return rcp( const_cast<SC&>(mv).CloneCopy() ); }
582 
583  static RCP<SC > CloneCopy( const SC& mv, const std::vector<int>& index )
584  { return rcp( const_cast<SC&>(mv).CloneCopy(index) ); }
585 
586  static ptrdiff_t GetGlobalLength( const SC& mv )
587  { return mv.GetGlobalLength(); }
588 
589  static int GetNumberVecs( const SC& mv )
590  { return mv.GetNumberVecs(); }
591 
592  static void MvTimesMatAddMv( ScalarType alpha, const SC& A,
594  ScalarType beta, SC& mv )
595  { mv.MvTimesMatAddMv(alpha, A, B, beta); }
596 
597  static void MvAddMv( ScalarType alpha, const SC& A, ScalarType beta, const SC& B, SC& mv )
598  { mv.MvAddMv(alpha, A, beta, B); }
599 
600  static void MvTransMv( ScalarType alpha, const SC& A, const SC& mv, Teuchos::SerialDenseMatrix<int,ScalarType>& B)
601  { mv.MvTransMv(alpha, A, B); }
602 
603  static void MvDot( const SC& mv, const SC& A, std::vector<ScalarType> & b)
604  { mv.MvDot( A, b); }
605 
606  static void MvScale ( SC& mv, ScalarType alpha )
607  { mv.MvScale( alpha ); }
608 
609  static void MvScale ( SC& mv, const std::vector<ScalarType>& alpha )
610  { mv.MvScale( alpha ); }
611 
612  static void MvNorm( const SC& mv, std::vector<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType> & normvec)
613  { mv.MvNorm(normvec); }
614 
615  static void SetBlock( const SC& A, const std::vector<int>& index, SC& mv )
616  { mv.SetBlock(A, index); }
617 
618  static void Assign( const SC& A, SC& mv )
619  { mv.Assign(A); }
620 
621  static void MvRandom( SC& mv )
622  { mv.MvRandom(); }
623 
624  static void MvInit( SC& mv, ScalarType alpha = Teuchos::ScalarTraits<ScalarType>::zero() )
625  { mv.MvInit(alpha); }
626 
627  static void MvPrint( const SC& mv, std::ostream& os )
628  { mv.MvPrint(os); }
629 };
630 
631 } // end namespace Anasazi
632 
633 #ifdef HAVE_ANASAZI_BELOS
634 namespace Belos
635 {
636 
637 template<class ScalarType, class MV >
638 class MultiVecTraits< ScalarType, Anasazi::Experimental::SaddleContainer<ScalarType,MV> >
639 {
640 typedef Anasazi::Experimental::SaddleContainer<ScalarType,MV> SC;
641 public:
642  static RCP<SC > Clone( const SC& mv, const int numvecs )
643  { return rcp( const_cast<SC&>(mv).Clone(numvecs) ); }
644 
645  static RCP<SC > CloneCopy( const SC& mv )
646  { return rcp( const_cast<SC&>(mv).CloneCopy() ); }
647 
648  static RCP<SC > CloneCopy( const SC& mv, const std::vector<int>& index )
649  { return rcp( const_cast<SC&>(mv).CloneCopy(index) ); }
650 
651  static RCP<SC> CloneViewNonConst( SC& mv, const std::vector<int>& index )
652  { return rcp( mv.CloneViewNonConst(index) ); }
653 
654  static RCP<SC> CloneViewNonConst( SC& mv, const Teuchos::Range1D& index )
655  { return rcp( mv.CloneViewNonConst(index) ); }
656 
657  static RCP<const SC> CloneView( const SC& mv, const std::vector<int> & index )
658  { return rcp( mv.CloneView(index) ); }
659 
660  static RCP<const SC> CloneView( const SC& mv, const Teuchos::Range1D& index )
661  { return rcp( mv.CloneView(index) ); }
662 
663  static ptrdiff_t GetGlobalLength( const SC& mv )
664  { return mv.GetGlobalLength(); }
665 
666  static int GetNumberVecs( const SC& mv )
667  { return mv.GetNumberVecs(); }
668 
669  static void MvTimesMatAddMv( ScalarType alpha, const SC& A,
671  ScalarType beta, SC& mv )
672  { mv.MvTimesMatAddMv(alpha, A, B, beta); }
673 
674  static void MvAddMv( ScalarType alpha, const SC& A, ScalarType beta, const SC& B, SC& mv )
675  { mv.MvAddMv(alpha, A, beta, B); }
676 
677  static void MvTransMv( ScalarType alpha, const SC& A, const SC& mv, Teuchos::SerialDenseMatrix<int,ScalarType>& B)
678  { mv.MvTransMv(alpha, A, B); }
679 
680  static void MvDot( const SC& mv, const SC& A, std::vector<ScalarType> & b)
681  { mv.MvDot( A, b); }
682 
683  static void MvScale ( SC& mv, ScalarType alpha )
684  { mv.MvScale( alpha ); }
685 
686  static void MvScale ( SC& mv, const std::vector<ScalarType>& alpha )
687  { mv.MvScale( alpha ); }
688 
689  // TODO: MAKE SURE TYPE == TWONORM
690  static void MvNorm( const SC& mv, std::vector<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType> & normvec, NormType type=TwoNorm)
691  { mv.MvNorm(normvec); }
692 
693  static void SetBlock( const SC& A, const std::vector<int>& index, SC& mv )
694  { mv.SetBlock(A, index); }
695 
696  static void Assign( const SC& A, SC& mv )
697  { mv.Assign(A); }
698 
699  static void MvRandom( SC& mv )
700  { mv.MvRandom(); }
701 
702  static void MvInit( SC& mv, ScalarType alpha = Teuchos::ScalarTraits<ScalarType>::zero() )
703  { mv.MvInit(alpha); }
704 
705  static void MvPrint( const SC& mv, std::ostream& os )
706  { mv.MvPrint(os); }
707 
708 #ifdef HAVE_BELOS_TSQR
709  typedef Belos::details::StubTsqrAdapter<SC> tsqr_adaptor_type;
722 #endif // HAVE_BELOS_TSQR
723 };
724 
725 } // end namespace Belos
726 #endif
727 
728 #endif
static ptrdiff_t GetGlobalLength(const MV &mv)
Return the number of rows in the given multivector mv.
static void MvInit(MV &mv, const ScalarType alpha=Teuchos::ScalarTraits< ScalarType >::zero())
Replace each element of the vectors in mv with alpha.
int multiply(ETransp transa, ETransp transb, ScalarType alpha, const SerialDenseMatrix< OrdinalType, ScalarType > &A, const SerialDenseMatrix< OrdinalType, ScalarType > &B, ScalarType beta)
Namespace Anasazi contains the classes, structs, enums and utilities used by the Anasazi package...
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
OrdinalType numRows() const
Traits class which defines basic operations on multivectors.
static RCP< FancyOStream > getDefaultOStream()
Anasazi header file which uses auto-configuration information to include necessary C++ headers...
Ordinal lbound() const
static int GetNumberVecs(const MV &mv)
Obtain the number of vectors in mv.
Ordinal size() const
static Teuchos::RCP< MV > Clone(const MV &mv, const int numvecs)
Creates a new empty MV containing numvecs columns.