Xpetra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Xpetra_BlockedVector_def.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 // Tobias Wiesner (tawiesn@sandia.gov)
42 // Ray Tuminaro (rstumin@sandia.gov)
43 //
44 // ***********************************************************************
45 //
46 // @HEADER
47 #ifndef XPETRA_BLOCKEDVECTOR_DEF_HPP
48 #define XPETRA_BLOCKEDVECTOR_DEF_HPP
49 
51 
52 #include "Xpetra_BlockedMultiVector.hpp"
53 #include "Xpetra_Exceptions.hpp"
54 
55 namespace Xpetra {
56 
57 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
59  BlockedVector(const Teuchos::RCP<const Xpetra::BlockedMap<LocalOrdinal, GlobalOrdinal, Node>>& map, bool zeroOut)
60  : Xpetra::BlockedMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>(map, 1, zeroOut) {}
61 
62 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
63 BlockedVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>::
64  BlockedVector(Teuchos::RCP<const Xpetra::BlockedMap<LocalOrdinal, GlobalOrdinal, Node>> bmap,
66  : Xpetra::BlockedMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>(bmap, v) {}
67 
68 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
72  : Xpetra::BlockedMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>(mapExtractor, v) {}
73 
74 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
77 
78 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
82  assign(rhs); // dispatch to protected virtual method
83  return *this;
84 }
85 
86 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
88  replaceGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar& value) {
89  BlockedMultiVector::replaceGlobalValue(globalRow, vectorIndex, value);
90 }
91 
92 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
94  sumIntoGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar& value) {
95  BlockedMultiVector::sumIntoGlobalValue(globalRow, vectorIndex, value);
96 }
97 
98 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
100  replaceLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar& value) {
101  BlockedMultiVector::replaceLocalValue(myRow, vectorIndex, value);
102 }
103 
104 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
106  sumIntoLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar& value) {
107  BlockedMultiVector::sumIntoLocalValue(myRow, vectorIndex, value);
108 }
109 
110 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
112  replaceGlobalValue(GlobalOrdinal globalRow, const Scalar& value) {
113  BlockedMultiVector::replaceGlobalValue(globalRow, 0, value);
114 }
115 
116 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
118  sumIntoGlobalValue(GlobalOrdinal globalRow, const Scalar& value) {
119  BlockedMultiVector::sumIntoGlobalValue(globalRow, 0, value);
120 }
121 
122 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
124  replaceLocalValue(LocalOrdinal myRow, const Scalar& value) {
125  BlockedMultiVector::replaceLocalValue(myRow, 0, value);
126 }
127 
128 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
130  sumIntoLocalValue(LocalOrdinal myRow, const Scalar& value) {
131  BlockedMultiVector::sumIntoLocalValue(myRow, 0, value);
132 }
133 
134 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
136  putScalar(const Scalar& value) {
138 }
139 
140 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
141 Teuchos::RCP<const Xpetra::Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
143  getVector(size_t j) const {
145 }
146 
147 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
148 Teuchos::RCP<Xpetra::Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
152 }
153 
154 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
155 Teuchos::ArrayRCP<const Scalar>
157  getData(size_t j) const {
158  return BlockedMultiVector::getData(j);
159 }
160 
161 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
162 Teuchos::ArrayRCP<Scalar>
164  getDataNonConst(size_t j) {
166 }
167 
168 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
170  dot(const MultiVector& A, const Teuchos::ArrayView<Scalar>& dots) const {
171  BlockedMultiVector::dot(A, dots);
172  return;
173 }
174 
175 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
176 Scalar
179  Teuchos::Array<Scalar> dots = Teuchos::Array<Scalar>(1);
180  BlockedMultiVector::dot(A, dots);
181  return dots[0];
182 }
183 
184 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
188  return;
189 }
190 
191 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
195  return;
196 }
197 
198 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
200  scale(const Scalar& alpha) {
202  return;
203 }
204 
205 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
207  scale(Teuchos::ArrayView<const Scalar> alpha) {
209  return;
210 }
211 
212 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
214  update(const Scalar& alpha,
216  const Scalar& beta) {
217  BlockedMultiVector::update(alpha, A, beta);
218  return;
219 }
220 
221 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
223  update(const Scalar& alpha,
225  const Scalar& beta,
227  const Scalar& gamma) {
228  BlockedMultiVector::update(alpha, A, beta, B, gamma);
229  return;
230 }
231 
232 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
233 typename Teuchos::ScalarTraits<Scalar>::magnitudeType
235  norm1() const {
236  using Array = Teuchos::Array<typename Teuchos::ScalarTraits<Scalar>::magnitudeType>;
237  Array norm = Array(1);
238  this->norm1(norm);
239  return norm[0];
240 }
241 
242 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
243 typename Teuchos::ScalarTraits<Scalar>::magnitudeType
245  norm2() const {
246  Teuchos::Array<typename Teuchos::ScalarTraits<Scalar>::magnitudeType> norm =
247  Teuchos::Array<typename Teuchos::ScalarTraits<Scalar>::magnitudeType>(1);
248  this->norm2(norm);
249  return norm[0];
250 }
251 
252 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
253 typename Teuchos::ScalarTraits<Scalar>::magnitudeType
255  normInf() const {
256  Teuchos::Array<typename Teuchos::ScalarTraits<Scalar>::magnitudeType>
257  norm = Teuchos::Array<typename Teuchos::ScalarTraits<Scalar>::magnitudeType>(1);
258  this->normInf(norm);
259  return norm[0];
260 }
261 
262 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
264  norm1(const Teuchos::ArrayView<typename Teuchos::ScalarTraits<Scalar>::magnitudeType>& norms) const {
266 }
267 
268 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
270  norm2(const Teuchos::ArrayView<typename Teuchos::ScalarTraits<Scalar>::magnitudeType>& norms) const {
272 }
273 
274 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
276  normInf(const Teuchos::ArrayView<typename Teuchos::ScalarTraits<Scalar>::magnitudeType>& norms) const {
278 }
279 
280 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
282  meanValue(const Teuchos::ArrayView<Scalar>& /* means */) const {
283  throw Xpetra::Exceptions::RuntimeError("BlockedVector::meanValue: Not (yet) supported by BlockedVector.");
284 }
285 
286 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
287 Scalar
289  meanValue() const {
290  throw Xpetra::Exceptions::RuntimeError("BlockedVector::meanValue: Not (yet) supported by BlockedVector.");
291 }
292 
293 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
295  multiply(Teuchos::ETransp /* transA */,
296  Teuchos::ETransp /* transB */,
297  const Scalar& /* alpha */,
300  const Scalar& /* beta */) {
301  throw Xpetra::Exceptions::RuntimeError("BlockedVector::multiply: Not (yet) supported by BlockedVector.");
302 }
303 
304 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
306  multiply(Teuchos::ETransp /* transA */,
307  Teuchos::ETransp /* transB */,
308  const Scalar& /* alpha */,
311  const Scalar& /* beta */) {
312  throw Xpetra::Exceptions::RuntimeError("BlockedVector::multiply: Not (yet) supported by BlockedVector.");
313 }
314 
315 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
317  elementWiseMultiply(Scalar /* scalarAB */,
320  Scalar /* scalarThis */) {
321  throw Xpetra::Exceptions::RuntimeError("BlockedVector::elementWiseMultiply: Not (yet) supported by BlockedVector.");
322 }
323 
324 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
326  elementWiseMultiply(Scalar /* scalarAB */,
329  Scalar /* scalarThis */) {
330  XPETRA_TEST_FOR_EXCEPTION(B.getMap()->isSameAs(*(this->getMap())) == false,
332  "BlockedVector::elementWiseMultipy: B must have same blocked map than this.");
333  TEUCHOS_TEST_FOR_EXCEPTION(A.getMap()->getLocalNumElements() != B.getMap()->getLocalNumElements(),
335  "BlockedVector::elementWiseMultipy: A has "
336  << A.getMap()->getLocalNumElements() << " elements, B has " << B.getMap()->getLocalNumElements()
337  << ".");
338  TEUCHOS_TEST_FOR_EXCEPTION(A.getMap()->getGlobalNumElements() != B.getMap()->getGlobalNumElements(),
340  "BlockedVector::elementWiseMultipy: A has " << A.getMap()->getGlobalNumElements()
341  << " elements, B has "
342  << B.getMap()->getGlobalNumElements() << ".");
343 
344  RCP<const BlockedMap> bmap = this->getBlockedMap();
345  RCP<const Xpetra::Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node>> rcpA = Teuchos::rcpFromRef(A);
346  RCP<const Xpetra::Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node>> bmvec = Teuchos::rcpFromRef(B);
347  RCP<const BlockedVector> bbmvec = Teuchos::rcp_dynamic_cast<const BlockedVector>(bmvec);
348  TEUCHOS_TEST_FOR_EXCEPTION(bbmvec.is_null() == true,
350  "BlockedVector::elementWiseMultipy: B must be a BlockedVector.");
351 
352  // TODO implement me
353  /*RCP<Xpetra::MapExtractor<Scalar,LocalOrdinal,GlobalOrdinal,Node> > me = Teuchos::rcp(new
354  Xpetra::MapExtractor<Scalar,LocalOrdinal,GlobalOrdinal,Node>(bmap));
355 
356  for(size_t m = 0; m < bmap->getNumMaps(); m++) {
357  // TODO introduce BlockedVector objects and "skip" this expensive ExtractVector call
358  RCP<const Xpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > pd = me->ExtractVector(rcpA,m,bmap->getThyraMode());
359  XPETRA_TEST_FOR_EXCEPTION(pd->getMap()->isSameAs(*(this->getBlockedMap()->getMap(m,bmap->getThyraMode())))==false,
360  Xpetra::Exceptions::RuntimeError, "BlockedVector::elementWiseMultipy: sub map of B does not fit with sub map of this.");
361  this->getMultiVector(m,bmap->getThyraMode())->elementWiseMultiply(scalarAB,*pd,*(bbmvec->getMultiVector(m,bmap->getThyraMode())),scalarThis);
362  }*/
363 }
364 
365 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
366 size_t
368  getNumVectors() const {
369  return 1;
370 }
371 
372 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
373 size_t
377  "BlockedVector::getLocalLength: routine not implemented. It has no value as one must iterate on the partial vectors.");
378  TEUCHOS_UNREACHABLE_RETURN(0);
379 }
380 
381 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
385  return this->getBlockedMap()->getFullMap()->getGlobalNumElements();
386 }
387 
388 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
392  "BlockedVector::isSameSize: routine not implemented. It has no value as one must iterate on the partial vectors.");
393  TEUCHOS_UNREACHABLE_RETURN(0);
394 }
395 
396 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
397 std::string
399  description() const {
400  return std::string("BlockedVector");
401 }
402 
403 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
405  describe(Teuchos::FancyOStream& out, const Teuchos::EVerbosityLevel verbLevel) const {
406  out << description() << std::endl;
407  for (size_t r = 0; r < this->getBlockedMap()->getNumMaps(); r++) {
408  getMultiVector(r)->describe(out, verbLevel);
409  }
410 }
411 
412 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
414  replaceMap(const RCP<const Map>& map) {
416 }
417 
418 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
421  const Import& /* importer */,
422  CombineMode /* CM */) {
423  throw Xpetra::Exceptions::RuntimeError("BlockedVector::doImport: Not supported by BlockedVector.");
424 }
425 
426 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
429  const Import& /* importer */,
430  CombineMode /* CM */) {
431  throw Xpetra::Exceptions::RuntimeError("BlockedVector::doExport: Not supported by BlockedVector.");
432 }
433 
434 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
437  const Export& /* exporter */,
438  CombineMode /* CM */) {
439  throw Xpetra::Exceptions::RuntimeError("BlockedVector::doImport: Not supported by BlockedVector.");
440 }
441 
442 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
445  const Export& /* exporter */,
446  CombineMode /* CM */) {
447  throw Xpetra::Exceptions::RuntimeError("BlockedVector::doExport: Not supported by BlockedVector.");
448 }
449 
450 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
452  setSeed(unsigned int seed) {
453  for (size_t r = 0; r < this->getBlockedMap()->getNumMaps(); ++r) {
454  getMultiVector(r)->setSeed(seed);
455  }
456 }
457 
458 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
460  randomize(bool bUseXpetraImplementation) {
461  for (size_t r = 0; r < this->getBlockedMap()->getNumMaps(); ++r) {
462  getMultiVector(r)->randomize(bUseXpetraImplementation);
463  }
464 }
465 
466 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
468  randomize(const Scalar& minVal, const Scalar& maxVal, bool bUseXpetraImplementation) {
469  for (size_t r = 0; r < this->getBlockedMap()->getNumMaps(); ++r) {
470  getMultiVector(r)->randomize(minVal, maxVal, bUseXpetraImplementation);
471  }
472 }
473 
474 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
477  {
479  }
480 }
481 
482 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
484  Xpetra_randomize(const Scalar& minVal, const Scalar& maxVal) {
485  {
487  }
488 }
489 
490 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
491 Teuchos::RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>>
493  getMap() const {
494  XPETRA_MONITOR("BlockedVector::getMap");
495  return this->getBlockedMap();
496 }
497 
498 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
499 Teuchos::RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
501  getMultiVector(size_t r) const {
503 }
504 
505 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
506 Teuchos::RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
508  getMultiVector(size_t r, bool bThyraMode) const {
509  return BlockedMultiVector::getMultiVector(r, bThyraMode);
510 }
511 
512 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
516  bool bThyraMode) {
517  BlockedMultiVector::setMultiVector(r, v, bThyraMode);
518  return;
519 }
520 
521 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
522 Teuchos::RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
524  Merge() const {
525  return BlockedMultiVector::Merge();
526 }
527 
528 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
532 }
533 
534 // template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
535 // virtual void BlockedVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::
536 // assign (const XpetrA::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>& rhs)
537 // {
538 // throw Xpetra::Exceptions::RuntimeError("BlockedVector::assign: Not supported by BlockedVector.");
539 // }
540 
541 } // namespace Xpetra
542 
543 #endif // XPETRA_BLOCKEDVECTOR_DEF_HPP
virtual Teuchos::RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getVectorNonConst(size_t j)
Return a Vector which is a nonconst view of column j.
virtual void reciprocal(const MultiVector &)
Put element-wise reciprocal values of input Multi-vector in target, this(i,j) = 1/A(i,j).
virtual Scalar meanValue() const
Compute mean (average) value of this Vector.
void setMultiVector(size_t r, Teuchos::RCP< const Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node >> v, bool bThyraMode)
set partial Vector associated with block row r
virtual void setSeed(unsigned int seed)
Set seed for Random function.
virtual void multiply(Teuchos::ETransp, Teuchos::ETransp, const Scalar &, const Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &, const Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &, const Scalar &)
Matrix-matrix multiplication: this = beta*this + alpha*op(A)*op(B).
virtual void sumIntoLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value)
Add value to existing value, using local (row) index.
virtual void update(const Scalar &alpha, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta)
Update multi-vector values with scaled values of A, this = beta*this + alpha*A.
virtual Teuchos::RCP< const Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getVector(size_t j) const
Return a Vector which is a const view of column j.
virtual void sumIntoGlobalValue(GlobalOrdinal, size_t, const Scalar &)
Add value to existing value, using global (row) index.
virtual Teuchos::RCP< const Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getVector(size_t j) const
Return a Vector which is a const view of column j.
virtual void elementWiseMultiply(Scalar, const Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &, Scalar)
Element-wise multiply of a Vector A with a MultiVector B.
virtual void scale(const Scalar &alpha)
Scale the current values of a vector, this = alpha*this.
virtual void replaceMap(const RCP< const Map > &map)
virtual size_t getNumVectors() const
Number of columns in the Vector.
virtual Teuchos::RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getVectorNonConst(size_t j)
Return a Vector which is a nonconst view of column j.
virtual void replaceGlobalValue(GlobalOrdinal, size_t, const Scalar &)
Replace value, using global (row) index.
virtual std::string description() const
A simple one-line description of this object.
Exception throws to report errors in the internal logical of the program.
virtual void Xpetra_randomize()
Set multi-vector values to random numbers. XPetra implementation.
virtual void abs(const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
Put element-wise absolute values of input vector in target: A = abs(this).
virtual void sumIntoLocalValue(LocalOrdinal, size_t, const Scalar &)
Add value to existing value, using local (row) index.
virtual void dot(const Xpetra::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]).
virtual void update(const Scalar &alpha, const MultiVector &A, const Scalar &beta)
Update multi-vector values with scaled values of A, this = beta*this + alpha*A.
virtual void replaceLocalValue(LocalOrdinal, size_t, const Scalar &)
Replace value, using local (row) index.
virtual 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::RCP< MultiVector > Merge() const
merge BlockedMultiVector blocks to a single MultiVector
virtual void reciprocal(const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
Put element-wise reciprocal values of input vector in target, this(i,j) = 1/A(i,j).
virtual void norm2(const Teuchos::ArrayView< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
virtual bool isSameSize(const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &) const
Local number of rows on the calling process.
virtual void assign(const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &rhs)
Implementation of the assignment operator (operator=); does a deep copy.
virtual void replaceMap(const RCP< const Map > &map)
Teuchos::RCP< const Map > getMap() const
Access function for the underlying Map this DistObject was constructed with.
virtual size_t getLocalLength() const
Local number of rows on the calling process.
virtual Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const =0
The Map describing the parallel distribution of this object.
virtual void randomize(bool bUseXpetraImplementation=false)
Set multi-vector values to random numbers.
virtual void scale(const Scalar &alpha)
Scale the current values of a multi-vector, this = alpha*this.
virtual ~BlockedVector()
Destructor.
virtual Teuchos::ScalarTraits< Scalar >::magnitudeType norm2() const
Compute 2-norm of vector.
void setMultiVector(size_t r, Teuchos::RCP< const MultiVector > v, bool bThyraMode)
set partial multivector associated with block row r
virtual void norm1(const Teuchos::ArrayView< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
Compute 1-norm of each vector in multi-vector.
virtual void normInf(const Teuchos::ArrayView< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
Compute Inf-norm of each vector in multi-vector.
size_t global_size_t
Global size_t object.
virtual Teuchos::ArrayRCP< const Scalar > getData(size_t j) const
Const view of the local values in a particular vector of this multivector.
virtual void abs(const MultiVector &)
Put element-wise absolute values of input Multi-vector in target: A = abs(this).
BlockedVector(const Teuchos::RCP< const BlockedMap > &map, bool zeroOut=true)
Constructor.
virtual void dot(const MultiVector &, const Teuchos::ArrayView< Scalar > &) const
Compute dot product of each corresponding pair of vectors, dots[i] = this[i].dot(A[i]).
Teuchos::RCP< MultiVector > getMultiVector(size_t r) const
return partial multivector associated with block row r
virtual void assign(const MultiVector &rhs)
Implementation of the assignment operator (operator=); does a deep copy.
Teuchos::RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Merge() const
merge BlockedVector blocks to a single Vector
virtual Teuchos::ArrayRCP< Scalar > getDataNonConst(size_t j)
View of the local values in a particular vector of this multivector.
virtual Teuchos::ArrayRCP< const Scalar > getData(size_t j) const
Const view of the local values in a particular vector of this vector.
#define XPETRA_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
virtual Teuchos::ScalarTraits< Scalar >::magnitudeType norm1() const
Compute 1-norm of vector.
virtual void putScalar(const Scalar &value)
Set all values in the multivector with the given value.
virtual Teuchos::ArrayRCP< Scalar > getDataNonConst(size_t j)
View of the local values in a particular vector of this vector.
virtual global_size_t getGlobalLength() const
Global number of rows in the Vector.
virtual void replaceGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value)
Replace value, using global (row) index.
BlockedVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & operator=(const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &rhs)
Assignment operator: Does a deep copy.
CombineMode
Xpetra::Combine Mode enumerable type.
#define XPETRA_MONITOR(funcName)
virtual void sumIntoGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value)
Add value to existing value, using global (row) index.
virtual Teuchos::ScalarTraits< Scalar >::magnitudeType normInf() const
Compute Inf-norm in vector.
virtual void doImport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &, const Import &, CombineMode)
Import.
virtual void doExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &, const Import &, CombineMode)
Export.
virtual void putScalar(const Scalar &value)
Set all values in the vector with the given value.
virtual void Xpetra_randomize()
Set vector values to random numbers. XPetra implementation.
Teuchos::RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getMultiVector(size_t r) const
return partial Vector associated with block row r
virtual void replaceLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value)
Replace value, using local (row) index.