MueLu  Version of the Day
BelosXpetraAdapterMultiVector.hpp
Go to the documentation of this file.
1 // @HEADER
2 //
3 // ***********************************************************************
4 //
5 // MueLu: A package for multigrid based preconditioning
6 // Copyright 2012 Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact
39 // Jonathan Hu (jhu@sandia.gov)
40 // Andrey Prokopenko (aprokop@sandia.gov)
41 // Ray Tuminaro (rstumin@sandia.gov)
42 //
43 // ***********************************************************************
44 //
45 // @HEADER
46 #ifndef BELOS_XPETRA_ADAPTER_MULTIVECTOR_HPP
47 #define BELOS_XPETRA_ADAPTER_MULTIVECTOR_HPP
48 
49 #include <Xpetra_ConfigDefs.hpp>
50 
51 #include <Xpetra_MultiVector.hpp>
52 
53 #ifdef HAVE_XPETRA_EPETRA
54 #include <Xpetra_EpetraMultiVector.hpp>
55 #endif
56 
57 #ifdef HAVE_XPETRA_TPETRA
58 #include <Xpetra_TpetraMultiVector.hpp>
59 #endif
60 
61 #include <BelosConfigDefs.hpp>
62 #include <BelosTypes.hpp>
63 #include <BelosMultiVecTraits.hpp>
64 #include <BelosOperatorTraits.hpp>
65 
66 #ifdef HAVE_XPETRA_EPETRA
67 #include <BelosEpetraAdapter.hpp>
68 #endif
69 
70 #ifdef HAVE_XPETRA_TPETRA
71 #include <BelosTpetraAdapter.hpp>
72 #endif
73 
74 namespace Belos { // should be moved to Belos or Xpetra?
75 
76  using Teuchos::RCP;
77  using Teuchos::rcp;
78 
80  //
81  // Implementation of the Belos::MultiVecTraits for Xpetra::MultiVector.
82  //
84 
89  template<class Scalar, class LO, class GO, class Node>
90  class MultiVecTraits<Scalar, Xpetra::MultiVector<Scalar,LO,GO,Node> >
91  {
92 
93  private:
94 #ifdef HAVE_XPETRA_TPETRA
95  typedef Xpetra::TpetraMultiVector<Scalar,LO,GO,Node> TpetraMultiVector;
96  typedef MultiVecTraits <Scalar,Tpetra::MultiVector<Scalar,LO,GO,Node> > MultiVecTraitsTpetra;
97 #endif
98 
99  public:
100 
101 #ifdef HAVE_BELOS_XPETRA_TIMERS
102  static RCP<Teuchos::Time> mvTimesMatAddMvTimer_, mvTransMvTimer_;
103 #endif
104 
105  static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> > Clone( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, const int numvecs )
106  {
107 
108 #ifdef HAVE_XPETRA_TPETRA
109  if (mv.getMap()->lib() == Xpetra::UseTpetra)
110  return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::Clone(toTpetra(mv), numvecs)));
111 #endif
112 
113  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
114  XPETRA_FACTORY_END;
115  }
116 
117  static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> > CloneCopy( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
118  {
119 
120 #ifdef HAVE_XPETRA_TPETRA
121  if (mv.getMap()->lib() == Xpetra::UseTpetra)
122  return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::CloneCopy(toTpetra(mv))));
123 #endif
124 
125  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
126  XPETRA_FACTORY_END;
127  }
128 
129  static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> > CloneCopy( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, const std::vector<int>& index )
130  {
131 
132 #ifdef HAVE_XPETRA_TPETRA
133  if (mv.getMap()->lib() == Xpetra::UseTpetra)
134  return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::CloneCopy(toTpetra(mv), index)));
135 #endif
136 
137  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
138  XPETRA_FACTORY_END;
139  }
140 
141  static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
142  CloneCopy (const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
143  const Teuchos::Range1D& index)
144  {
145 
146 #ifdef HAVE_XPETRA_TPETRA
147  if (mv.getMap()->lib() == Xpetra::UseTpetra)
148  return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::CloneCopy(toTpetra(mv), index)));
149 #endif
150 
151  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
152  XPETRA_FACTORY_END;
153  }
154 
155  static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> > CloneViewNonConst( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, const std::vector<int>& index )
156  {
157 
158 #ifdef HAVE_XPETRA_TPETRA
159  if (mv.getMap()->lib() == Xpetra::UseTpetra)
160  return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::CloneViewNonConst(toTpetra(mv), index)));
161 #endif
162 
163  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
164  XPETRA_FACTORY_END;
165  }
166 
167  static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
168  CloneViewNonConst(Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
169  const Teuchos::Range1D& index)
170  {
171 
172 #ifdef HAVE_XPETRA_TPETRA
173  if (mv.getMap()->lib() == Xpetra::UseTpetra)
174  return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::CloneViewNonConst(toTpetra(mv), index)));
175 #endif
176 
177  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
178  XPETRA_FACTORY_END;
179  }
180 
181  static RCP<const Xpetra::MultiVector<Scalar,LO,GO,Node> > CloneView(const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, const std::vector<int>& index )
182  {
183 
184 #ifdef HAVE_XPETRA_TPETRA
185  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
186  //TODO: double check if the const_cast is safe here.
187  RCP<const Tpetra::MultiVector<Scalar,LO,GO,Node> > r = MultiVecTraitsTpetra::CloneView(toTpetra(mv), index);
188  return rcp(new TpetraMultiVector(Teuchos::rcp_const_cast<Tpetra::MultiVector<Scalar,LO,GO,Node> >(r)));
189  }
190 #endif
191 
192  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
193  XPETRA_FACTORY_END;
194  }
195 
196  static RCP<const Xpetra::MultiVector<Scalar,LO,GO,Node> >
197  CloneView (const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
198  const Teuchos::Range1D& index)
199  {
200 
201 #ifdef HAVE_XPETRA_TPETRA
202  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
203  //TODO: double check if the const_cast is safe here.
204  RCP<const Tpetra::MultiVector<Scalar,LO,GO,Node> > r = MultiVecTraitsTpetra::CloneView(toTpetra(mv), index);
205  return rcp(new TpetraMultiVector(Teuchos::rcp_const_cast<Tpetra::MultiVector<Scalar,LO,GO,Node> >(r)));
206  }
207 #endif
208 
209  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
210  XPETRA_FACTORY_END;
211  }
212 
213  static ptrdiff_t GetGlobalLength( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
214  {
215 
216 #ifdef HAVE_XPETRA_TPETRA
217  if (mv.getMap()->lib() == Xpetra::UseTpetra)
218  return MultiVecTraitsTpetra::GetGlobalLength(toTpetra(mv));
219 #endif
220 
221  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
222  XPETRA_FACTORY_END;
223  }
224 
225  static int GetNumberVecs( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
226  {
227 
228 #ifdef HAVE_XPETRA_TPETRA
229  if (mv.getMap()->lib() == Xpetra::UseTpetra)
230  return MultiVecTraitsTpetra::GetNumberVecs(toTpetra(mv));
231 #endif
232 
233  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
234  XPETRA_FACTORY_END;
235  }
236 
237  static bool HasConstantStride( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
238  {
239 
240 #ifdef HAVE_XPETRA_TPETRA
241  if (mv.getMap()->lib() == Xpetra::UseTpetra)
242  return MultiVecTraitsTpetra::HasConstantStride(toTpetra(mv));
243 #endif
244 
245  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
246  XPETRA_FACTORY_END;
247  }
248 
249  static void MvTimesMatAddMv( Scalar alpha, const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
250  const Teuchos::SerialDenseMatrix<int,Scalar>& B,
251  Scalar beta, Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
252  {
253 
254 #ifdef HAVE_BELOS_XPETRA_TIMERS
255  Teuchos::TimeMonitor lcltimer(*mvTimesMatAddMvTimer_);
256 #endif
257 
258  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
259  XPETRA_FACTORY_END;
260 
261 
262 #ifdef HAVE_XPETRA_TPETRA
263  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
264  MultiVecTraitsTpetra::MvTimesMatAddMv(alpha, toTpetra(A), B, beta, toTpetra(mv));
265  return;
266  }
267 #endif
268 
269  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
270  XPETRA_FACTORY_END;
271  }
272 
273  static void MvAddMv( Scalar alpha, const Xpetra::MultiVector<Scalar,LO,GO,Node>& A, Scalar beta, const Xpetra::MultiVector<Scalar,LO,GO,Node>& B, Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
274  {
275 
276 #ifdef HAVE_XPETRA_TPETRA
277  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
278  MultiVecTraitsTpetra::MvAddMv(alpha, toTpetra(A), beta, toTpetra(B), toTpetra(mv));
279  return;
280  }
281 
282 #endif
283 
284  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
285  XPETRA_FACTORY_END;
286  }
287 
288  static void MvScale ( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, Scalar alpha )
289  {
290 
291 #ifdef HAVE_XPETRA_TPETRA
292  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
293  MultiVecTraitsTpetra::MvScale(toTpetra(mv), alpha);
294  return;
295  }
296 #endif
297 
298  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
299  XPETRA_FACTORY_END;
300  }
301 
302  static void MvScale ( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, const std::vector<Scalar>& alphas )
303  {
304 
305 #ifdef HAVE_XPETRA_TPETRA
306  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
307  MultiVecTraitsTpetra::MvScale(toTpetra(mv), alphas);
308  return;
309  }
310 #endif
311 
312  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
313  XPETRA_FACTORY_END;
314  }
315 
316  static void MvTransMv( Scalar alpha, const Xpetra::MultiVector<Scalar,LO,GO,Node>& A, const Xpetra::MultiVector<Scalar,LO,GO,Node>& B, Teuchos::SerialDenseMatrix<int,Scalar>& C)
317  {
318 
319 #ifdef HAVE_BELOS_XPETRA_TIMERS
320  Teuchos::TimeMonitor lcltimer(*mvTransMvTimer_);
321 #endif
322 
323 #ifdef HAVE_XPETRA_TPETRA
324  if (A.getMap()->lib() == Xpetra::UseTpetra) {
325  MultiVecTraitsTpetra::MvTransMv(alpha, toTpetra(A), toTpetra(B), C);
326  return;
327  }
328 #endif
329 
330  XPETRA_FACTORY_ERROR_IF_EPETRA(A.getMap()->lib());
331  XPETRA_FACTORY_END;
332  }
333 
334  static void MvDot( const Xpetra::MultiVector<Scalar,LO,GO,Node>& A, const Xpetra::MultiVector<Scalar,LO,GO,Node>& B, std::vector<Scalar> &dots)
335  {
336 
337 #ifdef HAVE_XPETRA_TPETRA
338  if (A.getMap()->lib() == Xpetra::UseTpetra) {
339  MultiVecTraitsTpetra::MvDot(toTpetra(A), toTpetra(B), dots);
340  return;
341  }
342 #endif
343 
344  XPETRA_FACTORY_ERROR_IF_EPETRA(A.getMap()->lib());
345  XPETRA_FACTORY_END;
346  }
347 
348  static void MvNorm(const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, std::vector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType> &normvec, NormType type=TwoNorm)
349  {
350 
351 #ifdef HAVE_XPETRA_TPETRA
352  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
353  MultiVecTraitsTpetra::MvNorm(toTpetra(mv), normvec, type);
354  return;
355  }
356 #endif
357 
358  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
359  XPETRA_FACTORY_END;
360  }
361 
362  static void SetBlock( const Xpetra::MultiVector<Scalar,LO,GO,Node>& A, const std::vector<int>& index, Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
363  {
364 
365 #ifdef HAVE_XPETRA_TPETRA
366  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
367  MultiVecTraitsTpetra::SetBlock(toTpetra(A), index, toTpetra(mv));
368  return;
369  }
370 #endif
371 
372  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
373  XPETRA_FACTORY_END;
374  }
375 
376  static void
377  SetBlock (const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
378  const Teuchos::Range1D& index,
379  Xpetra::MultiVector<Scalar,LO,GO,Node>& mv)
380  {
381 
382 #ifdef HAVE_XPETRA_TPETRA
383  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
384  MultiVecTraitsTpetra::SetBlock(toTpetra(A), index, toTpetra(mv));
385  return;
386  }
387 #endif
388 
389  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
390  XPETRA_FACTORY_END;
391  }
392 
393  static void
394  Assign (const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
395  Xpetra::MultiVector<Scalar,LO,GO,Node>& mv)
396  {
397 
398 #ifdef HAVE_XPETRA_TPETRA
399  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
400  MultiVecTraitsTpetra::Assign(toTpetra(A), toTpetra(mv));
401  return;
402  }
403 #endif
404 
405  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
406  XPETRA_FACTORY_END;
407  }
408 
409  static void MvRandom( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
410  {
411 
412 #ifdef HAVE_XPETRA_TPETRA
413  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
414  MultiVecTraitsTpetra::MvRandom(toTpetra(mv));
415  return;
416  }
417 #endif
418 
419  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
420  XPETRA_FACTORY_END;
421  }
422 
423  static void MvInit( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, Scalar alpha = Teuchos::ScalarTraits<Scalar>::zero() )
424  {
425 
426 #ifdef HAVE_XPETRA_TPETRA
427  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
428  MultiVecTraitsTpetra::MvInit(toTpetra(mv), alpha);
429  return;
430  }
431 #endif
432 
433  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
434  XPETRA_FACTORY_END;
435  }
436 
437  static void MvPrint( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, std::ostream& os )
438  {
439 
440 #ifdef HAVE_XPETRA_TPETRA
441  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
442  MultiVecTraitsTpetra::MvPrint(toTpetra(mv), os);
443  return;
444  }
445 #endif
446 
447  XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
448  XPETRA_FACTORY_END;
449  }
450 
451  };
452 
453 
454  template<>
455  class MultiVecTraits<double, Xpetra::MultiVector<double,int,int,KokkosClassic::DefaultNode::DefaultNodeType> >
456  {
457 
458  private:
459 
460  typedef double Scalar;
461  typedef int LO;
462  typedef int GO;
463  typedef KokkosClassic::DefaultNode::DefaultNodeType Node;
464 
465 #ifdef HAVE_XPETRA_TPETRA
466  typedef Xpetra::TpetraMultiVector<Scalar,LO,GO,Node> TpetraMultiVector;
467  typedef MultiVecTraits <Scalar, Tpetra::MultiVector<Scalar,LO,GO,Node> > MultiVecTraitsTpetra;
468 #endif
469 
470 #ifdef HAVE_XPETRA_EPETRA
471  typedef Xpetra::EpetraMultiVectorT<GO,Node> EpetraMultiVector;
472  typedef MultiVecTraits <Scalar, Epetra_MultiVector> MultiVecTraitsEpetra;
473 #endif
474 
475  public:
476 
477 #ifdef HAVE_BELOS_XPETRA_TIMERS
478  static RCP<Teuchos::Time> mvTimesMatAddMvTimer_, mvTransMvTimer_;
479 #endif
480 
481  static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> > Clone( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, const int numvecs )
482  {
483 
484 #ifdef HAVE_XPETRA_TPETRA
485  if (mv.getMap()->lib() == Xpetra::UseTpetra)
486  return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::Clone(toTpetra(mv), numvecs)));
487 #endif
488 
489 #ifdef HAVE_XPETRA_EPETRA
490  if (mv.getMap()->lib() == Xpetra::UseEpetra)
491  return rcp(new EpetraMultiVector(MultiVecTraitsEpetra::Clone(toEpetra(mv), numvecs)));
492 #endif
493 
494  XPETRA_FACTORY_END;
495  }
496 
497  static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> > CloneCopy( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
498  {
499 
500 #ifdef HAVE_XPETRA_TPETRA
501  if (mv.getMap()->lib() == Xpetra::UseTpetra)
502  return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::CloneCopy(toTpetra(mv))));
503 #endif
504 
505 #ifdef HAVE_XPETRA_EPETRA
506  if (mv.getMap()->lib() == Xpetra::UseEpetra)
507  return rcp(new EpetraMultiVector(MultiVecTraitsEpetra::CloneCopy(toEpetra(mv))));
508 #endif
509 
510  XPETRA_FACTORY_END;
511  }
512 
513  static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> > CloneCopy( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, const std::vector<int>& index )
514  {
515 
516 #ifdef HAVE_XPETRA_TPETRA
517  if (mv.getMap()->lib() == Xpetra::UseTpetra)
518  return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::CloneCopy(toTpetra(mv), index)));
519 #endif
520 
521 #ifdef HAVE_XPETRA_EPETRA
522  if (mv.getMap()->lib() == Xpetra::UseEpetra)
523  return rcp(new EpetraMultiVector(MultiVecTraitsEpetra::CloneCopy(toEpetra(mv), index)));
524 #endif
525 
526  XPETRA_FACTORY_END;
527  }
528 
529  static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
530  CloneCopy (const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
531  const Teuchos::Range1D& index)
532  {
533 
534 #ifdef HAVE_XPETRA_TPETRA
535  if (mv.getMap()->lib() == Xpetra::UseTpetra)
536  return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::CloneCopy(toTpetra(mv), index)));
537 #endif
538 
539 #ifdef HAVE_XPETRA_EPETRA
540  if (mv.getMap()->lib() == Xpetra::UseEpetra)
541  return rcp(new EpetraMultiVector(MultiVecTraitsEpetra::CloneCopy(toEpetra(mv), index)));
542 #endif
543 
544  XPETRA_FACTORY_END;
545  }
546 
547  static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> > CloneViewNonConst( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, const std::vector<int>& index )
548  {
549 
550 #ifdef HAVE_XPETRA_TPETRA
551  if (mv.getMap()->lib() == Xpetra::UseTpetra)
552  return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::CloneViewNonConst(toTpetra(mv), index)));
553 #endif
554 
555 #ifdef HAVE_XPETRA_EPETRA
556  if (mv.getMap()->lib() == Xpetra::UseEpetra)
557  return rcp(new EpetraMultiVector(MultiVecTraitsEpetra::CloneViewNonConst(toEpetra(mv), index)));
558 #endif
559 
560  XPETRA_FACTORY_END;
561  }
562 
563  static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
564  CloneViewNonConst(Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
565  const Teuchos::Range1D& index)
566  {
567 
568 #ifdef HAVE_XPETRA_TPETRA
569  if (mv.getMap()->lib() == Xpetra::UseTpetra)
570  return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::CloneViewNonConst(toTpetra(mv), index)));
571 #endif
572 
573 #ifdef HAVE_XPETRA_EPETRA
574  if (mv.getMap()->lib() == Xpetra::UseEpetra)
575  return rcp(new EpetraMultiVector(MultiVecTraitsEpetra::CloneViewNonConst(toEpetra(mv), index)));
576 #endif
577 
578  XPETRA_FACTORY_END;
579  }
580 
581  static RCP<const Xpetra::MultiVector<Scalar,LO,GO,Node> > CloneView(const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, const std::vector<int>& index )
582  {
583 
584 #ifdef HAVE_XPETRA_TPETRA
585  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
586  //TODO: double check if the const_cast is safe here.
587  RCP<const Tpetra::MultiVector<Scalar,LO,GO,Node> > r = MultiVecTraitsTpetra::CloneView(toTpetra(mv), index);
588  return rcp(new TpetraMultiVector(Teuchos::rcp_const_cast<Tpetra::MultiVector<Scalar,LO,GO,Node> >(r)));
589  }
590 #endif
591 
592 #ifdef HAVE_XPETRA_EPETRA
593  if (mv.getMap()->lib() == Xpetra::UseEpetra) {
594  //TODO: double check if the const_cast is safe here.
595  RCP<const Epetra_MultiVector > r = MultiVecTraitsEpetra::CloneView(toEpetra(mv), index);
596  return rcp(new EpetraMultiVector(Teuchos::rcp_const_cast<Epetra_MultiVector>(r)));
597  }
598 #endif
599 
600  XPETRA_FACTORY_END;
601  }
602 
603  static RCP<const Xpetra::MultiVector<Scalar,LO,GO,Node> >
604  CloneView (const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
605  const Teuchos::Range1D& index)
606  {
607 
608 #ifdef HAVE_XPETRA_TPETRA
609  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
610  //TODO: double check if the const_cast is safe here.
611  RCP<const Tpetra::MultiVector<Scalar,LO,GO,Node> > r = MultiVecTraitsTpetra::CloneView(toTpetra(mv), index);
612  return rcp(new TpetraMultiVector(Teuchos::rcp_const_cast<Tpetra::MultiVector<Scalar,LO,GO,Node> >(r)));
613  }
614 #endif
615 
616 #ifdef HAVE_XPETRA_EPETRA
617  if (mv.getMap()->lib() == Xpetra::UseEpetra) {
618  //TODO: double check if the const_cast is safe here.
619  RCP<const Epetra_MultiVector > r = MultiVecTraitsEpetra::CloneView(toEpetra(mv), index);
620  return rcp(new EpetraMultiVector(Teuchos::rcp_const_cast<Epetra_MultiVector>(r)));
621  }
622 #endif
623 
624  XPETRA_FACTORY_END;
625  }
626 
627  static ptrdiff_t GetGlobalLength( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
628  {
629 
630 #ifdef HAVE_XPETRA_TPETRA
631  if (mv.getMap()->lib() == Xpetra::UseTpetra)
632  return MultiVecTraitsTpetra::GetGlobalLength(toTpetra(mv));
633 #endif
634 
635 #ifdef HAVE_XPETRA_EPETRA
636  if (mv.getMap()->lib() == Xpetra::UseEpetra)
637  return MultiVecTraitsEpetra::GetGlobalLength(toEpetra(mv));
638 #endif
639 
640  XPETRA_FACTORY_END;
641  }
642 
643  static int GetNumberVecs( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
644  {
645 
646 #ifdef HAVE_XPETRA_TPETRA
647  if (mv.getMap()->lib() == Xpetra::UseTpetra)
648  return MultiVecTraitsTpetra::GetNumberVecs(toTpetra(mv));
649 #endif
650 
651 #ifdef HAVE_XPETRA_EPETRA
652  if (mv.getMap()->lib() == Xpetra::UseEpetra)
653  return MultiVecTraitsEpetra::GetNumberVecs(toEpetra(mv));
654 #endif
655 
656  XPETRA_FACTORY_END;
657  }
658 
659  static bool HasConstantStride( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
660  {
661 
662 #ifdef HAVE_XPETRA_TPETRA
663  if (mv.getMap()->lib() == Xpetra::UseTpetra)
664  return MultiVecTraitsTpetra::HasConstantStride(toTpetra(mv));
665 #endif
666 
667 #ifdef HAVE_XPETRA_EPETRA
668  if (mv.getMap()->lib() == Xpetra::UseEpetra)
669  return MultiVecTraitsEpetra::HasConstantStride(toEpetra(mv));
670 #endif
671 
672  XPETRA_FACTORY_END;
673  }
674 
675  static void MvTimesMatAddMv( Scalar alpha, const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
676  const Teuchos::SerialDenseMatrix<int,Scalar>& B,
677  Scalar beta, Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
678  {
679 
680 #ifdef HAVE_BELOS_XPETRA_TIMERS
681  Teuchos::TimeMonitor lcltimer(*mvTimesMatAddMvTimer_);
682 #endif
683 
684 #ifdef HAVE_XPETRA_TPETRA
685  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
686  MultiVecTraitsTpetra::MvTimesMatAddMv(alpha, toTpetra(A), B, beta, toTpetra(mv));
687  return;
688  }
689 #endif
690 
691 #ifdef HAVE_XPETRA_EPETRA
692  if (mv.getMap()->lib() == Xpetra::UseEpetra) {
693  MultiVecTraitsEpetra::MvTimesMatAddMv(alpha, toEpetra(A), B, beta, toEpetra(mv));
694  return;
695  }
696 #endif
697 
698  XPETRA_FACTORY_END;
699  }
700 
701  static void MvAddMv( Scalar alpha, const Xpetra::MultiVector<Scalar,LO,GO,Node>& A, Scalar beta, const Xpetra::MultiVector<Scalar,LO,GO,Node>& B, Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
702  {
703 
704 #ifdef HAVE_XPETRA_TPETRA
705  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
706  MultiVecTraitsTpetra::MvAddMv(alpha, toTpetra(A), beta, toTpetra(B), toTpetra(mv));
707  return;
708  }
709 #endif
710 
711 #ifdef HAVE_XPETRA_EPETRA
712  if (mv.getMap()->lib() == Xpetra::UseEpetra) {
713  MultiVecTraitsEpetra::MvAddMv(alpha, toEpetra(A), beta, toEpetra(B), toEpetra(mv));
714  return;
715  }
716 #endif
717 
718  XPETRA_FACTORY_END;
719  }
720 
721  static void MvScale ( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, Scalar alpha )
722  {
723 
724 #ifdef HAVE_XPETRA_TPETRA
725  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
726  MultiVecTraitsTpetra::MvScale(toTpetra(mv), alpha);
727  return;
728  }
729 #endif
730 
731 #ifdef HAVE_XPETRA_EPETRA
732  if (mv.getMap()->lib() == Xpetra::UseEpetra) {
733  MultiVecTraitsEpetra::MvScale(toEpetra(mv), alpha);
734  return;
735  }
736 #endif
737 
738  XPETRA_FACTORY_END;
739  }
740 
741  static void MvScale ( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, const std::vector<Scalar>& alphas )
742  {
743 
744 #ifdef HAVE_XPETRA_TPETRA
745  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
746  MultiVecTraitsTpetra::MvScale(toTpetra(mv), alphas);
747  return;
748  }
749 #endif
750 
751 #ifdef HAVE_XPETRA_EPETRA
752  if (mv.getMap()->lib() == Xpetra::UseEpetra) {
753  MultiVecTraitsEpetra::MvScale(toEpetra(mv), alphas);
754  return;
755  }
756 #endif
757 
758  XPETRA_FACTORY_END;
759  }
760 
761  static void MvTransMv( Scalar alpha, const Xpetra::MultiVector<Scalar,LO,GO,Node>& A, const Xpetra::MultiVector<Scalar,LO,GO,Node>& B, Teuchos::SerialDenseMatrix<int,Scalar>& C)
762  {
763 
764 #ifdef HAVE_BELOS_XPETRA_TIMERS
765  Teuchos::TimeMonitor lcltimer(*mvTransMvTimer_);
766 #endif
767 
768 #ifdef HAVE_XPETRA_TPETRA
769  if (A.getMap()->lib() == Xpetra::UseTpetra) {
770  MultiVecTraitsTpetra::MvTransMv(alpha, toTpetra(A), toTpetra(B), C);
771  return;
772  }
773 #endif
774 
775 #ifdef HAVE_XPETRA_EPETRA
776  if (A.getMap()->lib() == Xpetra::UseEpetra) {
777  MultiVecTraitsEpetra::MvTransMv(alpha, toEpetra(A), toEpetra(B), C);
778  return;
779  }
780 #endif
781 
782  XPETRA_FACTORY_END;
783  }
784 
785  static void MvDot( const Xpetra::MultiVector<Scalar,LO,GO,Node>& A, const Xpetra::MultiVector<Scalar,LO,GO,Node>& B, std::vector<Scalar> &dots)
786  {
787 
788 #ifdef HAVE_XPETRA_TPETRA
789  if (A.getMap()->lib() == Xpetra::UseTpetra) {
790  MultiVecTraitsTpetra::MvDot(toTpetra(A), toTpetra(B), dots);
791  return;
792  }
793 #endif
794 
795 #ifdef HAVE_XPETRA_EPETRA
796  if (A.getMap()->lib() == Xpetra::UseEpetra) {
797  MultiVecTraitsEpetra::MvDot(toEpetra(A), toEpetra(B), dots);
798  return;
799  }
800 #endif
801 
802  XPETRA_FACTORY_END;
803  }
804 
805  static void MvNorm(const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, std::vector<Teuchos::ScalarTraits<Scalar>::magnitudeType> &normvec, NormType type=TwoNorm)
806  {
807 
808 #ifdef HAVE_XPETRA_TPETRA
809  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
810  MultiVecTraitsTpetra::MvNorm(toTpetra(mv), normvec, type);
811  return;
812  }
813 #endif
814 
815 #ifdef HAVE_XPETRA_EPETRA
816  if (mv.getMap()->lib() == Xpetra::UseEpetra) {
817  MultiVecTraitsEpetra::MvNorm(toEpetra(mv), normvec, type);
818  return;
819  }
820 #endif
821 
822  XPETRA_FACTORY_END;
823  }
824 
825  static void SetBlock( const Xpetra::MultiVector<Scalar,LO,GO,Node>& A, const std::vector<int>& index, Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
826  {
827 
828 #ifdef HAVE_XPETRA_TPETRA
829  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
830  MultiVecTraitsTpetra::SetBlock(toTpetra(A), index, toTpetra(mv));
831  return;
832  }
833 #endif
834 
835 #ifdef HAVE_XPETRA_EPETRA
836  if (mv.getMap()->lib() == Xpetra::UseEpetra) {
837  MultiVecTraitsEpetra::SetBlock(toEpetra(A), index, toEpetra(mv));
838  return;
839  }
840 #endif
841 
842  XPETRA_FACTORY_END;
843  }
844 
845  static void
846  SetBlock (const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
847  const Teuchos::Range1D& index,
848  Xpetra::MultiVector<Scalar,LO,GO,Node>& mv)
849  {
850 
851 #ifdef HAVE_XPETRA_TPETRA
852  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
853  MultiVecTraitsTpetra::SetBlock(toTpetra(A), index, toTpetra(mv));
854  return;
855  }
856 #endif
857 
858 #ifdef HAVE_XPETRA_EPETRA
859  if (mv.getMap()->lib() == Xpetra::UseEpetra) {
860  MultiVecTraitsEpetra::SetBlock(toEpetra(A), index, toEpetra(mv));
861  return;
862  }
863 #endif
864 
865  XPETRA_FACTORY_END;
866  }
867 
868  static void
869  Assign (const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
870  Xpetra::MultiVector<Scalar,LO,GO,Node>& mv)
871  {
872 
873 #ifdef HAVE_XPETRA_TPETRA
874  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
875  MultiVecTraitsTpetra::Assign(toTpetra(A), toTpetra(mv));
876  return;
877  }
878 #endif
879 
880 #ifdef HAVE_XPETRA_EPETRA
881  if (mv.getMap()->lib() == Xpetra::UseEpetra) {
882  MultiVecTraitsEpetra::Assign(toEpetra(A), toEpetra(mv));
883  return;
884  }
885 #endif
886 
887  XPETRA_FACTORY_END;
888  }
889 
890  static void MvRandom( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
891  {
892 
893 #ifdef HAVE_XPETRA_TPETRA
894  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
895  MultiVecTraitsTpetra::MvRandom(toTpetra(mv));
896  return;
897  }
898 #endif
899 
900 #ifdef HAVE_XPETRA_EPETRA
901  if (mv.getMap()->lib() == Xpetra::UseEpetra) {
902  MultiVecTraitsEpetra::MvRandom(toEpetra(mv));
903  return;
904  }
905 #endif
906 
907  XPETRA_FACTORY_END;
908  }
909 
910  static void MvInit( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, Scalar alpha = Teuchos::ScalarTraits<Scalar>::zero() )
911  {
912 
913 #ifdef HAVE_XPETRA_TPETRA
914  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
915  MultiVecTraitsTpetra::MvInit(toTpetra(mv), alpha);
916  return;
917  }
918 #endif
919 
920 #ifdef HAVE_XPETRA_EPETRA
921  if (mv.getMap()->lib() == Xpetra::UseEpetra) {
922  MultiVecTraitsEpetra::MvInit(toEpetra(mv), alpha);
923  return;
924  }
925 #endif
926 
927  XPETRA_FACTORY_END;
928  }
929 
930  static void MvPrint( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, std::ostream& os )
931  {
932 
933 #ifdef HAVE_XPETRA_TPETRA
934  if (mv.getMap()->lib() == Xpetra::UseTpetra) {
935  MultiVecTraitsTpetra::MvPrint(toTpetra(mv), os);
936  return;
937  }
938 #endif
939 
940 #ifdef HAVE_XPETRA_EPETRA
941  if (mv.getMap()->lib() == Xpetra::UseEpetra) {
942  MultiVecTraitsEpetra::MvPrint(toEpetra(mv), os);
943  return;
944  }
945 #endif
946 
947  XPETRA_FACTORY_END;
948  }
949 
950  };
951 
952 } // end of Belos namespace
953 
954 #endif
static bool HasConstantStride(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneCopy(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneViewNonConst(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< int > &index)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneCopy(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvInit(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, Scalar alpha=Teuchos::ScalarTraits< Scalar >::zero())
static void MvTransMv(Scalar alpha, const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const Xpetra::MultiVector< Scalar, LO, GO, Node > &B, Teuchos::SerialDenseMatrix< int, Scalar > &C)
static RCP< const Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneView(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const Teuchos::Range1D &index)
static void MvTimesMatAddMv(Scalar alpha, const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const Teuchos::SerialDenseMatrix< int, Scalar > &B, Scalar beta, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > Clone(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const int numvecs)
static void MvPrint(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, std::ostream &os)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneCopy(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< int > &index)
static void MvInit(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, Scalar alpha=Teuchos::ScalarTraits< Scalar >::zero())
static void SetBlock(const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const std::vector< int > &index, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvTransMv(Scalar alpha, const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const Xpetra::MultiVector< Scalar, LO, GO, Node > &B, Teuchos::SerialDenseMatrix< int, Scalar > &C)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneViewNonConst(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< int > &index)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneViewNonConst(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const Teuchos::Range1D &index)
static void MvScale(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, Scalar alpha)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneCopy(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const Teuchos::Range1D &index)
static void MvDot(const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const Xpetra::MultiVector< Scalar, LO, GO, Node > &B, std::vector< Scalar > &dots)
static RCP< const Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneView(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const Teuchos::Range1D &index)
static void MvAddMv(Scalar alpha, const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, Scalar beta, const Xpetra::MultiVector< Scalar, LO, GO, Node > &B, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static RCP< const Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneView(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< int > &index)
static int GetNumberVecs(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvScale(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< Scalar > &alphas)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneCopy(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< int > &index)
static void SetBlock(const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const Teuchos::Range1D &index, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void SetBlock(const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const std::vector< int > &index, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvScale(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< Scalar > &alphas)
static ptrdiff_t GetGlobalLength(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvDot(const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const Xpetra::MultiVector< Scalar, LO, GO, Node > &B, std::vector< Scalar > &dots)
static RCP< const Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneView(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< int > &index)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneViewNonConst(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const Teuchos::Range1D &index)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > Clone(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const int numvecs)
static void Assign(const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void SetBlock(const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const Teuchos::Range1D &index, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvAddMv(Scalar alpha, const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, Scalar beta, const Xpetra::MultiVector< Scalar, LO, GO, Node > &B, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvTimesMatAddMv(Scalar alpha, const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const Teuchos::SerialDenseMatrix< int, Scalar > &B, Scalar beta, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvNorm(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, std::vector< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &normvec, NormType type=TwoNorm)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneCopy(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const Teuchos::Range1D &index)
static void Assign(const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvNorm(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, std::vector< Teuchos::ScalarTraits< Scalar >::magnitudeType > &normvec, NormType type=TwoNorm)