Xpetra Version of the Day
Loading...
Searching...
No Matches
Xpetra_EpetraMultiVector.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_EPETRAMULTIVECTOR_HPP
47#define XPETRA_EPETRAMULTIVECTOR_HPP
48
49/* this file is automatically generated - do not edit (see script/epetra.py) */
50
51#include <Kokkos_Core.hpp>
52#include <Kokkos_DualView.hpp>
53
54
56
57#include "Xpetra_MultiVector.hpp"
58#include "Xpetra_Vector.hpp"
59
60#include "Xpetra_EpetraMap.hpp"
62#include "Xpetra_Utils.hpp"
65#include "Xpetra_Exceptions.hpp"
66#include "Epetra_SerialComm.h"
67
68#include <Epetra_MultiVector.h>
69#include <Epetra_Vector.h>
70
71namespace Xpetra {
72
73 // TODO: move that elsewhere
74 template<class GlobalOrdinal, class Node>
75 const Epetra_MultiVector & toEpetra(const MultiVector<double,int,GlobalOrdinal,Node> &);
76 template<class GlobalOrdinal, class Node>
77 Epetra_MultiVector & toEpetra(MultiVector<double, int,GlobalOrdinal,Node> &);
78 template<class GlobalOrdinal, class Node>
79 RCP<MultiVector<double, int, GlobalOrdinal, Node> > toXpetra(RCP<Epetra_MultiVector> vec);
80
81 // we need this forward declaration
82#ifndef DOXYGEN_SHOULD_SKIP_THIS
83 template<class GlobalOrdinal, class Node> class EpetraVectorT;
84#endif
85
86 template<class EpetraGlobalOrdinal, class Node>
88 : public virtual MultiVector<double, int, EpetraGlobalOrdinal, Node>
89 {
90 typedef double Scalar;
91 typedef int LocalOrdinal;
92 typedef EpetraGlobalOrdinal GlobalOrdinal;
93
94 public:
95
97
98
100 EpetraMultiVectorT(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, size_t NumVectors, bool zeroOut=true) {
102 "Xpetra::EpetraMultiVector only available for GO=int or GO=long long with EpetraNode (Serial or OpenMP depending on configuration)");
103 }
104
108 "Xpetra::EpetraMultiVector only available for GO=int or GO=long long with EpetraNode (Serial or OpenMP depending on configuration)");
109 }
110
114 "Xpetra::EpetraMultiVector only available for GO=int or GO=long long with EpetraNode (Serial or OpenMP depending on configuration)");
115 }
116
119
121
123
124
126 void replaceGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value) { }
127
129 void sumIntoGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value) { }
130
132 void replaceLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value) { }
133
135 void sumIntoLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value) { }
136
138 void putScalar(const Scalar &value) { }
139
141
143
144
147 return Teuchos::null;
148 }
149
152 return Teuchos::null;
153 }
154
157 return ArrayRCP<const Scalar>();
158 }
159
162 return ArrayRCP<Scalar>();
163 }
164
166
168
169
172
175
178
180 void scale(const Scalar &alpha) { }
181
184
187
190
193
196
199
201 void meanValue(const Teuchos::ArrayView< Scalar > &means) const { }
202
205
208
210
212
213
215 size_t getNumVectors() const { return 0; }
216
218 size_t getLocalLength() const { return 0; }
219
221 global_size_t getGlobalLength() const { return 0; }
222
223 // \brief Checks to see if the local length, number of vectors and size of Scalar type match
224 bool isSameSize(const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> & vec) const { return false; }
225
227
229
230
232 std::string description() const { return std::string(""); }
233
236
238
240 void randomize(bool bUseXpetraImplementation = false) { }
241
243 void randomize(const Scalar& minVal, const Scalar& maxVal, bool bUseXpetraImplementation = false) { }
244
246 //{@
247
250
253
256
259
262
265
267
269
270
272 EpetraMultiVectorT(const RCP<Epetra_MultiVector> &vec) { //TODO removed const
274 "Xpetra::EpetraMultiVector only available for GO=int or GO=long long with EpetraNode (Serial or OpenMP depending on configuration)");
275 }
276
278 RCP<Epetra_MultiVector> getEpetra_MultiVector() const { return Teuchos::null; }
279
281 void setSeed(unsigned int seed) { }
282
284
285 protected:
288 virtual void
290
291 }; // EpetraMultiVectorT class
292
293 // specialization on GO=int and Node=EpetraNode
294#ifndef XPETRA_EPETRA_NO_32BIT_GLOBAL_INDICES
295 template<>
297 : public virtual MultiVector<double, int, int, EpetraNode>
298 {
299 typedef double Scalar;
300 typedef int LocalOrdinal;
301 typedef int GlobalOrdinal;
303
304 public:
305
307
308
310 EpetraMultiVectorT(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, size_t NumVectors, bool zeroOut=true)
311 : vec_(Teuchos::rcp(new Epetra_MultiVector(toEpetra<GlobalOrdinal,Node>(map), Teuchos::as<int>(NumVectors), zeroOut))) { }
312
315 if (copyOrView == Teuchos::Copy)
316 vec_ = Teuchos::rcp(new Epetra_MultiVector(toEpetra<GlobalOrdinal,Node>(source)));
317 else {
318 int* indices = new int [source.getNumVectors()];
319 for (size_t i = 0; i < source.getNumVectors(); i++)
320 indices[i] = i;
321 vec_ = Teuchos::rcp(new Epetra_MultiVector(View, toEpetra<GlobalOrdinal,Node>(source), indices, source.getNumVectors()));
322 delete [] indices;
323 }
324 }
325
328 //TODO: input argument 'NumVectors' is not necessary in both Xpetra and Tpetra interface. Should it be removed?
329
330 const std::string tfecfFuncName("MultiVector(ArrayOfPtrs)");
331 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(NumVectors < 1 || NumVectors != Teuchos::as<size_t>(ArrayOfPtrs.size()), std::runtime_error,
332 ": ArrayOfPtrs.size() must be strictly positive and as large as ArrayOfPtrs.");
333
334 #ifdef HAVE_XPETRA_DEBUG
335 // This cannot be tested by Epetra itself
336 {
337 size_t localLength = map->getLocalNumElements();
338 for(int j=0; j<ArrayOfPtrs.size(); j++) {
339 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(Teuchos::as<size_t>(ArrayOfPtrs[j].size()) != localLength, std::runtime_error,
340 ": ArrayOfPtrs[" << j << "].size() (== " << ArrayOfPtrs[j].size() <<
341 ") is not equal to getLocalLength() (== " << localLength);
342
343 }
344 }
345 #endif
346
347 // Convert Teuchos::ArrayView< const Teuchos::ArrayView< const Scalar > > to double**
348 Array<const double*> arrayOfRawPtrs(ArrayOfPtrs.size());
349 for(int i=0; i<ArrayOfPtrs.size(); i++) {
350 arrayOfRawPtrs[i] = ArrayOfPtrs[i].getRawPtr();
351 }
352 double** rawArrayOfRawPtrs = const_cast<double**>(arrayOfRawPtrs.getRawPtr()); // This const_cast should be fine, because Epetra_DataAccess=Copy.
353
354 vec_ = Teuchos::rcp(new Epetra_MultiVector(Copy, toEpetra<GlobalOrdinal,Node>(map), rawArrayOfRawPtrs, static_cast<int>(NumVectors)));
355 }
356
359
361
363
364
366 void replaceGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value) { XPETRA_MONITOR("EpetraMultiVectorT::replaceGlobalValue"); vec_->ReplaceGlobalValue(globalRow, Teuchos::as<int>(vectorIndex), value); }
367
369 void sumIntoGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value) { XPETRA_MONITOR("EpetraMultiVectorT::sumIntoGlobalValue"); vec_->SumIntoGlobalValue(globalRow, Teuchos::as<int>(vectorIndex), value); }
370
372 void replaceLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value) { XPETRA_MONITOR("EpetraMultiVectorT::replaceLocalValue"); vec_->ReplaceMyValue(myRow, Teuchos::as<int>(vectorIndex), value); }
373
375 void sumIntoLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value) { XPETRA_MONITOR("EpetraMultiVectorT::sumIntoLocalValue"); vec_->SumIntoMyValue(myRow, Teuchos::as<int>(vectorIndex), value); }
376
378 void putScalar(const Scalar &value) { XPETRA_MONITOR("EpetraMultiVectorT::putScalar"); vec_->PutScalar(value); }
379
381
383
384
387
390
393 XPETRA_MONITOR("EpetraMultiVectorT::getData");
394
395 double ** arrayOfPointers;
396
397 vec_->ExtractView(&arrayOfPointers);
398
399 double * data = arrayOfPointers[j];
400 int localLength = vec_->MyLength();
401
402 return ArrayRCP<double>(data, 0, localLength, false); // no ownership
403 }
404
407 XPETRA_MONITOR("EpetraMultiVectorT::getDataNonConst");
408
409 double ** arrayOfPointers;
410
411 vec_->ExtractView(&arrayOfPointers);
412
413 double * data = arrayOfPointers[j];
414 int localLength = vec_->MyLength();
415
416 return ArrayRCP<double>(data, 0, localLength, false); // no ownership
417 }
418
420
422
423
426 XPETRA_MONITOR("EpetraMultiVectorT::dot");
427
428 XPETRA_DYNAMIC_CAST(const EpetraMultiVectorT<GlobalOrdinal XPETRA_COMMA Node>, A, eA, "This Xpetra::EpetraMultiVectorT method only accept Xpetra::EpetraMultiVectorT as input arguments.");
429 vec_->Dot(*eA.getEpetra_MultiVector(), dots.getRawPtr());
430 }
431
433 void abs(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A) { XPETRA_MONITOR("EpetraMultiVectorT::abs"); vec_->Abs(toEpetra<GlobalOrdinal,Node>(A)); }
434
436 void reciprocal(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A) { XPETRA_MONITOR("EpetraMultiVectorT::reciprocal"); vec_->Reciprocal(toEpetra<GlobalOrdinal,Node>(A)); }
437
439 void scale(const Scalar &alpha) { XPETRA_MONITOR("EpetraMultiVectorT::scale"); vec_->Scale(alpha); }
440
443 XPETRA_MONITOR("EpetraMultiVectorT::scale");
444 // Epetra, unlike Tpetra, doesn't implement this version of
445 // scale(). Deal with this by scaling one column at a time.
446 const size_t numVecs = this->getNumVectors ();
447 for (size_t j = 0; j < numVecs; ++j) {
448 Epetra_Vector *v = (*vec_)(j);
449 v->Scale (alpha[j]);
450 }
451 }
452
454 void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta) { XPETRA_MONITOR("EpetraMultiVectorT::update"); vec_->Update(alpha, toEpetra<GlobalOrdinal,Node>(A), beta); }
455
457 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("EpetraMultiVectorT::update"); vec_->Update(alpha, toEpetra<GlobalOrdinal,Node>(A), beta, toEpetra<GlobalOrdinal,Node>(B), gamma); }
458
460 void norm1(const Teuchos::ArrayView< Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const { XPETRA_MONITOR("EpetraMultiVectorT::norm1"); vec_->Norm1(norms.getRawPtr()); }
461
463 void norm2(const Teuchos::ArrayView< Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const { XPETRA_MONITOR("EpetraMultiVectorT::norm2"); vec_->Norm2(norms.getRawPtr()); }
464
466 void normInf(const Teuchos::ArrayView< Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const { XPETRA_MONITOR("EpetraMultiVectorT::normInf"); vec_->NormInf(norms.getRawPtr()); }
467
469 void meanValue(const Teuchos::ArrayView< Scalar > &means) const { XPETRA_MONITOR("EpetraMultiVectorT::meanValue"); vec_->MeanValue(means.getRawPtr()); } //TODO: modify ArrayView size ??
470
472 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("EpetraMultiVectorT::multiply"); vec_->Multiply(toEpetra(transA), toEpetra(transB), alpha, toEpetra(A), toEpetra(B), beta); }
473
475 void elementWiseMultiply(Scalar scalarAB, const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, Scalar scalarThis) { XPETRA_MONITOR("EpetraMultiVectorT::elementWiseMultiply"); vec_->Multiply(scalarAB, toEpetra<GlobalOrdinal,Node>(A), toEpetra<GlobalOrdinal,Node>(B), scalarThis); }
476
478
480
481
483 size_t getNumVectors() const { XPETRA_MONITOR("EpetraMultiVectorT::getNumVectors"); return vec_->NumVectors(); }
484
486 size_t getLocalLength() const { XPETRA_MONITOR("EpetraMultiVectorT::getLocalLength"); return vec_->MyLength(); }
487
489 global_size_t getGlobalLength() const { XPETRA_MONITOR("EpetraMultiVectorT::getGlobalLength"); return vec_->GlobalLength64(); }
490
493 XPETRA_MONITOR("EpetraMultiVectorT::isSameSize");
494 auto vv = toEpetra<GlobalOrdinal,Node>(vec);
495 return ( (vec_->MyLength() == vv.MyLength()) && (vec_->NumVectors() == vv.NumVectors()));
496 }
497
499
501
502
504 std::string description() const {
505 XPETRA_MONITOR("EpetraMultiVectorT::description");
508 }
509
512 XPETRA_MONITOR("EpetraMultiVectorT::describe");
513 vec_->Print(out);
514 }
515
517
519 void randomize(bool bUseXpetraImplementation = false) {
520 XPETRA_MONITOR("EpetraMultiVectorT::randomize");
521
522 if (bUseXpetraImplementation)
524 else
525 vec_->Random();
526 }
527
529 void randomize(const Scalar& minVal, const Scalar& maxVal, bool bUseXpetraImplementation = false) {
530 XPETRA_MONITOR("EpetraMultiVectorT::randomize");
531
532 if (bUseXpetraImplementation)
534 else {
535 vec_->Random();
536 const size_t numVectors = getNumVectors();
537 for(size_t i = 0; i < numVectors; i++)
538 {
540
541 const size_t myLength = getLocalLength();
542 for(size_t j = 0; j < myLength; j++)
543 {
544 datai[ j ] = 0.5*(maxVal-minVal)*datai[ j ]+0.5*(maxVal+minVal);
545 }
546 }
547 }
548 }
549
551 //{@
552
554 Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const { XPETRA_MONITOR("EpetraMultiVectorT::getMap"); return toXpetra<GlobalOrdinal,Node>(vec_->Map()); }
555
558 XPETRA_MONITOR("EpetraMultiVectorT::doImport");
559
560 XPETRA_DYNAMIC_CAST(const EpetraMultiVectorT<GlobalOrdinal XPETRA_COMMA Node>, source, tSource, "Xpetra::EpetraMultiVectorT::doImport only accept Xpetra::EpetraMultiVectorT as input arguments.");
561 XPETRA_DYNAMIC_CAST(const EpetraImportT<GlobalOrdinal XPETRA_COMMA Node>, importer, tImporter, "Xpetra::EpetraMultiVectorT::doImport only accept Xpetra::EpetraImportT as input arguments.");
562
563 RCP<Epetra_MultiVector> v = tSource.getEpetra_MultiVector();
564 int err = this->getEpetra_MultiVector()->Import(*v, *tImporter.getEpetra_Import(), toEpetra(CM));
565 TEUCHOS_TEST_FOR_EXCEPTION(err != 0, std::runtime_error, "Catch error code returned by Epetra is " << err);
566 }
567
570 XPETRA_MONITOR("EpetraMultiVectorT::doExport");
571
572 XPETRA_DYNAMIC_CAST(const EpetraMultiVectorT<GlobalOrdinal XPETRA_COMMA Node>, dest, tDest, "Xpetra::EpetraMultiVectorT::doImport only accept Xpetra::EpetraMultiVectorT as input arguments.");
573 XPETRA_DYNAMIC_CAST(const EpetraImportT<GlobalOrdinal XPETRA_COMMA Node>, importer, tImporter, "Xpetra::EpetraMultiVectorT::doImport only accept Xpetra::EpetraImportT as input arguments.");
574
575 RCP<Epetra_MultiVector> v = tDest.getEpetra_MultiVector();
576 int err = this->getEpetra_MultiVector()->Export(*v, *tImporter.getEpetra_Import(), toEpetra(CM));
577 TEUCHOS_TEST_FOR_EXCEPTION(err != 0, std::runtime_error, "Catch error code returned by Epetra.");
578 }
579
582 XPETRA_MONITOR("EpetraMultiVectorT::doImport");
583
584 XPETRA_DYNAMIC_CAST(const EpetraMultiVectorT<GlobalOrdinal XPETRA_COMMA Node>, source, tSource, "Xpetra::EpetraMultiVectorT::doImport only accept Xpetra::EpetraMultiVectorT as input arguments.");
585 XPETRA_DYNAMIC_CAST(const EpetraExportT<GlobalOrdinal XPETRA_COMMA Node>, exporter, tExporter, "Xpetra::EpetraMultiVectorT::doImport only accept Xpetra::EpetraImportT as input arguments.");
586
587 RCP<Epetra_MultiVector> v = tSource.getEpetra_MultiVector();
588 int err = this->getEpetra_MultiVector()->Import(*v, *tExporter.getEpetra_Export(), toEpetra(CM));
589 TEUCHOS_TEST_FOR_EXCEPTION(err != 0, std::runtime_error, "Catch error code returned by Epetra.");
590 }
591
594 XPETRA_MONITOR("EpetraMultiVectorT::doExport");
595
596 XPETRA_DYNAMIC_CAST(const EpetraMultiVectorT<GlobalOrdinal XPETRA_COMMA Node>, dest, tDest, "Xpetra::EpetraMultiVectorT::doImport only accept Xpetra::EpetraMultiVectorT as input arguments.");
597 XPETRA_DYNAMIC_CAST(const EpetraExportT<GlobalOrdinal XPETRA_COMMA Node>, exporter, tExporter, "Xpetra::EpetraMultiVectorT::doImport only accept Xpetra::EpetraImportT as input arguments.");
598
599 RCP<Epetra_MultiVector> v = tDest.getEpetra_MultiVector();
600 int err = this->getEpetra_MultiVector()->Export(*v, *tExporter.getEpetra_Export(), toEpetra(CM));
601 TEUCHOS_TEST_FOR_EXCEPTION(err != 0, std::runtime_error, "Catch error code returned by Epetra.");
602 }
603
606 XPETRA_MONITOR("EpetraMultiVectorT::replaceMap");
607 int err = 0;
608 if (!map.is_null()) {
609 err = this->getEpetra_MultiVector()->ReplaceMap(toEpetra<GlobalOrdinal,Node>(map));
610
611 } else {
612 // Replace map with a dummy map to avoid potential hangs later
613 Epetra_SerialComm SComm;
614 Epetra_Map NewMap((GlobalOrdinal) vec_->MyLength(), (GlobalOrdinal) vec_->Map().IndexBase64(), SComm);
615 err = this->getEpetra_MultiVector()->ReplaceMap(NewMap);
616 }
617 TEUCHOS_TEST_FOR_EXCEPTION(err != 0, std::runtime_error, "Catch error code returned by Epetra.");
618 }
619
621
623
624
626 EpetraMultiVectorT(const RCP<Epetra_MultiVector> &vec) : vec_(vec) { } //TODO removed const
627
630
632 void setSeed(unsigned int seed) {
633 XPETRA_MONITOR("EpetraMultiVectorT::seedrandom");
634
636 vec_->SetSeed(seed);
637 }
638
640
641 typename dual_view_type::t_host_const_um getHostLocalView (Access::ReadOnlyStruct) const override { return getHostLocalView(Access::ReadWrite);}
642
643 typename dual_view_type::t_dev_const_um getDeviceLocalView(Access::ReadOnlyStruct) const override { return getDeviceLocalView(Access::ReadWrite);}
644
645 typename dual_view_type::t_host_um getHostLocalView (Access::OverwriteAllStruct) const override { return getHostLocalView(Access::ReadWrite);}
646
647 typename dual_view_type::t_dev_um getDeviceLocalView(Access::OverwriteAllStruct) const override { return getDeviceLocalView(Access::ReadWrite);}
648
649 typename dual_view_type::t_host_um getHostLocalView (Access::ReadWriteStruct) const override {
650 typedef Kokkos::View< typename dual_view_type::t_host::data_type ,
651 Kokkos::LayoutLeft,
652 typename dual_view_type::t_host::device_type ,
653 Kokkos::MemoryUnmanaged> epetra_view_type;
654
655 // access Epetra multivector data
656 double* data = NULL;
657 int myLDA;
658 vec_->ExtractView(&data, &myLDA);
659 int localLength = vec_->MyLength();
660 int numVectors = getNumVectors();
661
662 // create view
663 epetra_view_type test = epetra_view_type(data, localLength, numVectors);
664 typename dual_view_type::t_host_um ret = subview(test, Kokkos::ALL(), Kokkos::ALL());
665
666 return ret;
667 }
668
669 typename dual_view_type::t_dev_um getDeviceLocalView(Access::ReadWriteStruct) const override { return getHostLocalView(Access::ReadWrite); }
670
672
673 protected:
676 virtual void
678 typedef EpetraMultiVectorT this_type;
679 const this_type* rhsPtr = dynamic_cast<const this_type*> (&rhs);
681 rhsPtr == NULL, std::invalid_argument, "Xpetra::MultiVector::operator=: "
682 "The left-hand side (LHS) of the assignment has a different type than "
683 "the right-hand side (RHS). The LHS has type Xpetra::EpetraMultiVectorT "
684 "(which means it wraps an Epetra_MultiVector), but the RHS has some "
685 "other type. This probably means that the RHS wraps a Tpetra::Multi"
686 "Vector. Xpetra::MultiVector does not currently implement assignment "
687 "from a Tpetra object to an Epetra object, though this could be added "
688 "with sufficient interest.");
689
690 RCP<const Epetra_MultiVector> rhsImpl = rhsPtr->getEpetra_MultiVector ();
692
694 rhsImpl.is_null (), std::logic_error, "Xpetra::MultiVector::operator= "
695 "(in Xpetra::EpetraMultiVectorT::assign): *this (the right-hand side of "
696 "the assignment) has a null RCP<Epetra_MultiVector> inside. Please "
697 "report this bug to the Xpetra developers.");
699 lhsImpl.is_null (), std::logic_error, "Xpetra::MultiVector::operator= "
700 "(in Xpetra::EpetraMultiVectorT::assign): The left-hand side of the "
701 "assignment has a null RCP<Epetra_MultiVector> inside. Please report "
702 "this bug to the Xpetra developers.");
703
704 // Epetra_MultiVector's assignment operator does a deep copy.
705 *lhsImpl = *rhsImpl;
706 }
707
708 private:
711
712 }; // EpetraMultiVectorT class (specialization on GO=int, NO=EpetraNode
713#endif
714
715 // specialization on GO=long long and EpetraNode
716#ifndef XPETRA_EPETRA_NO_64BIT_GLOBAL_INDICES
717 template<>
719 : public virtual MultiVector<double, int, long long, EpetraNode>
720 {
721 typedef double Scalar;
722 typedef int LocalOrdinal;
723 typedef long long GlobalOrdinal;
725
726 public:
727
729
730
732 EpetraMultiVectorT(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, size_t NumVectors, bool zeroOut=true)
733 : vec_(Teuchos::rcp(new Epetra_MultiVector(toEpetra<GlobalOrdinal,Node>(map), Teuchos::as<int>(NumVectors), zeroOut))) { }
734
737 : vec_(Teuchos::rcp(new Epetra_MultiVector(toEpetra<GlobalOrdinal,Node>(source)))) { }
738
741 //TODO: input argument 'NumVectors' is not necessary in both Xpetra and Tpetra interface. Should it be removed?
742
743 const std::string tfecfFuncName("MultiVector(ArrayOfPtrs)");
744 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(NumVectors < 1 || NumVectors != Teuchos::as<size_t>(ArrayOfPtrs.size()), std::runtime_error,
745 ": ArrayOfPtrs.size() must be strictly positive and as large as ArrayOfPtrs.");
746
747 #ifdef HAVE_XPETRA_DEBUG
748 // This cannot be tested by Epetra itself
749 {
750 size_t localLength = map->getLocalNumElements();
751 for(int j=0; j<ArrayOfPtrs.size(); j++) {
752 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(Teuchos::as<size_t>(ArrayOfPtrs[j].size()) != localLength, std::runtime_error,
753 ": ArrayOfPtrs[" << j << "].size() (== " << ArrayOfPtrs[j].size() <<
754 ") is not equal to getLocalLength() (== " << localLength);
755
756 }
757 }
758 #endif
759
760 // Convert Teuchos::ArrayView< const Teuchos::ArrayView< const Scalar > > to double**
761 Array<const double*> arrayOfRawPtrs(ArrayOfPtrs.size());
762 for(int i=0; i<ArrayOfPtrs.size(); i++) {
763 arrayOfRawPtrs[i] = ArrayOfPtrs[i].getRawPtr();
764 }
765 double** rawArrayOfRawPtrs = const_cast<double**>(arrayOfRawPtrs.getRawPtr()); // This const_cast should be fine, because Epetra_DataAccess=Copy.
766
767 vec_ = Teuchos::rcp(new Epetra_MultiVector(Copy, toEpetra<GlobalOrdinal,Node>(map), rawArrayOfRawPtrs, NumVectors));
768 }
769
772
774
776
777
779 void replaceGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value) { XPETRA_MONITOR("EpetraMultiVectorT::replaceGlobalValue"); vec_->ReplaceGlobalValue(globalRow, Teuchos::as<int>(vectorIndex), value); }
780
782 void sumIntoGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value) { XPETRA_MONITOR("EpetraMultiVectorT::sumIntoGlobalValue"); vec_->SumIntoGlobalValue(globalRow, Teuchos::as<int>(vectorIndex), value); }
783
785 void replaceLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value) { XPETRA_MONITOR("EpetraMultiVectorT::replaceLocalValue"); vec_->ReplaceMyValue(myRow, Teuchos::as<int>(vectorIndex), value); }
786
788 void sumIntoLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value) { XPETRA_MONITOR("EpetraMultiVectorT::sumIntoLocalValue"); vec_->SumIntoMyValue(myRow, Teuchos::as<int>(vectorIndex), value); }
789
791 void putScalar(const Scalar &value) { XPETRA_MONITOR("EpetraMultiVectorT::putScalar"); vec_->PutScalar(value); }
792
794
796
797
800
803
806 XPETRA_MONITOR("EpetraMultiVectorT::getData");
807
808 double ** arrayOfPointers;
809
810 vec_->ExtractView(&arrayOfPointers);
811
812 double * data = arrayOfPointers[j];
813 int localLength = vec_->MyLength();
814
815 return ArrayRCP<double>(data, 0, localLength, false); // no ownership
816 }
817
820 XPETRA_MONITOR("EpetraMultiVectorT::getDataNonConst");
821
822 double ** arrayOfPointers;
823
824 vec_->ExtractView(&arrayOfPointers);
825
826 double * data = arrayOfPointers[j];
827 int localLength = vec_->MyLength();
828
829 return ArrayRCP<double>(data, 0, localLength, false); // no ownership
830 }
831
833
835
836
839 XPETRA_MONITOR("EpetraMultiVectorT::dot");
840
841 XPETRA_DYNAMIC_CAST(const EpetraMultiVectorT<GlobalOrdinal XPETRA_COMMA Node>, A, eA, "This Xpetra::EpetraMultiVectorT method only accept Xpetra::EpetraMultiVectorT as input arguments.");
842 vec_->Dot(*eA.getEpetra_MultiVector(), dots.getRawPtr());
843 }
844
846 void abs(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A) { XPETRA_MONITOR("EpetraMultiVectorT::abs"); vec_->Abs(toEpetra<GlobalOrdinal,Node>(A)); }
847
849 void reciprocal(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A) { XPETRA_MONITOR("EpetraMultiVectorT::reciprocal"); vec_->Reciprocal(toEpetra<GlobalOrdinal,Node>(A)); }
850
852 void scale(const Scalar &alpha) { XPETRA_MONITOR("EpetraMultiVectorT::scale"); vec_->Scale(alpha); }
853
856 XPETRA_MONITOR("EpetraMultiVectorT::scale");
857 // Epetra, unlike Tpetra, doesn't implement this version of
858 // scale(). Deal with this by scaling one column at a time.
859 const size_t numVecs = this->getNumVectors ();
860 for (size_t j = 0; j < numVecs; ++j) {
861 Epetra_Vector *v = (*vec_)(j);
862 v->Scale (alpha[j]);
863 }
864 }
865
867 void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta) { XPETRA_MONITOR("EpetraMultiVectorT::update"); vec_->Update(alpha, toEpetra<GlobalOrdinal,Node>(A), beta); }
868
870 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("EpetraMultiVectorT::update"); vec_->Update(alpha, toEpetra<GlobalOrdinal,Node>(A), beta, toEpetra<GlobalOrdinal,Node>(B), gamma); }
871
873 void norm1(const Teuchos::ArrayView< Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const { XPETRA_MONITOR("EpetraMultiVectorT::norm1"); vec_->Norm1(norms.getRawPtr()); }
874
876 void norm2(const Teuchos::ArrayView< Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const { XPETRA_MONITOR("EpetraMultiVectorT::norm2"); vec_->Norm2(norms.getRawPtr()); }
877
879 void normInf(const Teuchos::ArrayView< Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const { XPETRA_MONITOR("EpetraMultiVectorT::normInf"); vec_->NormInf(norms.getRawPtr()); }
880
882 void meanValue(const Teuchos::ArrayView< Scalar > &means) const { XPETRA_MONITOR("EpetraMultiVectorT::meanValue"); vec_->MeanValue(means.getRawPtr()); } //TODO: modify ArrayView size ??
883
885 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("EpetraMultiVectorT::multiply"); vec_->Multiply(toEpetra(transA), toEpetra(transB), alpha, toEpetra(A), toEpetra(B), beta); }
886
888 void elementWiseMultiply(Scalar scalarAB, const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, Scalar scalarThis) { XPETRA_MONITOR("EpetraMultiVectorT::elementWiseMultiply"); vec_->Multiply(scalarAB, toEpetra<GlobalOrdinal,Node>(A), toEpetra<GlobalOrdinal,Node>(B), scalarThis); }
889
891
893
894
896 size_t getNumVectors() const { XPETRA_MONITOR("EpetraMultiVectorT::getNumVectors"); return vec_->NumVectors(); }
897
899 size_t getLocalLength() const { XPETRA_MONITOR("EpetraMultiVectorT::getLocalLength"); return vec_->MyLength(); }
900
902 global_size_t getGlobalLength() const { XPETRA_MONITOR("EpetraMultiVectorT::getGlobalLength"); return vec_->GlobalLength64(); }
903
906 XPETRA_MONITOR("EpetraMultiVectorT::isSameSize");
907 auto vv = toEpetra<GlobalOrdinal,Node>(vec);
908 return ( (vec_->MyLength() == vv.MyLength()) && (vec_->NumVectors() == vv.NumVectors()));
909 }
910
912
914
915
917 std::string description() const {
918 XPETRA_MONITOR("EpetraMultiVectorT::description");
920 //return "TODO"; // unreachable
921 }
922
925 XPETRA_MONITOR("EpetraMultiVectorT::describe");
926 vec_->Print(out);
927 }
928
930
932 void randomize(bool bUseXpetraImplementation = false) {
933 XPETRA_MONITOR("EpetraMultiVectorT::randomize");
934
935 if (bUseXpetraImplementation)
937 else
938 vec_->Random();
939 }
940
942 void randomize(const Scalar& minVal, const Scalar& maxVal, bool bUseXpetraImplementation = false) {
943 XPETRA_MONITOR("EpetraMultiVectorT::randomize");
944
945 if (bUseXpetraImplementation)
947 else {
948 vec_->Random();
949 const size_t numVectors = getNumVectors();
950 for(size_t i = 0; i < numVectors; i++)
951 {
953
954 const size_t myLength = getLocalLength();
955 for(size_t j = 0; j < myLength; j++)
956 {
957 datai[ j ] = 0.5*(maxVal-minVal)*datai[ j ]+0.5*(maxVal+minVal);
958 }
959 }
960 }
961 }
962
964 //{@
965
967 Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const { XPETRA_MONITOR("EpetraMultiVectorT::getMap"); return toXpetra<GlobalOrdinal,Node>(vec_->Map()); }
968
971 XPETRA_MONITOR("EpetraMultiVectorT::doImport");
972
973 XPETRA_DYNAMIC_CAST(const EpetraMultiVectorT<GlobalOrdinal XPETRA_COMMA Node>, source, tSource, "Xpetra::EpetraMultiVectorT::doImport only accept Xpetra::EpetraMultiVectorT as input arguments.");
974 XPETRA_DYNAMIC_CAST(const EpetraImportT<GlobalOrdinal XPETRA_COMMA Node>, importer, tImporter, "Xpetra::EpetraMultiVectorT::doImport only accept Xpetra::EpetraImportT as input arguments.");
975
976 RCP<Epetra_MultiVector> v = tSource.getEpetra_MultiVector();
977 int err = this->getEpetra_MultiVector()->Import(*v, *tImporter.getEpetra_Import(), toEpetra(CM));
978 TEUCHOS_TEST_FOR_EXCEPTION(err != 0, std::runtime_error, "Catch error code returned by Epetra is " << err);
979 }
980
983 XPETRA_MONITOR("EpetraMultiVectorT::doExport");
984
985 XPETRA_DYNAMIC_CAST(const EpetraMultiVectorT<GlobalOrdinal XPETRA_COMMA Node>, dest, tDest, "Xpetra::EpetraMultiVectorT::doImport only accept Xpetra::EpetraMultiVectorT as input arguments.");
986 XPETRA_DYNAMIC_CAST(const EpetraImportT<GlobalOrdinal XPETRA_COMMA Node>, importer, tImporter, "Xpetra::EpetraMultiVectorT::doImport only accept Xpetra::EpetraImportT as input arguments.");
987
988 RCP<Epetra_MultiVector> v = tDest.getEpetra_MultiVector();
989 int err = this->getEpetra_MultiVector()->Export(*v, *tImporter.getEpetra_Import(), toEpetra(CM));
990 TEUCHOS_TEST_FOR_EXCEPTION(err != 0, std::runtime_error, "Catch error code returned by Epetra.");
991 }
992
995 XPETRA_MONITOR("EpetraMultiVectorT::doImport");
996
997 XPETRA_DYNAMIC_CAST(const EpetraMultiVectorT<GlobalOrdinal XPETRA_COMMA Node>, source, tSource, "Xpetra::EpetraMultiVectorT::doImport only accept Xpetra::EpetraMultiVectorT as input arguments.");
998 XPETRA_DYNAMIC_CAST(const EpetraExportT<GlobalOrdinal XPETRA_COMMA Node>, exporter, tExporter, "Xpetra::EpetraMultiVectorT::doImport only accept Xpetra::EpetraImportT as input arguments.");
999
1000 RCP<Epetra_MultiVector> v = tSource.getEpetra_MultiVector();
1001 int err = this->getEpetra_MultiVector()->Import(*v, *tExporter.getEpetra_Export(), toEpetra(CM));
1002 TEUCHOS_TEST_FOR_EXCEPTION(err != 0, std::runtime_error, "Catch error code returned by Epetra.");
1003 }
1004
1007 XPETRA_MONITOR("EpetraMultiVectorT::doExport");
1008
1009 XPETRA_DYNAMIC_CAST(const EpetraMultiVectorT<GlobalOrdinal XPETRA_COMMA Node>, dest, tDest, "Xpetra::EpetraMultiVectorT::doImport only accept Xpetra::EpetraMultiVectorT as input arguments.");
1010 XPETRA_DYNAMIC_CAST(const EpetraExportT<GlobalOrdinal XPETRA_COMMA Node>, exporter, tExporter, "Xpetra::EpetraMultiVectorT::doImport only accept Xpetra::EpetraImportT as input arguments.");
1011
1012 RCP<Epetra_MultiVector> v = tDest.getEpetra_MultiVector();
1013 int err = this->getEpetra_MultiVector()->Export(*v, *tExporter.getEpetra_Export(), toEpetra(CM));
1014 TEUCHOS_TEST_FOR_EXCEPTION(err != 0, std::runtime_error, "Catch error code returned by Epetra.");
1015 }
1016
1019 XPETRA_MONITOR("EpetraMultiVectorT::replaceMap");
1020 int err = 0;
1021 if (!map.is_null()) {
1022 err = this->getEpetra_MultiVector()->ReplaceMap(toEpetra<GlobalOrdinal,Node>(map));
1023
1024 } else {
1025 // Replace map with a dummy map to avoid potential hangs later
1026 Epetra_SerialComm SComm;
1027 Epetra_Map NewMap((GlobalOrdinal) vec_->MyLength(), (GlobalOrdinal) vec_->Map().IndexBase64(), SComm);
1028 err = this->getEpetra_MultiVector()->ReplaceMap(NewMap);
1029 }
1030 TEUCHOS_TEST_FOR_EXCEPTION(err != 0, std::runtime_error, "Catch error code returned by Epetra.");
1031 }
1032
1034
1036
1037
1039 EpetraMultiVectorT(const RCP<Epetra_MultiVector> &vec) : vec_(vec) { } //TODO removed const
1040
1043
1045 void setSeed(unsigned int seed) {
1046 XPETRA_MONITOR("EpetraMultiVectorT::seedrandom");
1047
1049 vec_->SetSeed(seed);
1050 }
1051
1053
1054 protected:
1057 virtual void
1059 typedef EpetraMultiVectorT this_type;
1060 const this_type* rhsPtr = dynamic_cast<const this_type*> (&rhs);
1062 rhsPtr == NULL, std::invalid_argument, "Xpetra::MultiVector::operator=: "
1063 "The left-hand side (LHS) of the assignment has a different type than "
1064 "the right-hand side (RHS). The LHS has type Xpetra::EpetraMultiVectorT "
1065 "(which means it wraps an Epetra_MultiVector), but the RHS has some "
1066 "other type. This probably means that the RHS wraps a Tpetra::Multi"
1067 "Vector. Xpetra::MultiVector does not currently implement assignment "
1068 "from a Tpetra object to an Epetra object, though this could be added "
1069 "with sufficient interest.");
1070
1071 RCP<const Epetra_MultiVector> rhsImpl = rhsPtr->getEpetra_MultiVector ();
1073
1075 rhsImpl.is_null (), std::logic_error, "Xpetra::MultiVector::operator= "
1076 "(in Xpetra::EpetraMultiVectorT::assign): *this (the right-hand side of "
1077 "the assignment) has a null RCP<Epetra_MultiVector> inside. Please "
1078 "report this bug to the Xpetra developers.");
1080 lhsImpl.is_null (), std::logic_error, "Xpetra::MultiVector::operator= "
1081 "(in Xpetra::EpetraMultiVectorT::assign): The left-hand side of the "
1082 "assignment has a null RCP<Epetra_MultiVector> inside. Please report "
1083 "this bug to the Xpetra developers.");
1084
1085 // Epetra_MultiVector's assignment operator does a deep copy.
1086 *lhsImpl = *rhsImpl;
1087 }
1088
1089 private:
1092
1093 }; // EpetraMultiVectorT class (specialization on GO=long long, NO=EpetraNode
1094#endif
1095
1096} // Xpetra namespace
1097
1098#include "Xpetra_EpetraVector.hpp" // to avoid incomplete type instantiated above in out-of-body functions.
1099
1100#endif // XPETRA_EPETRAMULTIVECTOR_HPP
#define XPETRA_MONITOR(funcName)
#define XPETRA_DYNAMIC_CAST(type, obj, newObj, exceptionMsg)
int Scale(double ScalarValue)
T * getRawPtr() const
static const EVerbosityLevel verbLevel_default
bool is_null() const
void doImport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Import (using an Exporter).
std::string description() const
A simple one-line description of this object.
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 replaceGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value)
Replace value, using global (row) index.
Teuchos::ArrayRCP< Scalar > getDataNonConst(size_t j)
View of the local values in a particular vector of this multivector.
dual_view_type::t_host_um getHostLocalView(Access::ReadWriteStruct) const override
void sumIntoGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value)
Add value to existing value, using global (row) index.
size_t getLocalLength() const
Local number of rows on the calling process.
void norm2(const Teuchos::ArrayView< Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
dual_view_type::t_dev_um getDeviceLocalView(Access::OverwriteAllStruct) const override
void replaceMap(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map)
Replace the underlying Map in place.
RCP< Epetra_MultiVector > getEpetra_MultiVector() const
Get the underlying Epetra multivector.
void dot(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Teuchos::ArrayView< Scalar > &dots) const
Compute the dot product of each corresponding pair of vectors (columns) in A and B.
void randomize(bool bUseXpetraImplementation=false)
Set multi-vector values to random numbers.
void abs(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
Put element-wise absolute values of input Multi-vector in target: A = abs(this).
void sumIntoLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value)
Add value to existing value, using local (row) index.
bool isSameSize(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec) const
Checks to see if the local length, number of vectors and size of Scalar type match.
void setSeed(unsigned int seed)
Set seed for Random function.
Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node >::dual_view_type dual_view_type
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: this = gamma*this + alpha*A + beta*B.
void replaceLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value)
Replace value, using local (row) index.
EpetraMultiVectorT(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Teuchos::DataAccess copyOrView=Teuchos::Copy)
MultiVector copy constructor.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel=Teuchos::Describable::verbLevel_default) const
Print the object with the given verbosity level to a FancyOStream.
void randomize(const Scalar &minVal, const Scalar &maxVal, bool bUseXpetraImplementation=false)
Set multi-vector values to random numbers.
void putScalar(const Scalar &value)
Set all values in the multivector with the given value.
dual_view_type::t_host_const_um getHostLocalView(Access::ReadOnlyStruct) const override
void scale(const Scalar &alpha)
Scale in place: 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,...
void normInf(const Teuchos::ArrayView< Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
Compute Inf-norm of each vector in multi-vector.
void elementWiseMultiply(Scalar scalarAB, const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, Scalar scalarThis)
Multiply a Vector A elementwise by a MultiVector B.
size_t getNumVectors() const
Number of columns in the multivector.
void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta)
Update: this = beta*this + alpha*A.
void norm1(const Teuchos::ArrayView< Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
Compute 1-norm of each vector in multi-vector.
void doImport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import.
RCP< Epetra_MultiVector > vec_
The Epetra_MultiVector which this class wraps.
dual_view_type::t_host_um getHostLocalView(Access::OverwriteAllStruct) const override
EpetraMultiVectorT(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, const Teuchos::ArrayView< const Teuchos::ArrayView< const Scalar > > &ArrayOfPtrs, size_t NumVectors)
Set multi-vector values from array of pointers using Teuchos memory management classes....
dual_view_type::t_dev_um getDeviceLocalView(Access::ReadWriteStruct) const override
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
dual_view_type::t_dev_const_um getDeviceLocalView(Access::ReadOnlyStruct) const override
EpetraMultiVectorT(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, size_t NumVectors, bool zeroOut=true)
Basic MultiVector constuctor.
void doExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Export (using an Importer).
void doExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
Teuchos::ArrayRCP< const Scalar > getData(size_t j) const
Const view of the local values in a particular vector of this multivector.
EpetraMultiVectorT(const RCP< Epetra_MultiVector > &vec)
EpetraMultiVectorT constructor to wrap a Epetra_MultiVector object.
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].
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...
global_size_t getGlobalLength() const
Global number of rows in the multivector.
void norm1(const Teuchos::ArrayView< Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
Compute 1-norm of each vector in multi-vector.
EpetraMultiVectorT(const RCP< Epetra_MultiVector > &vec)
EpetraMultiVectorT constructor to wrap a Epetra_MultiVector object.
Teuchos::RCP< const Vector< double, int, long long, EpetraNode > > getVector(size_t j) const
Return a Vector which is a const view of column j.
void normInf(const Teuchos::ArrayView< Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
Compute Inf-norm of each vector in multi-vector.
EpetraMultiVectorT(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, size_t NumVectors, bool zeroOut=true)
Basic MultiVector constuctor.
size_t getNumVectors() const
Number of columns in the multivector.
Teuchos::ArrayRCP< Scalar > getDataNonConst(size_t j)
View of the local values in a particular vector of this multivector.
void elementWiseMultiply(Scalar scalarAB, const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, Scalar scalarThis)
Multiply a Vector A elementwise by a MultiVector B.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel=Teuchos::Describable::verbLevel_default) const
Print the object with the given verbosity level to a FancyOStream.
void randomize(const Scalar &minVal, const Scalar &maxVal, bool bUseXpetraImplementation=false)
Set multi-vector values to random numbers.
EpetraMultiVectorT(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, const Teuchos::ArrayView< const Teuchos::ArrayView< const Scalar > > &ArrayOfPtrs, size_t NumVectors)
Set multi-vector values from array of pointers using Teuchos memory management classes....
void putScalar(const Scalar &value)
Set all values in the multivector with the given value.
void replaceGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value)
Replace value, using global (row) index.
void scale(const Scalar &alpha)
Scale in place: this = alpha*this.
Teuchos::RCP< Vector< double, int, long long, EpetraNode > > getVectorNonConst(size_t j)
Return a Vector which is a nonconst view of column j.
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,...
void scale(Teuchos::ArrayView< const Scalar > alpha)
Scale the current values of a multi-vector, this[j] = alpha[j]*this[j].
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: this = gamma*this + alpha*A + beta*B.
virtual void assign(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &rhs)
Implementation of the assignment operator (operator=); does a deep copy.
void sumIntoGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value)
Add value to existing value, using global (row) index.
RCP< Epetra_MultiVector > vec_
The Epetra_MultiVector which this class wraps.
void doExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
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).
EpetraMultiVectorT(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source)
MultiVector copy constructor.
size_t getLocalLength() const
Local number of rows on the calling process.
void replaceLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value)
Replace value, using local (row) index.
global_size_t getGlobalLength() const
Global number of rows in the multivector.
bool isSameSize(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec) const
Checks to see if the local length, number of vectors and size of Scalar type match.
void dot(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Teuchos::ArrayView< Scalar > &dots) const
Compute the dot product of each corresponding pair of vectors (columns) in A and B.
void doImport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import.
std::string description() const
A simple one-line description of this object.
void sumIntoLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value)
Add value to existing value, using local (row) index.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
void setSeed(unsigned int seed)
Set seed for Random function.
void randomize(bool bUseXpetraImplementation=false)
Set multi-vector values to random numbers.
void norm2(const Teuchos::ArrayView< Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
RCP< Epetra_MultiVector > getEpetra_MultiVector() const
Get the underlying Epetra multivector.
void doExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Export (using an Importer).
void abs(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
Put element-wise absolute values of input Multi-vector in target: A = abs(this).
Teuchos::ArrayRCP< const Scalar > getData(size_t j) const
Const view of the local values in a particular vector of this multivector.
void replaceMap(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map)
Replace the underlying Map in place.
void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta)
Update: this = beta*this + alpha*A.
void doImport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Import (using an Exporter).
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...
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 doImport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Import (using an Exporter).
void norm1(const Teuchos::ArrayView< Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
Compute 1-norm of each vector in multi-vector.
void doExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Export (using an Importer).
void dot(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Teuchos::ArrayView< Scalar > &dots) const
Compute the dot product of each corresponding pair of vectors (columns) in A and B.
void doExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
void putScalar(const Scalar &value)
Set all values in the multivector with the given value.
virtual void assign(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &rhs)
Implementation of the assignment operator (operator=); does a deep copy.
void norm2(const Teuchos::ArrayView< Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) 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.
void setSeed(unsigned int seed)
Set seed for Random function.
Teuchos::RCP< Vector< double, int, GlobalOrdinal, Node > > getVectorNonConst(size_t j)
Return a Vector which is a nonconst view of column j.
EpetraMultiVectorT(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Teuchos::DataAccess copyOrView=Teuchos::Copy)
MultiVector copy constructor.
Teuchos::RCP< const Vector< double, int, GlobalOrdinal, Node > > getVector(size_t j) const
Return a Vector which is a const view of column j.
void normInf(const Teuchos::ArrayView< Teuchos::ScalarTraits< Scalar >::magnitudeType > &norms) const
Compute Inf-norm of each vector in multi-vector.
void sumIntoGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value)
Add value to existing 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: 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 replaceMap(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map)
Replace the underlying Map in place.
std::string description() const
A simple one-line description of this object.
void elementWiseMultiply(Scalar scalarAB, const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, Scalar scalarThis)
Multiply a Vector A elementwise by a MultiVector B.
bool isSameSize(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec) const
void randomize(const Scalar &minVal, const Scalar &maxVal, bool bUseXpetraImplementation=false)
Set multi-vector values to random numbers.
void scale(Teuchos::ArrayView< const Scalar > alpha)
Scale the current values of a multi-vector, this[j] = alpha[j]*this[j].
EpetraMultiVectorT(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, size_t NumVectors, bool zeroOut=true)
Basic MultiVector constuctor.
Teuchos::ArrayRCP< const Scalar > getData(size_t j) const
Const view of the local values in a particular vector of this multivector.
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...
EpetraMultiVectorT(const RCP< Epetra_MultiVector > &vec)
EpetraMultiVectorT constructor to wrap a Epetra_MultiVector object.
void scale(const Scalar &alpha)
Scale in place: 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,...
EpetraMultiVectorT(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, const Teuchos::ArrayView< const Teuchos::ArrayView< const Scalar > > &ArrayOfPtrs, size_t NumVectors)
Set multi-vector values from array of pointers using Teuchos memory management classes....
void doImport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import.
size_t getLocalLength() const
Local number of rows on the calling process.
RCP< Epetra_MultiVector > getEpetra_MultiVector() const
Get the underlying Epetra multivector.
void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta)
Update: this = beta*this + alpha*A.
size_t getNumVectors() const
Number of columns in the multivector.
void replaceLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value)
Replace value, using local (row) index.
virtual ~EpetraMultiVectorT()
MultiVector destructor.
void randomize(bool bUseXpetraImplementation=false)
Set multi-vector values to random numbers.
global_size_t getGlobalLength() const
Global number of rows in the multivector.
void replaceGlobalValue(GlobalOrdinal globalRow, size_t vectorIndex, const Scalar &value)
Replace value, using global (row) index.
Teuchos::ArrayRCP< Scalar > getDataNonConst(size_t j)
View of the local values in a particular vector of this multivector.
void sumIntoLocalValue(LocalOrdinal myRow, size_t vectorIndex, const Scalar &value)
Add value to existing value, using local (row) index.
Exception throws when you call an unimplemented method of Xpetra.
Exception throws to report errors in the internal logical of the program.
virtual size_t getNumVectors() const =0
Number of columns in the multivector.
virtual void Xpetra_randomize()
Set multi-vector values to random numbers. XPetra implementation.
Kokkos::DualView< impl_scalar_type **, Kokkos::LayoutStride, typename node_type::device_type, Kokkos::MemoryUnmanaged > dual_view_type
virtual dual_view_type::t_host_const_um getHostLocalView(Access::ReadOnlyStruct) const
#define TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(throw_exception_test, Exception, msg)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
TypeTo as(const TypeFrom &t)
#define TEUCHOS_UNREACHABLE_RETURN(dummyReturnVal)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
View
constexpr struct ReadWriteStruct ReadWrite
Xpetra namespace
size_t global_size_t
Global size_t object.
const Epetra_CrsGraph & toEpetra(const RCP< const CrsGraph< int, GlobalOrdinal, Node > > &graph)
RCP< const CrsGraph< int, GlobalOrdinal, Node > > toXpetra(const Epetra_CrsGraph &g)
CombineMode
Xpetra::Combine Mode enumerable type.
static void seedrandom(unsigned int s)