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 "BelosMultiVecTraits.hpp"
44 #endif
45 
46 using Teuchos::RCP;
47 using Teuchos::rcp;
48 
49 namespace Anasazi {
50 namespace Experimental {
51 
52 template <class ScalarType, class MV>
53 class SaddleContainer //: public Anasazi::SaddleContainer<ScalarType, MV>
54 {
55  template <class Scalar_, class MV_, class OP_> friend class SaddleOperator;
56 
57 private:
59  typedef Teuchos::SerialDenseMatrix<int,ScalarType> SerialDenseMatrix;
60  const ScalarType ONE, ZERO;
61  RCP<MV> upper_;
62  RCP<SerialDenseMatrix> lowerRaw_;
63  std::vector<int> indices_;
64 
65  RCP<SerialDenseMatrix> getLower() const;
66  void setLower(const RCP<SerialDenseMatrix> lower);
67 
68 public:
69  // Constructors
70  SaddleContainer( ) : ONE(Teuchos::ScalarTraits<ScalarType>::one()), ZERO(Teuchos::ScalarTraits<ScalarType>::zero()) { };
71  SaddleContainer( const RCP<MV> X, bool eye=false );
72 
73  // Things that are necessary for compilation
74  // Returns a clone of the current vector
75  SaddleContainer<ScalarType, MV> * Clone(const int nvecs) const;
76  // Returns a duplicate of the current vector
77  SaddleContainer<ScalarType, MV> * CloneCopy() const;
78  // Returns a duplicate of the specified vectors
79  SaddleContainer<ScalarType, MV> * CloneCopy(const std::vector<int> &index) const;
80  SaddleContainer<ScalarType, MV> * CloneViewNonConst(const std::vector<int> &index) const;
81  const SaddleContainer<ScalarType, MV> * CloneView(const std::vector<int> &index) const;
82  ptrdiff_t GetGlobalLength() const { return MVT::GetGlobalLength(*upper_) + lowerRaw_->numRows(); };
83  int GetNumberVecs() const { return MVT::GetNumberVecs(*upper_); };
84  // Update *this with alpha * A * B + beta * (*this)
85  void MvTimesMatAddMv(ScalarType alpha, const SaddleContainer<ScalarType,MV> &A,
86  const Teuchos::SerialDenseMatrix<int, ScalarType> &B,
87  ScalarType beta);
88  // Replace *this with alpha * A + beta * B
89  void MvAddMv(ScalarType alpha, const SaddleContainer<ScalarType,MV>& A,
90  ScalarType beta, const SaddleContainer<ScalarType,MV>& B);
91  // Scale the vectors by alpha
92  void MvScale( ScalarType alpha );
93  // Scale the i-th vector by alpha[i]
94  void MvScale( const std::vector<ScalarType>& alpha );
95  // Compute a dense matrix B through the matrix-matrix multiply alpha * A^H * (*this)
96  void MvTransMv (ScalarType alpha, const SaddleContainer<ScalarType, MV>& A,
97  Teuchos::SerialDenseMatrix< int, ScalarType >& B) const;
98  // 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.
99  void MvDot (const SaddleContainer<ScalarType, MV>& A, std::vector<ScalarType> &b) const;
100  // Compute the 2-norm of each individual vector
101  void MvNorm ( std::vector<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType> &normvec) const;
102  // Copy the vectors in A to a set of vectors in *this. The numvecs vectors in
103  // A are copied to a subset of vectors in *this indicated by the indices given
104  // in index.
105  void SetBlock (const SaddleContainer<ScalarType, MV>& A, const std::vector<int> &index);
106  // Deep copy.
107  void Assign (const SaddleContainer<ScalarType, MV>&A);
108  // Fill the vectors in *this with random numbers.
109  void MvRandom ();
110  // Replace each element of the vectors in *this with alpha.
111  void MvInit (ScalarType alpha);
112  // Prints the multivector to an output stream
113  void MvPrint (std::ostream &os) const;
114 };
115 
116 
117 
118 // THIS IS NEW!
119 template <class ScalarType, class MV>
120 RCP< Teuchos::SerialDenseMatrix<int,ScalarType> > SaddleContainer<ScalarType, MV>::getLower() const
121 {
122  if(indices_.empty())
123  {
124  return lowerRaw_;
125  }
126 
127  int nrows = lowerRaw_->numRows();
128  int ncols = indices_.size();
129 
130  RCP<SerialDenseMatrix> lower = rcp(new SerialDenseMatrix(nrows,ncols,false));
131 
132  for(int r=0; r<nrows; r++)
133  {
134  for(int c=0; c<ncols; c++)
135  {
136  (*lower)(r,c) = (*lowerRaw_)(r,indices_[c]);
137  }
138  }
139 
140  return lower;
141 }
142 
143 
144 
145 // THIS IS NEW!
146 template <class ScalarType, class MV>
147 void SaddleContainer<ScalarType, MV>::setLower(const RCP<SerialDenseMatrix> lower)
148 {
149  // If the indices are empty, lower points to lowerRaw
150  if(indices_.empty())
151  {
152  return;
153  }
154 
155  int nrows = lowerRaw_->numRows();
156  int ncols = indices_.size();
157 
158  for(int r=0; r<nrows; r++)
159  {
160  for(int c=0; c<ncols; c++)
161  {
162  (*lowerRaw_)(r,indices_[c]) = (*lower)(r,c);
163  }
164  }
165 }
166 
167 
168 
169 // Constructor
170 template <class ScalarType, class MV>
171 SaddleContainer<ScalarType, MV>::SaddleContainer( const RCP<MV> X, bool eye )
172 : ONE(Teuchos::ScalarTraits<ScalarType>::one()), ZERO(Teuchos::ScalarTraits<ScalarType>::zero())
173 {
174  int nvecs = MVT::GetNumberVecs(*X);
175 
176  if(eye)
177  {
178  // Initialize upper_ as all 0s
179  upper_ = MVT::Clone(*X, nvecs);
180  MVT::MvInit(*upper_);
181 
182  // Initialize Y to be I
183  lowerRaw_ = rcp(new SerialDenseMatrix(nvecs,nvecs));
184  for(int i=0; i < nvecs; i++)
185  (*lowerRaw_)(i,i) = ONE;
186  }
187  else
188  {
189  // Point upper_ to X
190  upper_ = X;
191 
192  // Initialize Y to be 0
193  lowerRaw_ = rcp(new SerialDenseMatrix(nvecs,nvecs));
194  }
195 }
196 
197 
198 
199 // Returns a clone of the current vector
200 template <class ScalarType, class MV>
201 SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::Clone(const int nvecs) const
202 {
203  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
204 
205  newSC->upper_ = MVT::Clone(*upper_,nvecs);
206  newSC->lowerRaw_ = rcp(new SerialDenseMatrix(lowerRaw_->numRows(),nvecs));
207 
208  return newSC;
209 }
210 
211 
212 
213 // Returns a duplicate of the current vector
214 template <class ScalarType, class MV>
215 SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::CloneCopy() const
216 {
217  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
218 
219  newSC->upper_ = MVT::CloneCopy(*upper_);
220  newSC->lowerRaw_ = getLower();
221 
222  return newSC;
223 }
224 
225 
226 
227 // Returns a duplicate of the specified vectors
228 template <class ScalarType, class MV>
229 SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::CloneCopy(const std::vector< int > &index) const
230 {
231  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
232 
233  newSC->upper_ = MVT::CloneCopy(*upper_,index);
234 
235  int ncols = index.size();
236  int nrows = lowerRaw_->numRows();
237  RCP<SerialDenseMatrix> lower = getLower();
238  newSC->lowerRaw_ = rcp(new SerialDenseMatrix(nrows,ncols));
239  for(int c=0; c < ncols; c++)
240  {
241  for(int r=0; r < nrows; r++)
242  (*newSC->lowerRaw_)(r,c) = (*lower)(r,index[c]);
243  }
244 
245  return newSC;
246 }
247 
248 
249 
250 // THIS IS NEW!
251 template <class ScalarType, class MV>
252 SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::CloneViewNonConst(const std::vector<int> &index) const
253 {
254  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
255 
256  newSC->upper_ = MVT::CloneViewNonConst(*upper_,index);
257 
258  newSC->lowerRaw_ = lowerRaw_;
259 
260  if(!indices_.empty())
261  {
262  newSC->indices_.resize(index.size());
263  for(unsigned int i=0; i<index.size(); i++)
264  {
265  newSC->indices_[i] = indices_[index[i]];
266  }
267  }
268  else
269  {
270  newSC->indices_ = index;
271  }
272 
273  return newSC;
274 }
275 
276 
277 
278 // THIS IS NEW!
279 template <class ScalarType, class MV>
280 const SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::CloneView(const std::vector<int> &index) const
281 {
282  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
283 
284  newSC->upper_ = MVT::CloneViewNonConst(*upper_,index);
285 
286  newSC->lowerRaw_ = lowerRaw_;
287 
288  if(!indices_.empty())
289  {
290  newSC->indices_.resize(index.size());
291  for(unsigned int i=0; i<index.size(); i++)
292  {
293  newSC->indices_[i] = indices_[index[i]];
294  }
295  }
296  else
297  {
298  newSC->indices_ = index;
299  }
300 
301  return newSC;
302 }
303 
304 
305 
306 // Update *this with alpha * A * B + beta * (*this)
307 // THIS IS NEW!
308 template <class ScalarType, class MV>
309 void SaddleContainer<ScalarType, MV>::MvTimesMatAddMv(ScalarType alpha, const SaddleContainer<ScalarType,MV> &A,
310  const Teuchos::SerialDenseMatrix<int, ScalarType> &B,
311  ScalarType beta)
312 {
313  MVT::MvTimesMatAddMv(alpha,*(A.upper_),B,beta,*upper_);
314  RCP<SerialDenseMatrix> lower = getLower();
315  RCP<SerialDenseMatrix> Alower = A.getLower();
316  lower->multiply(Teuchos::NO_TRANS,Teuchos::NO_TRANS,alpha,*Alower,B,beta);
317  setLower(lower);
318 }
319 
320 
321 
322 // Replace *this with alpha * A + beta * B
323 template <class ScalarType, class MV>
324 void SaddleContainer<ScalarType, MV>::MvAddMv(ScalarType alpha, const SaddleContainer<ScalarType,MV>& A,
325  ScalarType beta, const SaddleContainer<ScalarType,MV>& B)
326 {
327  MVT::MvAddMv(alpha, *(A.upper_), beta, *(B.upper_), *upper_);
328 
329  RCP<SerialDenseMatrix> lower = getLower();
330  RCP<SerialDenseMatrix> Alower = A.getLower();
331  RCP<SerialDenseMatrix> Blower = B.getLower();
332 
333  //int ncolsA = Alower->numCols(); // unused
334  //int ncolsThis = lower->numCols(); // unused
335  //int nrows = lower->numRows(); // unused
336 
337  // Y = alpha A
338  lower->assign(*Alower);
339  if(alpha != ONE)
340  lower->scale(alpha);
341  // Y += beta B
342  if(beta == ONE)
343  *lower += *Blower;
344  else if(beta == -ONE)
345  *lower -= *Blower;
346  else if(beta != ZERO)
347  {
348  SerialDenseMatrix scaledB(*Blower);
349  scaledB.scale(beta);
350  *lower += *Blower;
351  }
352 
353  setLower(lower);
354 }
355 
356 
357 
358 // Scale the vectors by alpha
359 template <class ScalarType, class MV>
360 void SaddleContainer<ScalarType, MV>::MvScale( ScalarType alpha )
361 {
362  MVT::MvScale(*upper_, alpha);
363 
364 
365  RCP<SerialDenseMatrix> lower = getLower();
366  lower->scale(alpha);
367  setLower(lower);
368 }
369 
370 
371 
372 // Scale the i-th vector by alpha[i]
373 template <class ScalarType, class MV>
374 void SaddleContainer<ScalarType, MV>::MvScale( const std::vector<ScalarType>& alpha )
375 {
376  MVT::MvScale(*upper_, alpha);
377 
378  RCP<SerialDenseMatrix> lower = getLower();
379 
380  int nrows = lower->numRows();
381  int ncols = lower->numCols();
382 
383  for(int c=0; c<ncols; c++)
384  {
385  for(int r=0; r<nrows; r++)
386  (*lower)(r,c) *= alpha[c];
387  }
388 
389  setLower(lower);
390 }
391 
392 
393 
394 // Compute a dense matrix B through the matrix-matrix multiply alpha * A^H * (*this)
395 // THIS IS NEW!
396 template <class ScalarType, class MV>
397 void SaddleContainer<ScalarType, MV>::MvTransMv (ScalarType alpha, const SaddleContainer<ScalarType, MV>& A,
398  Teuchos::SerialDenseMatrix< int, ScalarType >& B) const
399 {
400  MVT::MvTransMv(alpha,*(A.upper_),*upper_,B);
401  RCP<SerialDenseMatrix> lower = getLower();
402  RCP<SerialDenseMatrix> Alower = A.getLower();
403  B.multiply(Teuchos::TRANS,Teuchos::NO_TRANS,alpha,*(Alower),*lower,ONE);
404 }
405 
406 
407 
408 // 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.
409 template <class ScalarType, class MV>
410 void SaddleContainer<ScalarType, MV>::MvDot (const SaddleContainer<ScalarType, MV>& A, std::vector<ScalarType> &b) const
411 {
412  MVT::MvDot(*upper_, *(A.upper_), b);
413 
414  RCP<SerialDenseMatrix> lower = getLower();
415  RCP<SerialDenseMatrix> Alower = A.getLower();
416 
417  int nrows = lower->numRows();
418  int ncols = lower->numCols();
419 
420  for(int c=0; c < ncols; c++)
421  {
422  for(int r=0; r < nrows; r++)
423  {
424  b[c] += ((*Alower)(r,c) * (*lower)(r,c));
425  }
426  }
427 }
428 
429 
430 
431 // Compute the 2-norm of each individual vector
432 // THIS IS NEW!
433 template <class ScalarType, class MV>
434 void SaddleContainer<ScalarType, MV>::MvNorm ( std::vector<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType> &normvec) const
435 {
436  // TODO: Make this better
437  MvDot(*this,normvec);
438  for(unsigned int i=0; i<normvec.size(); i++)
439  normvec[i] = sqrt(normvec[i]);
440 }
441 
442 
443 
444 // Copy the vectors in A to a set of vectors in *this. The numvecs vectors in
445 // A are copied to a subset of vectors in *this indicated by the indices given
446 // in index.
447 template <class ScalarType, class MV>
448 void SaddleContainer<ScalarType, MV>::SetBlock (const SaddleContainer<ScalarType, MV>& A, const std::vector<int> &index)
449 {
450  MVT::SetBlock(*(A.upper_), index, *upper_);
451 
452  RCP<SerialDenseMatrix> lower = getLower();
453  RCP<SerialDenseMatrix> Alower = A.getLower();
454 
455  int nrows = lower->numRows();
456 
457  int nvecs = index.size();
458  for(int c=0; c<nvecs; c++)
459  {
460  for(int r=0; r<nrows; r++)
461  (*lower)(r,index[c]) = (*Alower)(r,c);
462  }
463 
464  setLower(lower);
465 }
466 
467 
468 
469 // Deep copy.
470 template <class ScalarType, class MV>
471 void SaddleContainer<ScalarType, MV>::Assign (const SaddleContainer<ScalarType, MV>&A)
472 {
473  MVT::Assign(*(A.upper_),*(upper_));
474 
475  RCP<SerialDenseMatrix> lower = getLower();
476  RCP<SerialDenseMatrix> Alower = A.getLower();
477 
478  *lower = *Alower; // This is a well-defined operator for SerialDenseMatrix
479 
480  setLower(lower);
481 }
482 
483 
484 
485 // Fill the vectors in *this with random numbers.
486 // THIS IS NEW!
487 template <class ScalarType, class MV>
488 void SaddleContainer<ScalarType, MV>::MvRandom ()
489 {
490  MVT::MvRandom(*upper_);
491 
492  RCP<SerialDenseMatrix> lower = getLower();
493  lower->random();
494  setLower(lower);
495 }
496 
497 
498 
499 // Replace each element of the vectors in *this with alpha.
500 template <class ScalarType, class MV>
501 void SaddleContainer<ScalarType, MV>::MvInit (ScalarType alpha)
502 {
503  MVT::MvInit(*upper_,alpha);
504 
505  RCP<SerialDenseMatrix> lower = getLower();
506  lower->putScalar(alpha);
507  setLower(lower);
508 }
509 
510 
511 
512 // Prints the multivector to an output stream
513 template <class ScalarType, class MV>
514 void SaddleContainer<ScalarType, MV>::MvPrint (std::ostream &os) const
515 {
516  RCP<SerialDenseMatrix> lower = getLower();
517  //int nrows = lower->numRows(); // unused
518  //int ncols = lower->numCols(); // unused
519 
520  os << "Object SaddleContainer" << std::endl;
521  os << "X\n";
522  upper_->describe(*(Teuchos::VerboseObjectBase::getDefaultOStream()),Teuchos::VERB_EXTREME);
523 // os << "X\n" << *upper_ << std::endl;
524 
525  os << "Y\n" << *lower << std::endl;
526 }
527 
528 } // End namespace Experimental
529 
530 template<class ScalarType, class MV >
531 class MultiVecTraits<ScalarType,Experimental::SaddleContainer<ScalarType, MV> >
532 {
533 typedef Experimental::SaddleContainer<ScalarType,MV> SC;
534 
535 public:
536  static RCP<SC > Clone( const SC& mv, const int numvecs )
537  { return rcp( const_cast<SC&>(mv).Clone(numvecs) ); }
538 
539  static RCP<SC > CloneCopy( const SC& mv )
540  { return rcp( const_cast<SC&>(mv).CloneCopy() ); }
541 
542  static RCP<SC > CloneCopy( const SC& mv, const std::vector<int>& index )
543  { return rcp( const_cast<SC&>(mv).CloneCopy(index) ); }
544 
545  static ptrdiff_t GetGlobalLength( const SC& mv )
546  { return mv.GetGlobalLength(); }
547 
548  static int GetNumberVecs( const SC& mv )
549  { return mv.GetNumberVecs(); }
550 
551  static void MvTimesMatAddMv( ScalarType alpha, const SC& A,
552  const Teuchos::SerialDenseMatrix<int,ScalarType>& B,
553  ScalarType beta, SC& mv )
554  { mv.MvTimesMatAddMv(alpha, A, B, beta); }
555 
556  static void MvAddMv( ScalarType alpha, const SC& A, ScalarType beta, const SC& B, SC& mv )
557  { mv.MvAddMv(alpha, A, beta, B); }
558 
559  static void MvTransMv( ScalarType alpha, const SC& A, const SC& mv, Teuchos::SerialDenseMatrix<int,ScalarType>& B)
560  { mv.MvTransMv(alpha, A, B); }
561 
562  static void MvDot( const SC& mv, const SC& A, std::vector<ScalarType> & b)
563  { mv.MvDot( A, b); }
564 
565  static void MvScale ( SC& mv, ScalarType alpha )
566  { mv.MvScale( alpha ); }
567 
568  static void MvScale ( SC& mv, const std::vector<ScalarType>& alpha )
569  { mv.MvScale( alpha ); }
570 
571  static void MvNorm( const SC& mv, std::vector<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType> & normvec)
572  { mv.MvNorm(normvec); }
573 
574  static void SetBlock( const SC& A, const std::vector<int>& index, SC& mv )
575  { mv.SetBlock(A, index); }
576 
577  static void Assign( const SC& A, SC& mv )
578  { mv.Assign(A); }
579 
580  static void MvRandom( SC& mv )
581  { mv.MvRandom(); }
582 
583  static void MvInit( SC& mv, ScalarType alpha = Teuchos::ScalarTraits<ScalarType>::zero() )
584  { mv.MvInit(alpha); }
585 
586  static void MvPrint( const SC& mv, std::ostream& os )
587  { mv.MvPrint(os); }
588 };
589 
590 } // end namespace Anasazi
591 
592 #ifdef HAVE_ANASAZI_BELOS
593 namespace Belos
594 {
595 
596 template<class ScalarType, class MV >
597 class MultiVecTraits< ScalarType, Anasazi::Experimental::SaddleContainer<ScalarType,MV> >
598 {
599 typedef Anasazi::Experimental::SaddleContainer<ScalarType,MV> SC;
600 public:
601  static RCP<SC > Clone( const SC& mv, const int numvecs )
602  { return rcp( const_cast<SC&>(mv).Clone(numvecs) ); }
603 
604  static RCP<SC > CloneCopy( const SC& mv )
605  { return rcp( const_cast<SC&>(mv).CloneCopy() ); }
606 
607  static RCP<SC > CloneCopy( const SC& mv, const std::vector<int>& index )
608  { return rcp( const_cast<SC&>(mv).CloneCopy(index) ); }
609 
610  static RCP<SC> CloneViewNonConst( SC& mv, const std::vector<int>& index )
611  { return rcp( mv.CloneViewNonConst(index) ); }
612 
613  static RCP<const SC> CloneView( const SC& mv, const std::vector<int> & index )
614  { return rcp( mv.CloneView(index) ); }
615 
616  static ptrdiff_t GetGlobalLength( const SC& mv )
617  { return mv.GetGlobalLength(); }
618 
619  static int GetNumberVecs( const SC& mv )
620  { return mv.GetNumberVecs(); }
621 
622  static void MvTimesMatAddMv( ScalarType alpha, const SC& A,
623  const Teuchos::SerialDenseMatrix<int,ScalarType>& B,
624  ScalarType beta, SC& mv )
625  { mv.MvTimesMatAddMv(alpha, A, B, beta); }
626 
627  static void MvAddMv( ScalarType alpha, const SC& A, ScalarType beta, const SC& B, SC& mv )
628  { mv.MvAddMv(alpha, A, beta, B); }
629 
630  static void MvTransMv( ScalarType alpha, const SC& A, const SC& mv, Teuchos::SerialDenseMatrix<int,ScalarType>& B)
631  { mv.MvTransMv(alpha, A, B); }
632 
633  static void MvDot( const SC& mv, const SC& A, std::vector<ScalarType> & b)
634  { mv.MvDot( A, b); }
635 
636  static void MvScale ( SC& mv, ScalarType alpha )
637  { mv.MvScale( alpha ); }
638 
639  static void MvScale ( SC& mv, const std::vector<ScalarType>& alpha )
640  { mv.MvScale( alpha ); }
641 
642  // TODO: MAKE SURE TYPE == TWONORM
643  static void MvNorm( const SC& mv, std::vector<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType> & normvec, NormType type=TwoNorm)
644  { mv.MvNorm(normvec); }
645 
646  static void SetBlock( const SC& A, const std::vector<int>& index, SC& mv )
647  { mv.SetBlock(A, index); }
648 
649  static void Assign( const SC& A, SC& mv )
650  { mv.Assign(A); }
651 
652  static void MvRandom( SC& mv )
653  { mv.MvRandom(); }
654 
655  static void MvInit( SC& mv, ScalarType alpha = Teuchos::ScalarTraits<ScalarType>::zero() )
656  { mv.MvInit(alpha); }
657 
658  static void MvPrint( const SC& mv, std::ostream& os )
659  { mv.MvPrint(os); }
660 };
661 
662 } // end namespace Belos
663 #endif
664 
665 #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.
static void MvRandom(MV &mv)
Replace the vectors in mv with random vectors.
static void MvTimesMatAddMv(const ScalarType alpha, const MV &A, const Teuchos::SerialDenseMatrix< int, ScalarType > &B, const ScalarType beta, MV &mv)
Update mv with .
static void MvDot(const MV &mv, const MV &A, std::vector< ScalarType > &b)
Compute a vector b where the components are the individual dot-products of the i-th columns of A and ...
static void Assign(const MV &A, MV &mv)
mv := A
Namespace Anasazi contains the classes, structs, enums and utilities used by the Anasazi package...
Traits class which defines basic operations on multivectors.
static Teuchos::RCP< MV > CloneCopy(const MV &mv)
Creates a new MV and copies contents of mv into the new vector (deep copy).
static void MvTransMv(const ScalarType alpha, const MV &A, const MV &B, Teuchos::SerialDenseMatrix< int, ScalarType > &C)
Compute C := alpha * A^H B.
static void MvScale(MV &mv, const ScalarType alpha)
Scale each element of the vectors in mv with alpha.
Anasazi header file which uses auto-configuration information to include necessary C++ headers...
static Teuchos::RCP< MV > CloneViewNonConst(MV &mv, const std::vector< int > &index)
Creates a new MV that shares the selected contents of mv (shallow copy).
static void SetBlock(const MV &A, const std::vector< int > &index, MV &mv)
Copy the vectors in A to a set of vectors in mv indicated by the indices given in index...
static int GetNumberVecs(const MV &mv)
Obtain the number of vectors in mv.
static void MvAddMv(const ScalarType alpha, const MV &A, const ScalarType beta, const MV &B, MV &mv)
Replace mv with .
static Teuchos::RCP< MV > Clone(const MV &mv, const int numvecs)
Creates a new empty MV containing numvecs columns.