Xpetra_TpetraMultiVector.hpp
Go to the documentation of this file.
1 // @HEADER
2 //
3 // ***********************************************************************
4 //
5 // Xpetra: A linear algebra interface package
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 XPETRA_TPETRAMULTIVECTOR_HPP
47 #define XPETRA_TPETRAMULTIVECTOR_HPP
48 
49 /* this file is automatically generated - do not edit (see script/tpetra.py) */
50 
52 
53 #include "Xpetra_MultiVector.hpp"
54 
55 #include "Xpetra_TpetraMap.hpp" //TMP
56 #include "Xpetra_Utils.hpp"
57 #include "Xpetra_TpetraImport.hpp"
58 #include "Xpetra_TpetraExport.hpp"
59 
60 #include "Tpetra_MultiVector.hpp"
61 #include "Tpetra_Vector.hpp"
62 
63 namespace Xpetra {
64 
65  // TODO: move that elsewhere
66  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
67  const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node> & toTpetra(const MultiVector< Scalar,LocalOrdinal, GlobalOrdinal, Node> &);
68 
69  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
70  Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node> & toTpetra(MultiVector< Scalar,LocalOrdinal, GlobalOrdinal, Node> &);
71 
72 #ifndef DOXYGEN_SHOULD_SKIP_THIS
73  // forward declaration of TpetraVector, needed to prevent circular inclusions
74  template<class S, class LO, class GO, class N> class TpetraVector;
75 #endif
76 
77 
78  // Because we aren't including the header...
79  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
80  RCP<Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node > > toXpetra(RCP<Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec);
81 
82  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
83  RCP<const Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node > > toXpetra(RCP<const Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec);
84 
85 
86  template <class Scalar = MultiVector<>::scalar_type,
87  class LocalOrdinal = typename MultiVector<Scalar>::local_ordinal_type,
88  class GlobalOrdinal = typename MultiVector<Scalar, LocalOrdinal>::global_ordinal_type,
91  : public virtual MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node >
92  {
93 
94  // The following typedef are used by the XPETRA_DYNAMIC_CAST() macro.
96 
97  public:
98 
100 
101 
103  TpetraMultiVector(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, size_t NumVectors, bool zeroOut=true)
104  : vec_(Teuchos::rcp(new Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node >(toTpetra(map), NumVectors, zeroOut))) { }
105 
108  : vec_(Teuchos::rcp(new Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node >(Tpetra::createCopy(toTpetra(source))))) { }
109 
112  : vec_(Teuchos::rcp(new Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node >(toTpetra(map), A, LDA, NumVectors))) { }
113 
116  : vec_(Teuchos::rcp(new Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node >(toTpetra(map), ArrayOfPtrs, NumVectors))) { }
117 
118 
120  virtual ~TpetraMultiVector() { }
121 
123 
125 
126 
128  void replaceGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value) { XPETRA_MONITOR("TpetraMultiVector::replaceGlobalValue"); vec_->replaceGlobalValue(globalRow, vectorIndex, value); }
129 
131  void sumIntoGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value) { XPETRA_MONITOR("TpetraMultiVector::sumIntoGlobalValue"); vec_->sumIntoGlobalValue(globalRow, vectorIndex, value); }
132 
134  void replaceLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value) { XPETRA_MONITOR("TpetraMultiVector::replaceLocalValue"); vec_->replaceLocalValue(myRow, vectorIndex, value); }
135 
137  void sumIntoLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value) { XPETRA_MONITOR("TpetraMultiVector::sumIntoLocalValue"); vec_->sumIntoLocalValue(myRow, vectorIndex, value); }
138 
140  void putScalar(const Scalar &value) { XPETRA_MONITOR("TpetraMultiVector::putScalar"); vec_->putScalar(value); }
141 
143  void reduce() { XPETRA_MONITOR("TpetraMultiVector::reduce"); vec_->reduce(); }
144 
146 
148 
149 
151  Teuchos::RCP< const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getVector(size_t j) const { XPETRA_MONITOR("TpetraMultiVector::getVector"); return toXpetra(vec_->getVector(j)); }
152 
154  Teuchos::RCP< Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getVectorNonConst(size_t j) { XPETRA_MONITOR("TpetraMultiVector::getVectorNonConst"); return toXpetra(vec_->getVectorNonConst(j)); }
155 
157  Teuchos::ArrayRCP< const Scalar > getData(size_t j) const { XPETRA_MONITOR("TpetraMultiVector::getData"); return vec_->getData(j); }
158 
160  Teuchos::ArrayRCP< Scalar > getDataNonConst(size_t j) { XPETRA_MONITOR("TpetraMultiVector::getDataNonConst"); return vec_->getDataNonConst(j); }
161 
163  void get1dCopy(Teuchos::ArrayView< Scalar > A, size_t LDA) const { XPETRA_MONITOR("TpetraMultiVector::get1dCopy"); vec_->get1dCopy(A, LDA); }
164 
166  void get2dCopy(Teuchos::ArrayView< const Teuchos::ArrayView< Scalar > > ArrayOfPtrs) const { XPETRA_MONITOR("TpetraMultiVector::get2dCopy"); vec_->get2dCopy(ArrayOfPtrs); }
167 
169  Teuchos::ArrayRCP< const Scalar > get1dView() const { XPETRA_MONITOR("TpetraMultiVector::get1dView"); return vec_->get1dView(); }
170 
172  Teuchos::ArrayRCP< Teuchos::ArrayRCP< const Scalar > > get2dView() const { XPETRA_MONITOR("TpetraMultiVector::get2dView"); return vec_->get2dView(); }
173 
175  Teuchos::ArrayRCP< Scalar > get1dViewNonConst() { XPETRA_MONITOR("TpetraMultiVector::get1dViewNonConst"); return vec_->get1dViewNonConst(); }
176 
178  Teuchos::ArrayRCP< Teuchos::ArrayRCP< Scalar > > get2dViewNonConst() { XPETRA_MONITOR("TpetraMultiVector::get2dViewNonConst"); return vec_->get2dViewNonConst(); }
179 
181 
183 
184 
186  void dot(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Teuchos::ArrayView< Scalar > &dots) const { XPETRA_MONITOR("TpetraMultiVector::dot"); vec_->dot(toTpetra(A), dots); }
187 
189  void abs(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A) { XPETRA_MONITOR("TpetraMultiVector::abs"); vec_->abs(toTpetra(A)); }
190 
192  void reciprocal(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A) { XPETRA_MONITOR("TpetraMultiVector::reciprocal"); vec_->reciprocal(toTpetra(A)); }
193 
195  void scale(const Scalar &alpha) { XPETRA_MONITOR("TpetraMultiVector::scale"); vec_->scale(alpha); }
196 
198  void scale(Teuchos::ArrayView< const Scalar > alpha) { XPETRA_MONITOR("TpetraMultiVector::scale"); vec_->scale(alpha); }
199 
201  void scale(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A) { XPETRA_MONITOR("TpetraMultiVector::scale"); vec_->scale(alpha, toTpetra(A)); }
202 
204  void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta) { XPETRA_MONITOR("TpetraMultiVector::update"); vec_->update(alpha, toTpetra(A), beta); }
205 
207  void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, const Scalar &gamma) { XPETRA_MONITOR("TpetraMultiVector::update"); vec_->update(alpha, toTpetra(A), beta, toTpetra(B), gamma); }
208 
210  void norm1(const Teuchos::ArrayView< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const { XPETRA_MONITOR("TpetraMultiVector::norm1"); vec_->norm1(norms); }
211 
213  void norm2(const Teuchos::ArrayView< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const { XPETRA_MONITOR("TpetraMultiVector::norm2"); vec_->norm2(norms); }
214 
216  void normInf(const Teuchos::ArrayView< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const { XPETRA_MONITOR("TpetraMultiVector::normInf"); vec_->normInf(norms); }
217 
219  void meanValue(const Teuchos::ArrayView< Scalar > &means) const { XPETRA_MONITOR("TpetraMultiVector::meanValue"); vec_->meanValue(means); }
220 
222  void multiply(Teuchos::ETransp transA, Teuchos::ETransp transB, const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, const Scalar &beta) { XPETRA_MONITOR("TpetraMultiVector::multiply"); vec_->multiply(transA, transB, alpha, toTpetra(A), toTpetra(B), beta); }
223 
225 
227 
228 
230  size_t getNumVectors() const { XPETRA_MONITOR("TpetraMultiVector::getNumVectors"); return vec_->getNumVectors(); }
231 
233  size_t getLocalLength() const { XPETRA_MONITOR("TpetraMultiVector::getLocalLength"); return vec_->getLocalLength(); }
234 
236  global_size_t getGlobalLength() const { XPETRA_MONITOR("TpetraMultiVector::getGlobalLength"); return vec_->getGlobalLength(); }
237 
239 
241 
242 
244  std::string description() const { XPETRA_MONITOR("TpetraMultiVector::description"); return vec_->description(); }
245 
247  void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const { XPETRA_MONITOR("TpetraMultiVector::describe"); vec_->describe(out, verbLevel); }
248 
250 
252  void elementWiseMultiply(Scalar scalarAB, const Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> &A, const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> &B, Scalar scalarThis); // definition at the end of this file
253  //TODO: void elementWiseMultiply(Scalar scalarAB, const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, Scalar scalarThis){ vec_->elementWiseMultiply(scalarAB, toTpetra(A), toTpetra(B), scalarThis); }
254 
256  void randomize(bool bUseXpetraImplementation = false) {
257  XPETRA_MONITOR("TpetraMultiVector::randomize");
258 
259  if(bUseXpetraImplementation)
261  else
262  vec_->randomize();
263  }
264 
265  //{@
266  // Implements DistObject interface
267 
268  Teuchos::RCP< const Map<LocalOrdinal,GlobalOrdinal,Node> > getMap() const { XPETRA_MONITOR("TpetraMultiVector::getMap"); return toXpetra(vec_->getMap()); }
269 
271  XPETRA_MONITOR("TpetraMultiVector::doImport");
272 
273  XPETRA_DYNAMIC_CAST(const TpetraMultiVectorClass, source, tSource, "Xpetra::TpetraMultiVector::doImport only accept Xpetra::TpetraMultiVector as input arguments."); //TODO: remove and use toTpetra()
275  this->getTpetra_MultiVector()->doImport(*v, toTpetra(importer), toTpetra(CM));
276  }
277 
279  XPETRA_MONITOR("TpetraMultiVector::doExport");
280 
281  XPETRA_DYNAMIC_CAST(const TpetraMultiVectorClass, dest, tDest, "Xpetra::TpetraMultiVector::doImport only accept Xpetra::TpetraMultiVector as input arguments."); //TODO: remove and use toTpetra()
283  this->getTpetra_MultiVector()->doExport(*v, toTpetra(importer), toTpetra(CM));
284 
285  }
286 
288  XPETRA_MONITOR("TpetraMultiVector::doImport");
289 
290  XPETRA_DYNAMIC_CAST(const TpetraMultiVectorClass, source, tSource, "Xpetra::TpetraMultiVector::doImport only accept Xpetra::TpetraMultiVector as input arguments."); //TODO: remove and use toTpetra()
292  this->getTpetra_MultiVector()->doImport(*v, toTpetra(exporter), toTpetra(CM));
293 
294  }
295 
297  XPETRA_MONITOR("TpetraMultiVector::doExport");
298 
299  XPETRA_DYNAMIC_CAST(const TpetraMultiVectorClass, dest, tDest, "Xpetra::TpetraMultiVector::doImport only accept Xpetra::TpetraMultiVector as input arguments."); //TODO: remove and use toTpetra()
301  this->getTpetra_MultiVector()->doExport(*v, toTpetra(exporter), toTpetra(CM));
302 
303  }
304 
306  this->getTpetra_MultiVector()->replaceMap(toTpetra(map));
307  }
308 
309  template<class Node2>
312  //toXpetra(vec_->clone(node2));
313  }
314 
316 
318 
319 
321  TpetraMultiVector(const Teuchos::RCP<Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node> > &vec) : vec_(vec) { } //TODO removed const
322 
325 
327  void setSeed(unsigned int seed) { XPETRA_MONITOR("TpetraMultiVector::seedrandom"); Teuchos::ScalarTraits< Scalar >::seedrandom(seed); }
328 
329 
330 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
332  //typedef typename Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>::host_execution_space host_execution_space;
333  //typedef typename Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>::dev_execution_space dev_execution_space;
334 
345  template<class TargetDeviceType>
346  typename Kokkos::Impl::if_c<
347  Kokkos::Impl::is_same<
348  typename dual_view_type::t_dev_um::execution_space::memory_space,
349  typename TargetDeviceType::memory_space>::value,
350  typename dual_view_type::t_dev_um,
351  typename dual_view_type::t_host_um>::type
352  getLocalView () const {
353  return this->MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node >::template getLocalView<TargetDeviceType>();
354  }
355 
356  typename dual_view_type::t_host_um getHostLocalView () const {
357  return subview(vec_->getDualView().template view<typename dual_view_type::host_mirror_space> (),
358  Kokkos::ALL(), Kokkos::ALL());
359  }
360 
361  typename dual_view_type::t_dev_um getDeviceLocalView() const {
362  return subview(vec_->getDualView().template view<typename dual_view_type::t_dev_um::execution_space> (),
363  Kokkos::ALL(), Kokkos::ALL());
364  }
365 
366 #endif
367 
369 
370  protected:
373  virtual void
375  {
377  const this_type* rhsPtr = dynamic_cast<const this_type*> (&rhs);
379  rhsPtr == NULL, std::invalid_argument, "Xpetra::MultiVector::operator=:"
380  " The left-hand side (LHS) of the assignment has a different type than "
381  "the right-hand side (RHS). The LHS has type Xpetra::TpetraMultiVector"
382  " (which means it wraps a Tpetra::MultiVector), but the RHS has some "
383  "other type. This probably means that the RHS wraps an "
384  "Epetra_MultiVector. Xpetra::MultiVector does not currently implement "
385  "assignment from an Epetra object to a Tpetra object, though this could"
386  " be added with sufficient interest.");
387 
388  typedef Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> TMV;
389  RCP<const TMV> rhsImpl = rhsPtr->getTpetra_MultiVector ();
390  RCP<TMV> lhsImpl = this->getTpetra_MultiVector ();
391 
393  rhsImpl.is_null (), std::logic_error, "Xpetra::MultiVector::operator= "
394  "(in Xpetra::TpetraMultiVector::assign): *this (the right-hand side of "
395  "the assignment) has a null RCP<Tpetra::MultiVector> inside. Please "
396  "report this bug to the Xpetra developers.");
398  lhsImpl.is_null (), std::logic_error, "Xpetra::MultiVector::operator= "
399  "(in Xpetra::TpetraMultiVector::assign): The left-hand side of the "
400  "assignment has a null RCP<Tpetra::MultiVector> inside. Please report "
401  "this bug to the Xpetra developers.");
402 
403  Tpetra::deep_copy (*lhsImpl, *rhsImpl);
404  }
405 
406  private:
409 
410  }; // TpetraMultiVector class
411 
412  // TODO: move that elsewhere
413  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
414  const Tpetra::MultiVector< Scalar,LocalOrdinal, GlobalOrdinal, Node> & toTpetra(const MultiVector< Scalar,LocalOrdinal, GlobalOrdinal, Node> &x) {
416  XPETRA_DYNAMIC_CAST(const TpetraMultiVectorClass, x, tX, "toTpetra");
417  return *tX.getTpetra_MultiVector();
418  }
419 
420  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
421  Tpetra::MultiVector< Scalar,LocalOrdinal, GlobalOrdinal, Node> & toTpetra(MultiVector< Scalar,LocalOrdinal, GlobalOrdinal, Node> &x) {
423  XPETRA_DYNAMIC_CAST( TpetraMultiVectorClass, x, tX, "toTpetra");
424  return *tX.getTpetra_MultiVector();
425  }
426  //
427 
428 
429  // Things we actually need
430  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
431  RCP<MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node > > toXpetra(RCP<Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec) {
432  if (!vec.is_null())
434 
435  return Teuchos::null;
436  }
437 
438  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
439  RCP<const MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node > > toXpetra(RCP<const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec) {
440  if (!vec.is_null())
442 
443  return Teuchos::null;
444  }
445 
446 #ifdef HAVE_XPETRA_EPETRA
447 
448 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
449  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
450 
451  // specialization for TpetraMultiVector on EpetraNode and GO=int
452  template <class Scalar>
453  class TpetraMultiVector<Scalar,int,int,EpetraNode>
454  : public virtual MultiVector< Scalar, int, int, EpetraNode >
455  {
456  typedef int LocalOrdinal;
457  typedef int GlobalOrdinal;
458  typedef EpetraNode Node;
459 
460  // The following typedef are used by the XPETRA_DYNAMIC_CAST() macro.
462 
463  public:
464 
466 
467 
469  TpetraMultiVector(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, size_t NumVectors, bool zeroOut=true) {
470  XPETRA_TPETRA_ETI_EXCEPTION("TpetraMultiVector<int,int>", "TpetraMultiVector<int,int>", "int");
471  }
472 
475  XPETRA_TPETRA_ETI_EXCEPTION("TpetraMultiVector<int,int>", "TpetraMultiVector<int,int>", "int");
476  }
477 
480  XPETRA_TPETRA_ETI_EXCEPTION("TpetraMultiVector<int,int>", "TpetraMultiVector<int,int>", "int");
481  }
482 
485  XPETRA_TPETRA_ETI_EXCEPTION("TpetraMultiVector<int,int>", "TpetraMultiVector<int,int>", "int");
486  }
487 
488 
490  virtual ~TpetraMultiVector() { }
491 
493 
495 
496 
498  void replaceGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value) { }
499 
501  void sumIntoGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value) { }
502 
504  void replaceLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value) { }
505 
507  void sumIntoLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value) { }
508 
510  void putScalar(const Scalar &value) { }
511 
513  void reduce() { }
514 
516 
518 
519 
522 
525 
528 
531 
533  void get1dCopy(Teuchos::ArrayView< Scalar > A, size_t LDA) const { }
534 
536  void get2dCopy(Teuchos::ArrayView< const Teuchos::ArrayView< Scalar > > ArrayOfPtrs) const { }
537 
540 
543 
546 
549 
551 
553 
554 
557 
560 
563 
565  void scale(const Scalar &alpha) { }
566 
569 
571  void scale(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A) { }
572 
574  void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta) { }
575 
577  void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, const Scalar &gamma) { }
578 
581 
584 
587 
589  void meanValue(const Teuchos::ArrayView< Scalar > &means) const { }
590 
593 
595 
597 
598 
600  size_t getNumVectors() const { return 0; }
601 
603  size_t getLocalLength() const { return 0; }
604 
606  global_size_t getGlobalLength() const { return 0; }
607 
609 
611 
612 
614  std::string description() const { return std::string(""); }
615 
618 
620 
623 
625  void randomize(bool bUseXpetraImplementation = false) { }
626 
627  //{@
628  // Implements DistObject interface
629 
631 
633 
635 
637 
639 
641 
642  template<class Node2>
643  RCP<MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node2> > clone(const RCP<Node2> &node2) const { return Teuchos::null; }
644 
646 
648 
649 
651  TpetraMultiVector(const Teuchos::RCP<Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node> > &vec) { }
652 
655 
657  void setSeed(unsigned int seed) { }
658 
659 
660 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
662  //typedef typename Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>::host_execution_space host_execution_space;
663  //typedef typename Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>::dev_execution_space dev_execution_space;
664 
675  template<class TargetDeviceType>
676  typename Kokkos::Impl::if_c<
677  Kokkos::Impl::is_same<
678  typename dual_view_type::t_dev_um::execution_space::memory_space,
679  typename TargetDeviceType::memory_space>::value,
680  typename dual_view_type::t_dev_um,
681  typename dual_view_type::t_host_um>::type
682  getLocalView () const {
683  typename Kokkos::Impl::if_c<
684  Kokkos::Impl::is_same<
685  typename dual_view_type::t_dev_um::execution_space::memory_space,
686  typename TargetDeviceType::memory_space>::value,
687  typename dual_view_type::t_dev_um,
688  typename dual_view_type::t_host_um>::type dummy;
689  return dummy;
690  }
691 
692  typename dual_view_type::t_host_um getHostLocalView () const {
693  //return subview(vec_->getDualView().template view<typename dual_view_type::host_mirror_space> (),
694  // Kokkos::ALL(), Kokkos::ALL());
695  return typename dual_view_type::t_host_um();
696  }
697 
698  typename dual_view_type::t_dev_um getDeviceLocalView() const {
699  //return subview(vec_->getDualView().template view<typename dual_view_type::t_dev_um::execution_space> (),
700  // Kokkos::ALL(), Kokkos::ALL());
701  return typename dual_view_type::t_dev_um();
702  }
703 
704 #endif
705 
707 
708  protected:
711  virtual void
713  { }
714  }; // TpetraMultiVector class (specialization GO=int, NO=EpetraNode)
715 #endif
716 
717 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_LONG_LONG))) || \
718  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_LONG_LONG))))
719 
720  // specialization for TpetraMultiVector on EpetraNode and GO=long long
721  template <class Scalar>
722  class TpetraMultiVector<Scalar,int,long long,EpetraNode>
723  : public virtual MultiVector< Scalar, int, long long, EpetraNode >
724  {
725  typedef int LocalOrdinal;
726  typedef long long GlobalOrdinal;
727  typedef EpetraNode Node;
728 
729  // The following typedef are used by the XPETRA_DYNAMIC_CAST() macro.
731 
732  public:
733 
735 
736 
738  TpetraMultiVector(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, size_t NumVectors, bool zeroOut=true) {
739  XPETRA_TPETRA_ETI_EXCEPTION("TpetraMultiVector<int,int>", "TpetraMultiVector<int,int>", "int");
740  }
741 
744  XPETRA_TPETRA_ETI_EXCEPTION("TpetraMultiVector<int,int>", "TpetraMultiVector<int,int>", "int");
745  }
746 
749  XPETRA_TPETRA_ETI_EXCEPTION("TpetraMultiVector<int,int>", "TpetraMultiVector<int,int>", "int");
750  }
751 
754  XPETRA_TPETRA_ETI_EXCEPTION("TpetraMultiVector<int,int>", "TpetraMultiVector<int,int>", "int");
755  }
756 
757 
759  virtual ~TpetraMultiVector() { }
760 
762 
764 
765 
767  void replaceGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value) { }
768 
770  void sumIntoGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value) { }
771 
773  void replaceLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value) { }
774 
776  void sumIntoLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value) { }
777 
779  void putScalar(const Scalar &value) { }
780 
782  void reduce() { }
783 
785 
787 
788 
791 
794 
797 
800 
802  void get1dCopy(Teuchos::ArrayView< Scalar > A, size_t LDA) const { }
803 
805  void get2dCopy(Teuchos::ArrayView< const Teuchos::ArrayView< Scalar > > ArrayOfPtrs) const { }
806 
809 
812 
815 
818 
820 
822 
823 
826 
829 
832 
834  void scale(const Scalar &alpha) { }
835 
838 
840  void scale(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A) { }
841 
843  void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta) { }
844 
846  void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, const Scalar &gamma) { }
847 
850 
853 
856 
858  void meanValue(const Teuchos::ArrayView< Scalar > &means) const { }
859 
862 
864 
866 
867 
869  size_t getNumVectors() const { return 0; }
870 
872  size_t getLocalLength() const { return 0; }
873 
875  global_size_t getGlobalLength() const { return 0; }
876 
878 
880 
881 
883  std::string description() const { return std::string(""); }
884 
887 
889 
892 
894  void randomize(bool bUseXpetraImplementation = false) { }
895 
896  //{@
897  // Implements DistObject interface
898 
900 
902 
904 
906 
908 
910 
911  template<class Node2>
912  RCP<MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node2> > clone(const RCP<Node2> &node2) const { return Teuchos::null; }
913 
915 
917 
918 
920  TpetraMultiVector(const Teuchos::RCP<Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node> > &vec) { }
921 
924 
926  void setSeed(unsigned int seed) { }
927 
928 
929 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
931  //typedef typename Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>::host_execution_space host_execution_space;
932  //typedef typename Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>::dev_execution_space dev_execution_space;
933 
944  template<class TargetDeviceType>
945  typename Kokkos::Impl::if_c<
946  Kokkos::Impl::is_same<
947  typename dual_view_type::t_dev_um::execution_space::memory_space,
948  typename TargetDeviceType::memory_space>::value,
949  typename dual_view_type::t_dev_um,
950  typename dual_view_type::t_host_um>::type
951  getLocalView () const {
952  typename Kokkos::Impl::if_c<
953  Kokkos::Impl::is_same<
954  typename dual_view_type::t_dev_um::execution_space::memory_space,
955  typename TargetDeviceType::memory_space>::value,
956  typename dual_view_type::t_dev_um,
957  typename dual_view_type::t_host_um>::type dummy;
958  return dummy;
959  }
960 
961  typename dual_view_type::t_host_um getHostLocalView () const {
962  //return subview(vec_->getDualView().template view<typename dual_view_type::host_mirror_space> (),
963  // Kokkos::ALL(), Kokkos::ALL());
964  return typename dual_view_type::t_host_um();
965  }
966 
967  typename dual_view_type::t_dev_um getDeviceLocalView() const {
968  //return subview(vec_->getDualView().template view<typename dual_view_type::t_dev_um::execution_space> (),
969  // Kokkos::ALL(), Kokkos::ALL());
970  return typename dual_view_type::t_dev_um();
971  }
972 
973 #endif
974 
976 
977  protected:
980  virtual void
982  { }
983  }; // TpetraMultiVector class (specialization GO=int, NO=EpetraNode)
984 
985 #endif // TpetraMultiVector class (specialization GO=long long, NO=EpetraNode)
986 
987 #endif // HAVE_XPETRA_EPETRA
988 
989 } // Xpetra namespace
990 
991 // Following header file inculsion is needed for the dynamic_cast to TpetraVector in elementWiseMultiply (because we cannot dynamic_cast if target is not a complete type)
992 // It is included here to avoid circular dependency between Vector and MultiVector
993 // TODO: there is certainly a more elegant solution...
994 #include "Xpetra_TpetraVector.hpp"
995 
996 namespace Xpetra {
997  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
999  XPETRA_MONITOR("TpetraMultiVector::elementWiseMultiply");
1000 
1001  // XPETRA_DYNAMIC_CAST won't take TpetraVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>
1002  // as an argument, hence the following typedef.
1004  XPETRA_DYNAMIC_CAST(const tpv, A, tA, "Xpetra::TpetraMultiVectorMatrix->multiply() only accept Xpetra::TpetraMultiVector as input arguments.");
1005  XPETRA_DYNAMIC_CAST(const TpetraMultiVector, B, tB, "Xpetra::TpetraMultiVectorMatrix->multiply() only accept Xpetra::TpetraMultiVector as input arguments.");
1006  vec_->elementWiseMultiply(scalarAB, *tA.getTpetra_Vector(), *tB.getTpetra_MultiVector(), scalarThis);
1007  }
1008 
1009 } // Xpetra namespace
1010 
1011 #define XPETRA_TPETRAMULTIVECTOR_SHORT
1012 #endif // XPETRA_TPETRAMULTIVECTOR_HPP
void replaceMap(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map)
size_t getLocalLength() const
Local number of rows on the calling process.
void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta)
Update multi-vector values with scaled values of A, this = beta*this + alpha*A.
void norm1(const Teuchos::ArrayView< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
Compute 1-norm of each vector in multi-vector.
Teuchos::ArrayRCP< const Scalar > get1dView() const
Const persisting (1-D) view of this multivector&#39;s local values.
RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_MultiVector() const
Get the underlying Tpetra multivector.
Teuchos::ArrayRCP< Teuchos::ArrayRCP< const Scalar > > get2dView() const
Return const persisting pointers to values.
void doExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
void dot(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Teuchos::ArrayView< Scalar > &dots) const
Compute dot product of each corresponding pair of vectors, dots[i] = this[i].dot(A[i]).
void scale(const Scalar &alpha)
Scale the current values of a multi-vector, this = alpha*this.
TpetraMultiVector(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, size_t NumVectors, bool zeroOut=true)
Basic constuctor.
void replaceGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value)
Replace value, using global (row) index.
void setSeed(unsigned int seed)
Set seed for Random function.
TpetraMultiVector(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, const Teuchos::ArrayView< const Teuchos::ArrayView< const Scalar > > &ArrayOfPtrs, size_t NumVectors)
Create multivector by copying array of views of local data.
RCP< MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node2 > > clone(const RCP< Node2 > &node2) const
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with the given verbosity level to a FancyOStream.
Teuchos::ArrayRCP< Teuchos::ArrayRCP< Scalar > > get2dViewNonConst()
Return non-const persisting pointers to values.
void reduce()
Sum values of a locally replicated multivector across all processes.
void dot(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Teuchos::ArrayView< Scalar > &dots) const
Compute dot product of each corresponding pair of vectors, dots[i] = this[i].dot(A[i]).
void normInf(const Teuchos::ArrayView< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
Compute Inf-norm of each vector in multi-vector.
void scale(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
Replace multi-vector values with scaled values of A, this = alpha*A.
std::string description() const
A simple one-line description of this object.
virtual ~TpetraMultiVector()
Destructor (virtual for memory safety of derived classes).
TpetraMultiVector(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, const Teuchos::ArrayView< const Teuchos::ArrayView< const Scalar > > &ArrayOfPtrs, size_t NumVectors)
Create multivector by copying array of views of local data.
void reciprocal(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
Put element-wise reciprocal values of input Multi-vector in target, this(i,j) = 1/A(i,j).
Teuchos::ArrayRCP< Scalar > getDataNonConst(size_t j)
View of the local values in a particular vector of this multivector.
void doExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Teuchos::RCP< Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getVectorNonConst(size_t j)
Return a Vector which is a nonconst view of column j.
virtual ~TpetraMultiVector()
Destructor (virtual for memory safety of derived classes).
void setSeed(unsigned int seed)
Set seed for Random function.
void elementWiseMultiply(Scalar scalarAB, const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, Scalar scalarThis)
Element-wise multiply of a Vector A with a TpetraMultiVector B.
TpetraMultiVector(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, size_t NumVectors, bool zeroOut=true)
Basic constuctor.
void doExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
TpetraMultiVector(const Teuchos::RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &vec)
TpetraMultiVector constructor to wrap a Tpetra::MultiVector object.
size_t getNumVectors() const
Number of columns in the multivector.
RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_MultiVector() const
Get the underlying Tpetra multivector.
void norm2(const Teuchos::ArrayView< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
The Map describing the parallel distribution of this object.
void elementWiseMultiply(Scalar scalarAB, const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, Scalar scalarThis)
Element-wise multiply of a Vector A with a TpetraMultiVector B.
Teuchos::RCP< Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getVectorNonConst(size_t j)
Return a Vector which is a nonconst view of column j.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
virtual void Xpetra_randomize()
Set multi-vector values to random numbers. XPetra implementation.
RCP< MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node2 > > clone(const RCP< Node2 > &node2) const
void doImport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Import data into this object using an Export object ("reverse mode").
std::string description() const
A simple one-line description of this object.
Xpetra namespace
void normInf(const Teuchos::ArrayView< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
Compute Inf-norm of each vector in multi-vector.
size_t getNumVectors() const
Number of columns in the multivector.
void setSeed(unsigned int seed)
Set seed for Random function.
Teuchos::ArrayRCP< Scalar > getDataNonConst(size_t j)
View of the local values in a particular vector of this multivector.
void scale(const Scalar &alpha)
Scale the current values of a multi-vector, this = alpha*this.
TpetraMultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraMultiVectorClass
void doImport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
void randomize(bool bUseXpetraImplementation=false)
Set multi-vector values to random numbers.
void scale(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
Replace multi-vector values with scaled values of A, this = alpha*A.
void replaceMap(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map)
Teuchos::ArrayRCP< Teuchos::ArrayRCP< const Scalar > > get2dView() const
Return const persisting pointers to values.
void reduce()
Sum values of a locally replicated multivector across all processes.
void replaceGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value)
Replace value, using global (row) index.
void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, const Scalar &gamma)
Update multi-vector with scaled values of A and B, this = gamma*this + alpha*A + beta*B.
Teuchos::ArrayRCP< const Scalar > getData(size_t j) const
Const view of the local values in a particular vector of this multivector.
TpetraMultiVector(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source)
Copy constructor (performs a deep copy).
virtual void assign(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &rhs)
Implementation of the assignment operator (operator=); does a deep copy.
void replaceLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value)
Replace value, using local (row) index.
void get2dCopy(Teuchos::ArrayView< const Teuchos::ArrayView< Scalar > > ArrayOfPtrs) const
Fill the given array with a copy of this multivector&#39;s local values.
Teuchos::RCP< const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getVector(size_t j) const
Return a Vector which is a const view of column j.
void replaceGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value)
Replace value, using global (row) index.
void doImport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import data into this object using an Import object ("forward mode").
void meanValue(const Teuchos::ArrayView< Scalar > &means) const
Compute mean (average) value of each vector in multi-vector. The outcome of this routine is undefined...
TpetraMultiVector(const Teuchos::RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &vec)
TpetraMultiVector constructor to wrap a Tpetra::MultiVector object.
void get1dCopy(Teuchos::ArrayView< Scalar > A, size_t LDA) const
Fill the given array with a copy of this multivector&#39;s local values.
void elementWiseMultiply(Scalar scalarAB, const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, Scalar scalarThis)
Element-wise multiply of a Vector A with a TpetraMultiVector B.
void norm2(const Teuchos::ArrayView< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
void multiply(Teuchos::ETransp transA, Teuchos::ETransp transB, const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, const Scalar &beta)
Matrix-matrix multiplication: this = beta*this + alpha*op(A)*op(B).
Teuchos::ArrayRCP< Scalar > get1dViewNonConst()
Nonconst persisting (1-D) view of this multivector&#39;s local values.
TpetraMultiVector(const Teuchos::RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &vec)
TpetraMultiVector constructor to wrap a Tpetra::MultiVector object.
Teuchos::ArrayRCP< Teuchos::ArrayRCP< Scalar > > get2dViewNonConst()
Return non-const persisting pointers to values.
void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, const Scalar &gamma)
Update multi-vector with scaled values of A and B, this = gamma*this + alpha*A + beta*B.
void abs(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
Put element-wise absolute values of input Multi-vector in target: A = abs(this).
void putScalar(const Scalar &value)
Set all values in the multivector with the given value.
size_t getLocalLength() const
Local number of rows on the calling process.
void multiply(Teuchos::ETransp transA, Teuchos::ETransp transB, const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, const Scalar &beta)
Matrix-matrix multiplication: this = beta*this + alpha*op(A)*op(B).
void putScalar(const Scalar &value)
Set all values in the multivector with the given value.
TpetraMultiVector(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, const Teuchos::ArrayView< const Scalar > &A, size_t LDA, size_t NumVectors)
Create multivector by copying two-dimensional array of local data.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
void sumIntoGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value)
Add value to existing value, using global (row) index.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
The Map describing the parallel distribution of this object.
void get2dCopy(Teuchos::ArrayView< const Teuchos::ArrayView< Scalar > > ArrayOfPtrs) const
Fill the given array with a copy of this multivector&#39;s local values.
void multiply(Teuchos::ETransp transA, Teuchos::ETransp transB, const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, const Scalar &beta)
Matrix-matrix multiplication: this = beta*this + alpha*op(A)*op(B).
void norm1(const Teuchos::ArrayView< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
Compute 1-norm of each vector in multi-vector.
std::string description() const
A simple one-line description of this object.
void scale(Teuchos::ArrayView< const Scalar > alpha)
Scale the current values of a multi-vector, this[j] = alpha[j]*this[j].
void dot(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Teuchos::ArrayView< Scalar > &dots) const
Compute dot product of each corresponding pair of vectors, dots[i] = this[i].dot(A[i]).
size_t getNumVectors() const
Number of columns in the multivector.
Teuchos::ArrayRCP< Scalar > get1dViewNonConst()
Nonconst persisting (1-D) view of this multivector&#39;s local values.
void scale(Teuchos::ArrayView< const Scalar > alpha)
Scale the current values of a multi-vector, this[j] = alpha[j]*this[j].
void replaceLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value)
Replace value, using local (row) index.
virtual void assign(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &rhs)
Implementation of the assignment operator (operator=); does a deep copy.
void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta)
Update multi-vector values with scaled values of A, this = beta*this + alpha*A.
void reduce()
Sum values of a locally replicated multivector across all processes.
Teuchos::ArrayRCP< Scalar > getDataNonConst(size_t j)
View of the local values in a particular vector of this multivector.
void get2dCopy(Teuchos::ArrayView< const Teuchos::ArrayView< Scalar > > ArrayOfPtrs) const
Fill the given array with a copy of this multivector&#39;s local values.
#define XPETRA_TPETRA_ETI_EXCEPTION(cl, obj, go)
void sumIntoLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value)
Add value to existing value, using local (row) index.
void doImport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
void sumIntoGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value)
Add value to existing value, using global (row) index.
void scale(const Scalar &alpha)
Scale the current values of a multi-vector, this = alpha*this.
void reciprocal(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
Put element-wise reciprocal values of input Multi-vector in target, this(i,j) = 1/A(i,j).
virtual ~TpetraMultiVector()
Destructor (virtual for memory safety of derived classes).
LocalOrdinal local_ordinal_type
static void seedrandom(unsigned int s)
#define XPETRA_DYNAMIC_CAST(type, obj, newObj, exceptionMsg)
size_t global_size_t
Global size_t object.
TpetraMultiVector(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, size_t NumVectors, bool zeroOut=true)
Basic constuctor.
void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, const Scalar &gamma)
Update multi-vector with scaled values of A and B, this = gamma*this + alpha*A + beta*B.
void doExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
static const EVerbosityLevel verbLevel_default
void abs(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
Put element-wise absolute values of input Multi-vector in target: A = abs(this).
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
The Map describing the parallel distribution of this object.
Teuchos::RCP< const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getVector(size_t j) const
Return a Vector which is a const view of column j.
void doImport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
void get1dCopy(Teuchos::ArrayView< Scalar > A, size_t LDA) const
Fill the given array with a copy of this multivector&#39;s local values.
void get1dCopy(Teuchos::ArrayView< Scalar > A, size_t LDA) const
Fill the given array with a copy of this multivector&#39;s local values.
void reciprocal(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
Put element-wise reciprocal values of input Multi-vector in target, this(i,j) = 1/A(i,j).
global_size_t getGlobalLength() const
Global number of rows in the multivector.
RCP< const Tpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > toTpetra(const RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph)
GlobalOrdinal global_ordinal_type
void norm2(const Teuchos::ArrayView< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
TpetraMultiVector(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, const Teuchos::ArrayView< const Teuchos::ArrayView< const Scalar > > &ArrayOfPtrs, size_t NumVectors)
Create multivector by copying array of views of local data.
RCP< const CrsGraph< int, GlobalOrdinal, Node > > toXpetra(const Epetra_CrsGraph &g)
Teuchos::RCP< Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getVectorNonConst(size_t j)
Return a Vector which is a nonconst view of column j.
Teuchos::RCP< const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getVector(size_t j) const
Return a Vector which is a const view of column j.
void putScalar(const Scalar &value)
Set all values in the multivector with the given value.
size_t getLocalLength() const
Local number of rows on the calling process.
void scale(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
Replace multi-vector values with scaled values of A, this = alpha*A.
TpetraMultiVector(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, const Teuchos::ArrayView< const Scalar > &A, size_t LDA, size_t NumVectors)
Create multivector by copying two-dimensional array of local data.
TpetraMultiVector(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source)
Copy constructor (performs a deep copy).
void sumIntoGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value)
Add value to existing value, using global (row) index.
TpetraMultiVector(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, const Teuchos::ArrayView< const Scalar > &A, size_t LDA, size_t NumVectors)
Create multivector by copying two-dimensional array of local data.
void meanValue(const Teuchos::ArrayView< Scalar > &means) const
Compute mean (average) value of each vector in multi-vector. The outcome of this routine is undefined...
Teuchos::ArrayRCP< const Scalar > get1dView() const
Const persisting (1-D) view of this multivector&#39;s local values.
Teuchos::ArrayRCP< const Scalar > get1dView() const
Const persisting (1-D) view of this multivector&#39;s local values.
void normInf(const Teuchos::ArrayView< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
Compute Inf-norm of each vector in multi-vector.
virtual void assign(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &rhs)
Implementation of the assignment operator (operator=); does a deep copy.
void scale(Teuchos::ArrayView< const Scalar > alpha)
Scale the current values of a multi-vector, this[j] = alpha[j]*this[j].
RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec_
The Tpetra::MultiVector which this class wraps.
void doExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export data into this object using an Import object ("reverse mode").
CombineMode
Xpetra::Combine Mode enumerable type.
void doImport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
TpetraMultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraMultiVectorClass
#define XPETRA_MONITOR(funcName)
void meanValue(const Teuchos::ArrayView< Scalar > &means) const
Compute mean (average) value of each vector in multi-vector. The outcome of this routine is undefined...
Teuchos::ArrayRCP< const Scalar > getData(size_t j) const
Const view of the local values in a particular vector of this multivector.
Teuchos::ArrayRCP< Scalar > get1dViewNonConst()
Nonconst persisting (1-D) view of this multivector&#39;s local values.
TpetraMultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraMultiVectorClass
global_size_t getGlobalLength() const
Global number of rows in the multivector.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with the given verbosity level to a FancyOStream.
RCP< MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node2 > > clone(const RCP< Node2 > &node2) const
void abs(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
Put element-wise absolute values of input Multi-vector in target: A = abs(this).
void randomize(bool bUseXpetraImplementation=false)
Set multi-vector values to random numbers.
Teuchos::ArrayRCP< Teuchos::ArrayRCP< Scalar > > get2dViewNonConst()
Return non-const persisting pointers to values.
void sumIntoLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value)
Add value to existing value, using local (row) index.
void doExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Export data into this object using an Export object ("forward mode").
TpetraMultiVector(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source)
Copy constructor (performs a deep copy).
void norm1(const Teuchos::ArrayView< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
Compute 1-norm of each vector in multi-vector.
Teuchos::ArrayRCP< Teuchos::ArrayRCP< const Scalar > > get2dView() const
Return const persisting pointers to values.
void randomize(bool bUseXpetraImplementation=false)
Set multi-vector values to random numbers.
Teuchos::ArrayRCP< const Scalar > getData(size_t j) const
Const view of the local values in a particular vector of this multivector.
void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta)
Update multi-vector values with scaled values of A, this = beta*this + alpha*A.
void replaceLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value)
Replace value, using local (row) index.
void sumIntoLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value)
Add value to existing value, using local (row) index.
global_size_t getGlobalLength() const
Global number of rows in the multivector.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with the given verbosity level to a FancyOStream.
void replaceMap(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map)
bool is_null() const
RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_MultiVector() const
Get the underlying Tpetra multivector.