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 with Sandia Corporation,
8 // the U.S. Government retains certain rights in this software.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ***********************************************************************
40 // @HEADER
41 
49 #ifndef ANASAZI_SADDLE_CONTAINER_HPP
50 #define ANASAZI_SADDLE_CONTAINER_HPP
51 
52 #include "AnasaziConfigDefs.hpp"
53 #include "Teuchos_VerboseObject.hpp"
54 
55 #ifdef HAVE_ANASAZI_BELOS
56 #include "BelosConfigDefs.hpp"
57 #include "BelosMultiVecTraits.hpp"
58 #endif
59 
60 using Teuchos::RCP;
61 using Teuchos::rcp;
62 
63 namespace Anasazi {
64 namespace Experimental {
65 
66 template <class ScalarType, class MV>
67 class SaddleContainer //: public Anasazi::SaddleContainer<ScalarType, MV>
68 {
69  template <class Scalar_, class MV_, class OP_> friend class SaddleOperator;
70 
71 private:
73  typedef Teuchos::SerialDenseMatrix<int,ScalarType> SerialDenseMatrix;
74  const ScalarType ONE, ZERO;
75  RCP<MV> upper_;
76  RCP<SerialDenseMatrix> lowerRaw_;
77  std::vector<int> indices_;
78 
79  RCP<SerialDenseMatrix> getLower() const;
80  void setLower(const RCP<SerialDenseMatrix> lower);
81 
82 public:
83  // Constructors
84  SaddleContainer( ) : ONE(Teuchos::ScalarTraits<ScalarType>::one()), ZERO(Teuchos::ScalarTraits<ScalarType>::zero()) { };
85  SaddleContainer( const RCP<MV> X, bool eye=false );
86 
87  // Things that are necessary for compilation
88  // Returns a clone of the current vector
89  SaddleContainer<ScalarType, MV> * Clone(const int nvecs) const;
90  // Returns a duplicate of the current vector
91  SaddleContainer<ScalarType, MV> * CloneCopy() const;
92  // Returns a duplicate of the specified vectors
93  SaddleContainer<ScalarType, MV> * CloneCopy(const std::vector<int> &index) const;
94  SaddleContainer<ScalarType, MV> * CloneViewNonConst(const std::vector<int> &index) const;
95  SaddleContainer<ScalarType, MV> * CloneViewNonConst(const Teuchos::Range1D& index) const;
96  const SaddleContainer<ScalarType, MV> * CloneView(const std::vector<int> &index) const;
97  const SaddleContainer<ScalarType, MV> * CloneView(const Teuchos::Range1D& index) const;
98  ptrdiff_t GetGlobalLength() const { return MVT::GetGlobalLength(*upper_) + lowerRaw_->numRows(); };
99  int GetNumberVecs() const { return MVT::GetNumberVecs(*upper_); };
100  // Update *this with alpha * A * B + beta * (*this)
101  void MvTimesMatAddMv(ScalarType alpha, const SaddleContainer<ScalarType,MV> &A,
102  const Teuchos::SerialDenseMatrix<int, ScalarType> &B,
103  ScalarType beta);
104  // Replace *this with alpha * A + beta * B
105  void MvAddMv(ScalarType alpha, const SaddleContainer<ScalarType,MV>& A,
106  ScalarType beta, const SaddleContainer<ScalarType,MV>& B);
107  // Scale the vectors by alpha
108  void MvScale( ScalarType alpha );
109  // Scale the i-th vector by alpha[i]
110  void MvScale( const std::vector<ScalarType>& alpha );
111  // Compute a dense matrix B through the matrix-matrix multiply alpha * A^H * (*this)
112  void MvTransMv (ScalarType alpha, const SaddleContainer<ScalarType, MV>& A,
113  Teuchos::SerialDenseMatrix< int, ScalarType >& B) const;
114  // 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.
115  void MvDot (const SaddleContainer<ScalarType, MV>& A, std::vector<ScalarType> &b) const;
116  // Compute the 2-norm of each individual vector
117  void MvNorm ( std::vector<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType> &normvec) const;
118  // Copy the vectors in A to a set of vectors in *this. The numvecs vectors in
119  // A are copied to a subset of vectors in *this indicated by the indices given
120  // in index.
121  void SetBlock (const SaddleContainer<ScalarType, MV>& A, const std::vector<int> &index);
122  // Deep copy.
123  void Assign (const SaddleContainer<ScalarType, MV>&A);
124  // Fill the vectors in *this with random numbers.
125  void MvRandom ();
126  // Replace each element of the vectors in *this with alpha.
127  void MvInit (ScalarType alpha);
128  // Prints the multivector to an output stream
129  void MvPrint (std::ostream &os) const;
130 };
131 
132 
133 
134 // THIS IS NEW!
135 template <class ScalarType, class MV>
136 RCP< Teuchos::SerialDenseMatrix<int,ScalarType> > SaddleContainer<ScalarType, MV>::getLower() const
137 {
138  if(indices_.empty())
139  {
140  return lowerRaw_;
141  }
142 
143  int nrows = lowerRaw_->numRows();
144  int ncols = indices_.size();
145 
146  RCP<SerialDenseMatrix> lower = rcp(new SerialDenseMatrix(nrows,ncols,false));
147 
148  for(int r=0; r<nrows; r++)
149  {
150  for(int c=0; c<ncols; c++)
151  {
152  (*lower)(r,c) = (*lowerRaw_)(r,indices_[c]);
153  }
154  }
155 
156  return lower;
157 }
158 
159 
160 
161 // THIS IS NEW!
162 template <class ScalarType, class MV>
163 void SaddleContainer<ScalarType, MV>::setLower(const RCP<SerialDenseMatrix> lower)
164 {
165  // If the indices are empty, lower points to lowerRaw
166  if(indices_.empty())
167  {
168  return;
169  }
170 
171  int nrows = lowerRaw_->numRows();
172  int ncols = indices_.size();
173 
174  for(int r=0; r<nrows; r++)
175  {
176  for(int c=0; c<ncols; c++)
177  {
178  (*lowerRaw_)(r,indices_[c]) = (*lower)(r,c);
179  }
180  }
181 }
182 
183 
184 
185 // Constructor
186 template <class ScalarType, class MV>
187 SaddleContainer<ScalarType, MV>::SaddleContainer( const RCP<MV> X, bool eye )
188 : ONE(Teuchos::ScalarTraits<ScalarType>::one()), ZERO(Teuchos::ScalarTraits<ScalarType>::zero())
189 {
190  int nvecs = MVT::GetNumberVecs(*X);
191 
192  if(eye)
193  {
194  // Initialize upper_ as all 0s
195  upper_ = MVT::Clone(*X, nvecs);
196  MVT::MvInit(*upper_);
197 
198  // Initialize Y to be I
199  lowerRaw_ = rcp(new SerialDenseMatrix(nvecs,nvecs));
200  for(int i=0; i < nvecs; i++)
201  (*lowerRaw_)(i,i) = ONE;
202  }
203  else
204  {
205  // Point upper_ to X
206  upper_ = X;
207 
208  // Initialize Y to be 0
209  lowerRaw_ = rcp(new SerialDenseMatrix(nvecs,nvecs));
210  }
211 }
212 
213 
214 
215 // Returns a clone of the current vector
216 template <class ScalarType, class MV>
217 SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::Clone(const int nvecs) const
218 {
219  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
220 
221  newSC->upper_ = MVT::Clone(*upper_,nvecs);
222  newSC->lowerRaw_ = rcp(new SerialDenseMatrix(lowerRaw_->numRows(),nvecs));
223 
224  return newSC;
225 }
226 
227 
228 
229 // Returns a duplicate of the current vector
230 template <class ScalarType, class MV>
231 SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::CloneCopy() const
232 {
233  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
234 
235  newSC->upper_ = MVT::CloneCopy(*upper_);
236  newSC->lowerRaw_ = getLower();
237 
238  return newSC;
239 }
240 
241 
242 
243 // Returns a duplicate of the specified vectors
244 template <class ScalarType, class MV>
245 SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::CloneCopy(const std::vector< int > &index) const
246 {
247  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
248 
249  newSC->upper_ = MVT::CloneCopy(*upper_,index);
250 
251  int ncols = index.size();
252  int nrows = lowerRaw_->numRows();
253  RCP<SerialDenseMatrix> lower = getLower();
254  newSC->lowerRaw_ = rcp(new SerialDenseMatrix(nrows,ncols));
255  for(int c=0; c < ncols; c++)
256  {
257  for(int r=0; r < nrows; r++)
258  (*newSC->lowerRaw_)(r,c) = (*lower)(r,index[c]);
259  }
260 
261  return newSC;
262 }
263 
264 
265 
266 // THIS IS NEW!
267 template <class ScalarType, class MV>
268 SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::CloneViewNonConst(const std::vector<int> &index) const
269 {
270  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
271 
272  newSC->upper_ = MVT::CloneViewNonConst(*upper_,index);
273 
274  newSC->lowerRaw_ = lowerRaw_;
275 
276  if(!indices_.empty())
277  {
278  newSC->indices_.resize(index.size());
279  for(unsigned int i=0; i<index.size(); i++)
280  {
281  newSC->indices_[i] = indices_[index[i]];
282  }
283  }
284  else
285  {
286  newSC->indices_ = index;
287  }
288 
289  return newSC;
290 }
291 
292 
293 template <class ScalarType, class MV>
294 SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::CloneViewNonConst(const Teuchos::Range1D& index) const
295 {
296  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
297 
298  newSC->upper_ = MVT::CloneViewNonConst(*upper_,index);
299 
300  newSC->lowerRaw_ = lowerRaw_;
301 
302  newSC->indices_.resize(index.size());
303  for(unsigned int i=0; i<index.size(); i++)
304  {
305  newSC->indices_[i] = indices_[index.lbound()+i];
306  }
307 
308  return newSC;
309 }
310 
311 
312 
313 // THIS IS NEW!
314 template <class ScalarType, class MV>
315 const SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::CloneView(const std::vector<int> &index) const
316 {
317  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
318 
319  newSC->upper_ = MVT::CloneViewNonConst(*upper_,index);
320 
321  newSC->lowerRaw_ = lowerRaw_;
322 
323  if(!indices_.empty())
324  {
325  newSC->indices_.resize(index.size());
326  for(unsigned int i=0; i<index.size(); i++)
327  {
328  newSC->indices_[i] = indices_[index[i]];
329  }
330  }
331  else
332  {
333  newSC->indices_ = index;
334  }
335 
336  return newSC;
337 }
338 
339 
340 template <class ScalarType, class MV>
341 const SaddleContainer<ScalarType, MV> * SaddleContainer<ScalarType, MV>::CloneView(const Teuchos::Range1D& index) const
342 {
343  SaddleContainer<ScalarType, MV> * newSC = new SaddleContainer<ScalarType, MV>();
344 
345  newSC->upper_ = MVT::CloneViewNonConst(*upper_,index);
346 
347  newSC->lowerRaw_ = lowerRaw_;
348 
349  newSC->indices_.resize(index.size());
350  for(unsigned int i=0; i<index.size(); i++)
351  {
352  newSC->indices_[i] = indices_[index.lbound()+i];
353  }
354 
355  return newSC;
356 }
357 
358 
359 
360 // Update *this with alpha * A * B + beta * (*this)
361 // THIS IS NEW!
362 template <class ScalarType, class MV>
363 void SaddleContainer<ScalarType, MV>::MvTimesMatAddMv(ScalarType alpha, const SaddleContainer<ScalarType,MV> &A,
364  const Teuchos::SerialDenseMatrix<int, ScalarType> &B,
365  ScalarType beta)
366 {
367  MVT::MvTimesMatAddMv(alpha,*(A.upper_),B,beta,*upper_);
368  RCP<SerialDenseMatrix> lower = getLower();
369  RCP<SerialDenseMatrix> Alower = A.getLower();
370  lower->multiply(Teuchos::NO_TRANS,Teuchos::NO_TRANS,alpha,*Alower,B,beta);
371  setLower(lower);
372 }
373 
374 
375 
376 // Replace *this with alpha * A + beta * B
377 template <class ScalarType, class MV>
378 void SaddleContainer<ScalarType, MV>::MvAddMv(ScalarType alpha, const SaddleContainer<ScalarType,MV>& A,
379  ScalarType beta, const SaddleContainer<ScalarType,MV>& B)
380 {
381  MVT::MvAddMv(alpha, *(A.upper_), beta, *(B.upper_), *upper_);
382 
383  RCP<SerialDenseMatrix> lower = getLower();
384  RCP<SerialDenseMatrix> Alower = A.getLower();
385  RCP<SerialDenseMatrix> Blower = B.getLower();
386 
387  //int ncolsA = Alower->numCols(); // unused
388  //int ncolsThis = lower->numCols(); // unused
389  //int nrows = lower->numRows(); // unused
390 
391  // Y = alpha A
392  lower->assign(*Alower);
393  if(alpha != ONE)
394  lower->scale(alpha);
395  // Y += beta B
396  if(beta == ONE)
397  *lower += *Blower;
398  else if(beta == -ONE)
399  *lower -= *Blower;
400  else if(beta != ZERO)
401  {
402  SerialDenseMatrix scaledB(*Blower);
403  scaledB.scale(beta);
404  *lower += *Blower;
405  }
406 
407  setLower(lower);
408 }
409 
410 
411 
412 // Scale the vectors by alpha
413 template <class ScalarType, class MV>
414 void SaddleContainer<ScalarType, MV>::MvScale( ScalarType alpha )
415 {
416  MVT::MvScale(*upper_, alpha);
417 
418 
419  RCP<SerialDenseMatrix> lower = getLower();
420  lower->scale(alpha);
421  setLower(lower);
422 }
423 
424 
425 
426 // Scale the i-th vector by alpha[i]
427 template <class ScalarType, class MV>
428 void SaddleContainer<ScalarType, MV>::MvScale( const std::vector<ScalarType>& alpha )
429 {
430  MVT::MvScale(*upper_, alpha);
431 
432  RCP<SerialDenseMatrix> lower = getLower();
433 
434  int nrows = lower->numRows();
435  int ncols = lower->numCols();
436 
437  for(int c=0; c<ncols; c++)
438  {
439  for(int r=0; r<nrows; r++)
440  (*lower)(r,c) *= alpha[c];
441  }
442 
443  setLower(lower);
444 }
445 
446 
447 
448 // Compute a dense matrix B through the matrix-matrix multiply alpha * A^H * (*this)
449 // THIS IS NEW!
450 template <class ScalarType, class MV>
451 void SaddleContainer<ScalarType, MV>::MvTransMv (ScalarType alpha, const SaddleContainer<ScalarType, MV>& A,
452  Teuchos::SerialDenseMatrix< int, ScalarType >& B) const
453 {
454  MVT::MvTransMv(alpha,*(A.upper_),*upper_,B);
455  RCP<SerialDenseMatrix> lower = getLower();
456  RCP<SerialDenseMatrix> Alower = A.getLower();
457  B.multiply(Teuchos::TRANS,Teuchos::NO_TRANS,alpha,*(Alower),*lower,ONE);
458 }
459 
460 
461 
462 // 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.
463 template <class ScalarType, class MV>
464 void SaddleContainer<ScalarType, MV>::MvDot (const SaddleContainer<ScalarType, MV>& A, std::vector<ScalarType> &b) const
465 {
466  MVT::MvDot(*upper_, *(A.upper_), b);
467 
468  RCP<SerialDenseMatrix> lower = getLower();
469  RCP<SerialDenseMatrix> Alower = A.getLower();
470 
471  int nrows = lower->numRows();
472  int ncols = lower->numCols();
473 
474  for(int c=0; c < ncols; c++)
475  {
476  for(int r=0; r < nrows; r++)
477  {
478  b[c] += ((*Alower)(r,c) * (*lower)(r,c));
479  }
480  }
481 }
482 
483 
484 
485 // Compute the 2-norm of each individual vector
486 // THIS IS NEW!
487 template <class ScalarType, class MV>
488 void SaddleContainer<ScalarType, MV>::MvNorm ( std::vector<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType> &normvec) const
489 {
490  // TODO: Make this better
491  MvDot(*this,normvec);
492  for(unsigned int i=0; i<normvec.size(); i++)
493  normvec[i] = sqrt(normvec[i]);
494 }
495 
496 
497 
498 // Copy the vectors in A to a set of vectors in *this. The numvecs vectors in
499 // A are copied to a subset of vectors in *this indicated by the indices given
500 // in index.
501 template <class ScalarType, class MV>
502 void SaddleContainer<ScalarType, MV>::SetBlock (const SaddleContainer<ScalarType, MV>& A, const std::vector<int> &index)
503 {
504  MVT::SetBlock(*(A.upper_), index, *upper_);
505 
506  RCP<SerialDenseMatrix> lower = getLower();
507  RCP<SerialDenseMatrix> Alower = A.getLower();
508 
509  int nrows = lower->numRows();
510 
511  int nvecs = index.size();
512  for(int c=0; c<nvecs; c++)
513  {
514  for(int r=0; r<nrows; r++)
515  (*lower)(r,index[c]) = (*Alower)(r,c);
516  }
517 
518  setLower(lower);
519 }
520 
521 
522 
523 // Deep copy.
524 template <class ScalarType, class MV>
525 void SaddleContainer<ScalarType, MV>::Assign (const SaddleContainer<ScalarType, MV>&A)
526 {
527  MVT::Assign(*(A.upper_),*(upper_));
528 
529  RCP<SerialDenseMatrix> lower = getLower();
530  RCP<SerialDenseMatrix> Alower = A.getLower();
531 
532  *lower = *Alower; // This is a well-defined operator for SerialDenseMatrix
533 
534  setLower(lower);
535 }
536 
537 
538 
539 // Fill the vectors in *this with random numbers.
540 // THIS IS NEW!
541 template <class ScalarType, class MV>
542 void SaddleContainer<ScalarType, MV>::MvRandom ()
543 {
544  MVT::MvRandom(*upper_);
545 
546  RCP<SerialDenseMatrix> lower = getLower();
547  lower->random();
548  setLower(lower);
549 }
550 
551 
552 
553 // Replace each element of the vectors in *this with alpha.
554 template <class ScalarType, class MV>
555 void SaddleContainer<ScalarType, MV>::MvInit (ScalarType alpha)
556 {
557  MVT::MvInit(*upper_,alpha);
558 
559  RCP<SerialDenseMatrix> lower = getLower();
560  lower->putScalar(alpha);
561  setLower(lower);
562 }
563 
564 
565 
566 // Prints the multivector to an output stream
567 template <class ScalarType, class MV>
568 void SaddleContainer<ScalarType, MV>::MvPrint (std::ostream &os) const
569 {
570  RCP<SerialDenseMatrix> lower = getLower();
571  //int nrows = lower->numRows(); // unused
572  //int ncols = lower->numCols(); // unused
573 
574  os << "Object SaddleContainer" << std::endl;
575  os << "X\n";
576  upper_->describe(*(Teuchos::VerboseObjectBase::getDefaultOStream()),Teuchos::VERB_EXTREME);
577 // os << "X\n" << *upper_ << std::endl;
578 
579  os << "Y\n" << *lower << std::endl;
580 }
581 
582 } // End namespace Experimental
583 
584 template<class ScalarType, class MV >
585 class MultiVecTraits<ScalarType,Experimental::SaddleContainer<ScalarType, MV> >
586 {
587 typedef Experimental::SaddleContainer<ScalarType,MV> SC;
588 
589 public:
590  static RCP<SC > Clone( const SC& mv, const int numvecs )
591  { return rcp( const_cast<SC&>(mv).Clone(numvecs) ); }
592 
593  static RCP<SC > CloneCopy( const SC& mv )
594  { return rcp( const_cast<SC&>(mv).CloneCopy() ); }
595 
596  static RCP<SC > CloneCopy( const SC& mv, const std::vector<int>& index )
597  { return rcp( const_cast<SC&>(mv).CloneCopy(index) ); }
598 
599  static ptrdiff_t GetGlobalLength( const SC& mv )
600  { return mv.GetGlobalLength(); }
601 
602  static int GetNumberVecs( const SC& mv )
603  { return mv.GetNumberVecs(); }
604 
605  static void MvTimesMatAddMv( ScalarType alpha, const SC& A,
606  const Teuchos::SerialDenseMatrix<int,ScalarType>& B,
607  ScalarType beta, SC& mv )
608  { mv.MvTimesMatAddMv(alpha, A, B, beta); }
609 
610  static void MvAddMv( ScalarType alpha, const SC& A, ScalarType beta, const SC& B, SC& mv )
611  { mv.MvAddMv(alpha, A, beta, B); }
612 
613  static void MvTransMv( ScalarType alpha, const SC& A, const SC& mv, Teuchos::SerialDenseMatrix<int,ScalarType>& B)
614  { mv.MvTransMv(alpha, A, B); }
615 
616  static void MvDot( const SC& mv, const SC& A, std::vector<ScalarType> & b)
617  { mv.MvDot( A, b); }
618 
619  static void MvScale ( SC& mv, ScalarType alpha )
620  { mv.MvScale( alpha ); }
621 
622  static void MvScale ( SC& mv, const std::vector<ScalarType>& alpha )
623  { mv.MvScale( alpha ); }
624 
625  static void MvNorm( const SC& mv, std::vector<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType> & normvec)
626  { mv.MvNorm(normvec); }
627 
628  static void SetBlock( const SC& A, const std::vector<int>& index, SC& mv )
629  { mv.SetBlock(A, index); }
630 
631  static void Assign( const SC& A, SC& mv )
632  { mv.Assign(A); }
633 
634  static void MvRandom( SC& mv )
635  { mv.MvRandom(); }
636 
637  static void MvInit( SC& mv, ScalarType alpha = Teuchos::ScalarTraits<ScalarType>::zero() )
638  { mv.MvInit(alpha); }
639 
640  static void MvPrint( const SC& mv, std::ostream& os )
641  { mv.MvPrint(os); }
642 };
643 
644 } // end namespace Anasazi
645 
646 #ifdef HAVE_ANASAZI_BELOS
647 namespace Belos
648 {
649 
650 template<class ScalarType, class MV >
651 class MultiVecTraits< ScalarType, Anasazi::Experimental::SaddleContainer<ScalarType,MV> >
652 {
653 typedef Anasazi::Experimental::SaddleContainer<ScalarType,MV> SC;
654 public:
655  static RCP<SC > Clone( const SC& mv, const int numvecs )
656  { return rcp( const_cast<SC&>(mv).Clone(numvecs) ); }
657 
658  static RCP<SC > CloneCopy( const SC& mv )
659  { return rcp( const_cast<SC&>(mv).CloneCopy() ); }
660 
661  static RCP<SC > CloneCopy( const SC& mv, const std::vector<int>& index )
662  { return rcp( const_cast<SC&>(mv).CloneCopy(index) ); }
663 
664  static RCP<SC> CloneViewNonConst( SC& mv, const std::vector<int>& index )
665  { return rcp( mv.CloneViewNonConst(index) ); }
666 
667  static RCP<SC> CloneViewNonConst( SC& mv, const Teuchos::Range1D& index )
668  { return rcp( mv.CloneViewNonConst(index) ); }
669 
670  static RCP<const SC> CloneView( const SC& mv, const std::vector<int> & index )
671  { return rcp( mv.CloneView(index) ); }
672 
673  static RCP<const SC> CloneView( const SC& mv, const Teuchos::Range1D& index )
674  { return rcp( mv.CloneView(index) ); }
675 
676  static ptrdiff_t GetGlobalLength( const SC& mv )
677  { return mv.GetGlobalLength(); }
678 
679  static int GetNumberVecs( const SC& mv )
680  { return mv.GetNumberVecs(); }
681 
682  static void MvTimesMatAddMv( ScalarType alpha, const SC& A,
683  const Teuchos::SerialDenseMatrix<int,ScalarType>& B,
684  ScalarType beta, SC& mv )
685  { mv.MvTimesMatAddMv(alpha, A, B, beta); }
686 
687  static void MvAddMv( ScalarType alpha, const SC& A, ScalarType beta, const SC& B, SC& mv )
688  { mv.MvAddMv(alpha, A, beta, B); }
689 
690  static void MvTransMv( ScalarType alpha, const SC& A, const SC& mv, Teuchos::SerialDenseMatrix<int,ScalarType>& B)
691  { mv.MvTransMv(alpha, A, B); }
692 
693  static void MvDot( const SC& mv, const SC& A, std::vector<ScalarType> & b)
694  { mv.MvDot( A, b); }
695 
696  static void MvScale ( SC& mv, ScalarType alpha )
697  { mv.MvScale( alpha ); }
698 
699  static void MvScale ( SC& mv, const std::vector<ScalarType>& alpha )
700  { mv.MvScale( alpha ); }
701 
702  // TODO: MAKE SURE TYPE == TWONORM
703  static void MvNorm( const SC& mv, std::vector<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType> & normvec, NormType type=TwoNorm)
704  { mv.MvNorm(normvec); }
705 
706  static void SetBlock( const SC& A, const std::vector<int>& index, SC& mv )
707  { mv.SetBlock(A, index); }
708 
709  static void Assign( const SC& A, SC& mv )
710  { mv.Assign(A); }
711 
712  static void MvRandom( SC& mv )
713  { mv.MvRandom(); }
714 
715  static void MvInit( SC& mv, ScalarType alpha = Teuchos::ScalarTraits<ScalarType>::zero() )
716  { mv.MvInit(alpha); }
717 
718  static void MvPrint( const SC& mv, std::ostream& os )
719  { mv.MvPrint(os); }
720 
721 #ifdef HAVE_BELOS_TSQR
722  typedef Belos::details::StubTsqrAdapter<SC> tsqr_adaptor_type;
735 #endif // HAVE_BELOS_TSQR
736 };
737 
738 } // end namespace Belos
739 #endif
740 
741 #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.
Namespace Anasazi contains the classes, structs, enums and utilities used by the Anasazi package...
Traits class which defines basic operations on multivectors.
Anasazi header file which uses auto-configuration information to include necessary C++ headers...
static int GetNumberVecs(const MV &mv)
Obtain the number of vectors in mv.
static Teuchos::RCP< MV > Clone(const MV &mv, const int numvecs)
Creates a new empty MV containing numvecs columns.