Xpetra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Xpetra_TpetraCrsMatrix_decl.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Xpetra: A linear algebra interface package
4 //
5 // Copyright 2012 NTESS and the Xpetra contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 #ifndef XPETRA_TPETRACRSMATRIX_DECL_HPP
11 #define XPETRA_TPETRACRSMATRIX_DECL_HPP
12 
13 /* this file is automatically generated - do not edit (see scripts/tpetra.py) */
14 
15 // FIXME (mfh 03 Sep 2014) The above is probably not true anymore.
16 // Furthermore, I don't think anyone maintains the scripts.
17 // Feel free to correct this comment if I'm wrong.
18 
20 
21 #include "Tpetra_CrsMatrix.hpp"
22 #include "Tpetra_replaceDiagonalCrsMatrix.hpp"
23 
24 #include "Xpetra_CrsMatrix.hpp"
30 #include "Xpetra_TpetraCrsGraph.hpp"
31 #include "Xpetra_Exceptions.hpp"
32 
33 namespace Xpetra {
34 
35 template <class Scalar,
36  class LocalOrdinal,
37  class GlobalOrdinal,
38  class Node = Tpetra::KokkosClassic::DefaultNode::DefaultNodeType>
40  : public CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> //, public TpetraRowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>
41 {
42 #undef XPETRA_TPETRACRSMATRIX_SHORT
44 
45  // The following typedef are used by the XPETRA_DYNAMIC_CAST() macro.
50 
51  // The following typedefs are used by the Kokkos interface
52 #ifdef HAVE_XPETRA_TPETRA
54 #endif
55 
56  public:
58 
59 
61  TpetraCrsMatrix(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap, size_t maxNumEntriesPerRow, const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null);
62 
64  TpetraCrsMatrix(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap, const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc, const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null);
65 
67  TpetraCrsMatrix(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap, const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &colMap, size_t maxNumEntriesPerRow, const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null);
68 
70  TpetraCrsMatrix(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap, const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &colMap, const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc, const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null);
71 
73  TpetraCrsMatrix(const Teuchos::RCP<const CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > &graph, const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null);
74 
76  TpetraCrsMatrix(const Teuchos::RCP<const CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > &graph, typename Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::local_matrix_type::values_type &values, const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null);
77 
79  TpetraCrsMatrix(const Teuchos::RCP<const CrsMatrix> &sourceMatrix,
81  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap = Teuchos::null,
82  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap = Teuchos::null,
83  const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null);
84 
86  TpetraCrsMatrix(const Teuchos::RCP<const CrsMatrix> &sourceMatrix,
88  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap = Teuchos::null,
89  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap = Teuchos::null,
90  const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null);
91 
93  TpetraCrsMatrix(const Teuchos::RCP<const CrsMatrix> &sourceMatrix,
95  const Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Node> > DomainImporter,
96  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap,
97  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap,
98  const Teuchos::RCP<Teuchos::ParameterList> &params);
99 
101  TpetraCrsMatrix(const Teuchos::RCP<const CrsMatrix> &sourceMatrix,
103  const Teuchos::RCP<const Export<LocalOrdinal, GlobalOrdinal, Node> > DomainExporter,
104  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap,
105  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap,
106  const Teuchos::RCP<Teuchos::ParameterList> &params);
107 
108 #ifdef HAVE_XPETRA_TPETRA
109  TpetraCrsMatrix(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap,
130  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &colMap,
131  const local_matrix_type &lclMatrix,
132  const Teuchos::RCP<Teuchos::ParameterList> &params = null);
133 
136  const local_matrix_type &lclMatrix,
137  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap,
138  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &colMap,
139  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap = Teuchos::null,
140  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap = Teuchos::null,
141  const Teuchos::RCP<Teuchos::ParameterList> &params = null);
142 
145  const local_matrix_type &lclMatrix,
146  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap,
147  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &colMap,
148  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap,
149  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap,
150  const Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Node> > &importer,
151  const Teuchos::RCP<const Export<LocalOrdinal, GlobalOrdinal, Node> > &exporter,
152  const Teuchos::RCP<Teuchos::ParameterList> &params = null);
153 #endif
154 
156  virtual ~TpetraCrsMatrix();
157 
159 
161 
162 
164  void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView<const GlobalOrdinal> &cols, const ArrayView<const Scalar> &vals);
165 
167  void insertLocalValues(LocalOrdinal localRow, const ArrayView<const LocalOrdinal> &cols, const ArrayView<const Scalar> &vals);
168 
170  void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView<const GlobalOrdinal> &cols, const ArrayView<const Scalar> &vals);
171 
173  void
174  replaceLocalValues(LocalOrdinal localRow,
175  const ArrayView<const LocalOrdinal> &cols,
176  const ArrayView<const Scalar> &vals);
177 
179  void setAllToScalar(const Scalar &alpha);
180 
182  void scale(const Scalar &alpha);
183 
185  //** \warning This is an expert-only routine and should not be called from user code. */
186  void allocateAllValues(size_t numNonZeros, ArrayRCP<size_t> &rowptr, ArrayRCP<LocalOrdinal> &colind, ArrayRCP<Scalar> &values);
187 
189  void setAllValues(const ArrayRCP<size_t> &rowptr, const ArrayRCP<LocalOrdinal> &colind, const ArrayRCP<Scalar> &values);
190 
192  void getAllValues(ArrayRCP<const size_t> &rowptr, ArrayRCP<const LocalOrdinal> &colind, ArrayRCP<const Scalar> &values) const;
193 
195  void getAllValues(ArrayRCP<Scalar> &values);
196 
197  bool haveGlobalConstants() const;
198 
200 
202 
203 
205  void resumeFill(const RCP<ParameterList> &params = null);
206 
208  void fillComplete(const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap, const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap, const RCP<ParameterList> &params = null);
209 
211  void fillComplete(const RCP<ParameterList> &params = null);
212 
214  void replaceDomainMapAndImporter(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &newDomainMap, Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Node> > &newImporter);
215 
217  void expertStaticFillComplete(const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap,
218  const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap,
219  const RCP<const Import<LocalOrdinal, GlobalOrdinal, Node> > &importer = Teuchos::null,
220  const RCP<const Export<LocalOrdinal, GlobalOrdinal, Node> > &exporter = Teuchos::null,
221  const RCP<ParameterList> &params = Teuchos::null);
222 
224 
226 
227 
229  const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > getRowMap() const;
230 
232  const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > getColMap() const;
233 
235  RCP<const CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > getCrsGraph() const;
236 
239 
242 
244  size_t getLocalNumRows() const;
245 
247  size_t getLocalNumCols() const;
248 
251 
253  size_t getLocalNumEntries() const;
254 
256  size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const;
257 
259  size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const;
260 
262  size_t getGlobalMaxNumRowEntries() const;
263 
265  size_t getLocalMaxNumRowEntries() const;
266 
268  bool isLocallyIndexed() const;
269 
271  bool isGloballyIndexed() const;
272 
274  bool isFillComplete() const;
275 
277  bool isFillActive() const;
278 
280  typename ScalarTraits<Scalar>::magnitudeType getFrobeniusNorm() const;
281 
283  bool supportsRowViews() const;
284 
286  void getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView<LocalOrdinal> &Indices, const ArrayView<Scalar> &Values, size_t &NumEntries) const;
287 
289  void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView<const GlobalOrdinal> &indices, ArrayView<const Scalar> &values) const;
290 
292  void getGlobalRowCopy(GlobalOrdinal GlobalRow, const ArrayView<GlobalOrdinal> &indices, const ArrayView<Scalar> &values, size_t &numEntries) const;
293 
295  void getLocalRowView(LocalOrdinal LocalRow, ArrayView<const LocalOrdinal> &indices, ArrayView<const Scalar> &values) const;
296 
298 
300 
301 
303  void apply(const MultiVector &X, MultiVector &Y, Teuchos::ETransp mode = Teuchos::NO_TRANS, Scalar alpha = ScalarTraits<Scalar>::one(), Scalar beta = ScalarTraits<Scalar>::zero()) const;
304 
306  void apply(const MultiVector &X, MultiVector &Y, Teuchos::ETransp mode, Scalar alpha, Scalar beta, bool sumInterfaceValues, const RCP<Xpetra::Import<LocalOrdinal, GlobalOrdinal, Node> > &regionInterfaceImporter, const Teuchos::ArrayRCP<LocalOrdinal> &regionInterfaceLIDs) const;
307 
309  const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > getDomainMap() const;
310 
312  const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > getRangeMap() const;
313 
315 
317 
318 
320  std::string description() const;
321 
323  void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel = Teuchos::Describable::verbLevel_default) const;
324 
326 
328 
329  void setObjectLabel(const std::string &objectLabel);
331 
333  TpetraCrsMatrix(const TpetraCrsMatrix &matrix);
334 
336  void getLocalDiagCopy(Vector &diag) const;
337 
339  void getLocalDiagOffsets(Teuchos::ArrayRCP<size_t> &offsets) const;
340 
342  void getLocalDiagCopy(Vector &diag, const Teuchos::ArrayView<const size_t> &offsets) const;
343 
345  void getLocalDiagCopy(Vector &diag, const Kokkos::View<const size_t *, typename Node::device_type, Kokkos::MemoryUnmanaged> &offsets) const;
346 
348  void replaceDiag(const Vector &diag);
349 
351  void leftScale(const Vector &x);
352 
354  void rightScale(const Vector &x);
355 
357  //{@
358 
360  Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > getMap() const;
361 
365 
369 
373 
377 
378  void removeEmptyProcessesInPlace(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &newMap);
379 
380  // @}
382 
383 
385  bool hasMatrix() const;
386 
388  TpetraCrsMatrix(const Teuchos::RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &mtx);
389 
391  RCP<const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > getTpetra_CrsMatrix() const;
392 
394  RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > getTpetra_CrsMatrixNonConst() const; // TODO: remove
395 
396 #ifdef HAVE_XPETRA_TPETRA
397  typename local_matrix_type::HostMirror getLocalMatrixHost() const {
399  return getTpetra_CrsMatrixNonConst()->getLocalMatrixHost();
400  }
403  return getTpetra_CrsMatrixNonConst()->getLocalMatrixDevice();
404  }
405 
406  void setAllValues(const typename local_matrix_type::row_map_type &ptr,
407  const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type &ind,
408  const typename local_matrix_type::values_type &val) {
409  getTpetra_CrsMatrixNonConst()->setAllValues(ptr, ind, val);
410  }
411 #endif
412 
414  LocalOrdinal GetStorageBlockSize() const { return 1; }
415 
417  void residual(const MultiVector &X,
418  const MultiVector &B,
419  MultiVector &R) const;
420 
422 
423  private:
424  RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > mtx_;
425 }; // TpetraCrsMatrix class
426 
427 #ifdef HAVE_XPETRA_EPETRA
428 
429 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
430  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
431 
432 // specialization of TpetraCrsMatrix for GO=LO=int
433 template <class Scalar>
434 class TpetraCrsMatrix<Scalar, int, int, EpetraNode>
435  : public CrsMatrix<Scalar, int, int, EpetraNode> //, public TpetraRowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>
436 {
437  typedef int LocalOrdinal;
438  typedef int GlobalOrdinal;
439  typedef EpetraNode Node;
440  // The following typedef are used by the XPETRA_DYNAMIC_CAST() macro.
445 
446  // The following typedefs are used by the Kokkos interface
447 #ifdef HAVE_XPETRA_TPETRA
449 #endif
450 
451  public:
453 
454 
456  TpetraCrsMatrix(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap, size_t maxNumEntriesPerRow, const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null) {
458  }
459 
461  TpetraCrsMatrix(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap, const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc, const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null) {
463  }
464 
466  TpetraCrsMatrix(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap, const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &colMap, size_t maxNumEntriesPerRow, const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null) {
468  }
469 
471  TpetraCrsMatrix(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap, const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &colMap, const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc, const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null) {
473  }
474 
476  TpetraCrsMatrix(const Teuchos::RCP<const CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > &graph, const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null) {
478  }
479 
483  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap = Teuchos::null,
484  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap = Teuchos::null,
485  const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null) {
487  }
488 
492  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap = Teuchos::null,
493  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap = Teuchos::null,
494  const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null) {
496  }
497 
501  const Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Node> > DomainImporter,
502  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap,
503  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap,
504  const Teuchos::RCP<Teuchos::ParameterList> &params) {
506  }
507 
511  const Teuchos::RCP<const Export<LocalOrdinal, GlobalOrdinal, Node> > DomainExporter,
512  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap,
513  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap,
514  const Teuchos::RCP<Teuchos::ParameterList> &params) {
516  }
517 
518 #ifdef HAVE_XPETRA_TPETRA
519  TpetraCrsMatrix(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap,
540  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &colMap,
541  const local_matrix_type &lclMatrix,
542  const Teuchos::RCP<Teuchos::ParameterList> &params = null) {
543  XPETRA_TPETRA_ETI_EXCEPTION(typeid(TpetraCrsMatrix<Scalar, int, int, EpetraNode>).name(), "TpetraCrsMatrix<int,int>", "int", typeid(EpetraNode).name());
544  }
545 
548  const local_matrix_type &lclMatrix,
549  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap,
550  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &colMap,
551  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap = Teuchos::null,
552  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap = Teuchos::null,
553  const Teuchos::RCP<Teuchos::ParameterList> &params = null) {
554  XPETRA_TPETRA_ETI_EXCEPTION(typeid(TpetraCrsMatrix<Scalar, int, int, EpetraNode>).name(), "TpetraCrsMatrix<int,int>", "int", typeid(EpetraNode).name());
555  }
556 #endif
557 
559  virtual ~TpetraCrsMatrix() {}
560 
562 
564 
565 
567  void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView<const GlobalOrdinal> &cols, const ArrayView<const Scalar> &vals) {}
568 
570  void insertLocalValues(LocalOrdinal localRow, const ArrayView<const LocalOrdinal> &cols, const ArrayView<const Scalar> &vals) {}
571 
573  void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView<const GlobalOrdinal> &cols, const ArrayView<const Scalar> &vals) {}
574 
576  void
578  const ArrayView<const LocalOrdinal> &cols,
579  const ArrayView<const Scalar> &vals) {}
580 
582  void setAllToScalar(const Scalar &alpha) {}
583 
585  void scale(const Scalar &alpha) {}
586 
588  //** \warning This is an expert-only routine and should not be called from user code. */
589  void allocateAllValues(size_t numNonZeros, ArrayRCP<size_t> &rowptr, ArrayRCP<LocalOrdinal> &colind, ArrayRCP<Scalar> &values) {}
590 
592  void setAllValues(const ArrayRCP<size_t> &rowptr, const ArrayRCP<LocalOrdinal> &colind, const ArrayRCP<Scalar> &values) {}
593 
595  void getAllValues(ArrayRCP<const size_t> &rowptr, ArrayRCP<const LocalOrdinal> &colind, ArrayRCP<const Scalar> &values) const {}
596 
598  void getAllValues(ArrayRCP<Scalar> &values) {}
599 
600  bool haveGlobalConstants() const { return false; }
601 
603 
605 
606 
608  void resumeFill(const RCP<ParameterList> &params = null) {}
609 
611  void fillComplete(const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap, const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap, const RCP<ParameterList> &params = null) {}
612 
614  void fillComplete(const RCP<ParameterList> &params = null) {}
615 
617  void replaceDomainMapAndImporter(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &newDomainMap, Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Node> > &newImporter) {}
618 
621  const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap,
622  const RCP<const Import<LocalOrdinal, GlobalOrdinal, Node> > &importer = Teuchos::null,
623  const RCP<const Export<LocalOrdinal, GlobalOrdinal, Node> > &exporter = Teuchos::null,
624  const RCP<ParameterList> &params = Teuchos::null) {}
625 
627 
629 
630 
632  const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > getRowMap() const { return Teuchos::null; }
633 
635  const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > getColMap() const { return Teuchos::null; }
636 
638  RCP<const CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > getCrsGraph() const { return Teuchos::null; }
639 
641  global_size_t getGlobalNumRows() const { return 0; }
642 
644  global_size_t getGlobalNumCols() const { return 0; }
645 
647  size_t getLocalNumRows() const { return 0; }
648 
650  size_t getLocalNumCols() const { return 0; }
651 
653  global_size_t getGlobalNumEntries() const { return 0; }
654 
656  size_t getLocalNumEntries() const { return 0; }
657 
659  size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const { return 0; }
660 
662  size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const { return 0; }
663 
665  size_t getGlobalMaxNumRowEntries() const { return 0; }
666 
668  size_t getLocalMaxNumRowEntries() const { return 0; }
669 
671  bool isLocallyIndexed() const { return false; }
672 
674  bool isGloballyIndexed() const { return false; }
675 
677  bool isFillComplete() const { return false; }
678 
680  bool isFillActive() const { return false; }
681 
683  typename ScalarTraits<Scalar>::magnitudeType getFrobeniusNorm() const { return ScalarTraits<Scalar>::magnitude(ScalarTraits<Scalar>::zero()); }
684 
686  bool supportsRowViews() const { return false; }
687 
689  void getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView<LocalOrdinal> &Indices, const ArrayView<Scalar> &Values, size_t &NumEntries) const {}
690 
692  void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView<const GlobalOrdinal> &indices, ArrayView<const Scalar> &values) const {}
693 
695  void getGlobalRowCopy(GlobalOrdinal GlobalRow, const ArrayView<GlobalOrdinal> &indices, const ArrayView<Scalar> &values, size_t &numEntries) const {}
696 
698  void getLocalRowView(LocalOrdinal LocalRow, ArrayView<const LocalOrdinal> &indices, ArrayView<const Scalar> &values) const {}
699 
701 
703 
704 
706  void apply(const MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> &X, MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> &Y, Teuchos::ETransp mode = Teuchos::NO_TRANS, Scalar alpha = ScalarTraits<Scalar>::one(), Scalar beta = ScalarTraits<Scalar>::zero()) const {}
707 
709  void apply(const MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> &X, MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> &Y, Teuchos::ETransp mode, Scalar alpha, Scalar beta, bool sumInterfaceValues, const RCP<Xpetra::Import<LocalOrdinal, GlobalOrdinal, Node> > &regionInterfaceImporter, const Teuchos::ArrayRCP<LocalOrdinal> &regionInterfaceLIDs) const {}
710 
712  const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > getDomainMap() const { return Teuchos::null; }
713 
715  const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > getRangeMap() const { return Teuchos::null; }
716 
718 
720 
721 
723  std::string description() const { return std::string(""); }
724 
726  void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel = Teuchos::Describable::verbLevel_default) const {}
727 
729 
731 
732  void setObjectLabel(const std::string &objectLabel) {}
734 
737 
740 
742  void getLocalDiagOffsets(Teuchos::ArrayRCP<size_t> &offsets) const {}
743 
745  void getLocalDiagCopy(Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node> &diag, const Teuchos::ArrayView<const size_t> &offsets) const {}
746 
749 
752 
755 
757  //{@
758 
760  Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > getMap() const { return Teuchos::null; }
761 
765 
769 
773 
777 
778  void removeEmptyProcessesInPlace(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &newMap) {}
779 
780  // @}
782 
783 
785  bool hasMatrix() const { return false; }
786 
788  TpetraCrsMatrix(const Teuchos::RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &mtx) {
790  }
791 
793  RCP<const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > getTpetra_CrsMatrix() const { return Teuchos::null; }
794 
796  RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > getTpetra_CrsMatrixNonConst() const { return Teuchos::null; } // TODO: remove
797 
798 #ifdef HAVE_XPETRA_TPETRA
799  local_matrix_type getLocalMatrix() const {
801  TEUCHOS_UNREACHABLE_RETURN(local_matrix_type());
802  }
803 
804  void setAllValues(const typename local_matrix_type::row_map_type &ptr,
805  const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type &ind,
806  const typename local_matrix_type::values_type &val) {}
807 #endif
808 
810  LocalOrdinal GetStorageBlockSize() const { return 1; }
811 
816 
818 }; // TpetraCrsMatrix class (specialization for GO=int, NO=EpetraNode)
819 #endif
820 
821 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_LONG_LONG))) || \
822  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_LONG_LONG))))
823 
824 // specialization of TpetraCrsMatrix for GO=long long, NO=EpetraNode
825 template <class Scalar>
826 class TpetraCrsMatrix<Scalar, int, long long, EpetraNode>
827  : public CrsMatrix<Scalar, int, long long, EpetraNode> //, public TpetraRowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>
828 {
829  typedef int LocalOrdinal;
830  typedef long long GlobalOrdinal;
831  typedef EpetraNode Node;
832  // The following typedef are used by the XPETRA_DYNAMIC_CAST() macro.
837 
838  // The following typedefs are used by the Kokkos interface
839 #ifdef HAVE_XPETRA_TPETRA
841 #endif
842 
843  public:
845 
846 
848  TpetraCrsMatrix(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap, size_t maxNumEntriesPerRow, const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null) {
850  }
851 
853  TpetraCrsMatrix(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap, const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc, const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null) {
855  }
856 
858  TpetraCrsMatrix(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap, const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &colMap, size_t maxNumEntriesPerRow, const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null) {
860  }
861 
863  TpetraCrsMatrix(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap, const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &colMap, const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc, const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null) {
865  }
866 
868  TpetraCrsMatrix(const Teuchos::RCP<const CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > &graph, const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null) {
870  }
871 
875  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap = Teuchos::null,
876  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap = Teuchos::null,
877  const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null) {
879  }
880 
884  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap = Teuchos::null,
885  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap = Teuchos::null,
886  const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null) {
888  }
889 
893  const Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Node> > DomainImporter,
894  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap,
895  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap,
896  const Teuchos::RCP<Teuchos::ParameterList> &params) {
898  }
899 
903  const Teuchos::RCP<const Export<LocalOrdinal, GlobalOrdinal, Node> > DomainExporter,
904  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap,
905  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap,
906  const Teuchos::RCP<Teuchos::ParameterList> &params) {
908  }
909 
910 #ifdef HAVE_XPETRA_TPETRA
911  TpetraCrsMatrix(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap,
932  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &colMap,
933  const local_matrix_type &lclMatrix,
934  const Teuchos::RCP<Teuchos::ParameterList> &params = null) {
936  }
937 
940  const local_matrix_type &lclMatrix,
941  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rowMap,
942  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &colMap,
943  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap = Teuchos::null,
944  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap = Teuchos::null,
945  const Teuchos::RCP<Teuchos::ParameterList> &params = null) {
947  }
948 #endif
949 
951  virtual ~TpetraCrsMatrix() {}
952 
954 
956 
957 
959  void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView<const GlobalOrdinal> &cols, const ArrayView<const Scalar> &vals) {}
960 
962  void insertLocalValues(LocalOrdinal localRow, const ArrayView<const LocalOrdinal> &cols, const ArrayView<const Scalar> &vals) {}
963 
965  void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView<const GlobalOrdinal> &cols, const ArrayView<const Scalar> &vals) {}
966 
968  void
970  const ArrayView<const LocalOrdinal> &cols,
971  const ArrayView<const Scalar> &vals) {}
972 
974  void setAllToScalar(const Scalar &alpha) {}
975 
977  void scale(const Scalar &alpha) {}
978 
980  //** \warning This is an expert-only routine and should not be called from user code. */
981  void allocateAllValues(size_t numNonZeros, ArrayRCP<size_t> &rowptr, ArrayRCP<LocalOrdinal> &colind, ArrayRCP<Scalar> &values) {}
982 
984  void setAllValues(const ArrayRCP<size_t> &rowptr, const ArrayRCP<LocalOrdinal> &colind, const ArrayRCP<Scalar> &values) {}
985 
987  void getAllValues(ArrayRCP<const size_t> &rowptr, ArrayRCP<const LocalOrdinal> &colind, ArrayRCP<const Scalar> &values) const {}
988 
990  void getAllValues(ArrayRCP<Scalar> &values) {}
991 
992  bool haveGlobalConstants() const { return false; }
993 
995 
997 
998 
1000  void resumeFill(const RCP<ParameterList> &params = null) {}
1001 
1003  void fillComplete(const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap, const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap, const RCP<ParameterList> &params = null) {}
1004 
1006  void fillComplete(const RCP<ParameterList> &params = null) {}
1007 
1009  void replaceDomainMapAndImporter(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &newDomainMap, Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Node> > &newImporter) {}
1010 
1013  const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap,
1014  const RCP<const Import<LocalOrdinal, GlobalOrdinal, Node> > &importer = Teuchos::null,
1015  const RCP<const Export<LocalOrdinal, GlobalOrdinal, Node> > &exporter = Teuchos::null,
1016  const RCP<ParameterList> &params = Teuchos::null) {}
1017 
1019 
1021 
1022 
1024  const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > getRowMap() const { return Teuchos::null; }
1025 
1027  const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > getColMap() const { return Teuchos::null; }
1028 
1030  RCP<const CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > getCrsGraph() const { return Teuchos::null; }
1031 
1033  global_size_t getGlobalNumRows() const { return 0; }
1034 
1036  global_size_t getGlobalNumCols() const { return 0; }
1037 
1039  size_t getLocalNumRows() const { return 0; }
1040 
1042  size_t getLocalNumCols() const { return 0; }
1043 
1045  global_size_t getGlobalNumEntries() const { return 0; }
1046 
1048  size_t getLocalNumEntries() const { return 0; }
1049 
1051  size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const { return 0; }
1052 
1054  size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const { return 0; }
1055 
1057  size_t getGlobalMaxNumRowEntries() const { return 0; }
1058 
1060  size_t getLocalMaxNumRowEntries() const { return 0; }
1061 
1063  bool isLocallyIndexed() const { return false; }
1064 
1066  bool isGloballyIndexed() const { return false; }
1067 
1069  bool isFillComplete() const { return false; }
1070 
1072  bool isFillActive() const { return false; }
1073 
1075  typename ScalarTraits<Scalar>::magnitudeType getFrobeniusNorm() const { return ScalarTraits<Scalar>::magnitude(ScalarTraits<Scalar>::zero()); }
1076 
1078  bool supportsRowViews() const { return false; }
1079 
1081  void getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView<LocalOrdinal> &Indices, const ArrayView<Scalar> &Values, size_t &NumEntries) const {}
1082 
1084  void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView<const GlobalOrdinal> &indices, ArrayView<const Scalar> &values) const {}
1085 
1087  void getGlobalRowCopy(GlobalOrdinal GlobalRow, const ArrayView<GlobalOrdinal> &indices, const ArrayView<Scalar> &values, size_t &numEntries) const {}
1088 
1090  void getLocalRowView(LocalOrdinal LocalRow, ArrayView<const LocalOrdinal> &indices, ArrayView<const Scalar> &values) const {}
1091 
1093 
1095 
1096 
1098  void apply(const MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> &X, MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> &Y, Teuchos::ETransp mode = Teuchos::NO_TRANS, Scalar alpha = ScalarTraits<Scalar>::one(), Scalar beta = ScalarTraits<Scalar>::zero()) const {}
1099 
1101  void apply(const MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> &X, MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> &Y, Teuchos::ETransp mode, Scalar alpha, Scalar beta, bool sumInterfaceValues, const RCP<Xpetra::Import<LocalOrdinal, GlobalOrdinal, Node> > &regionInterfaceImporter, const Teuchos::ArrayRCP<LocalOrdinal> &regionInterfaceLIDs) const {}
1102 
1104  const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > getDomainMap() const { return Teuchos::null; }
1105 
1107  const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > getRangeMap() const { return Teuchos::null; }
1108 
1110 
1112 
1113 
1115  std::string description() const { return std::string(""); }
1116 
1118  void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel = Teuchos::Describable::verbLevel_default) const {}
1119 
1121 
1123 
1124  void setObjectLabel(const std::string &objectLabel) {}
1126 
1129 
1132 
1134  void getLocalDiagOffsets(Teuchos::ArrayRCP<size_t> &offsets) const {}
1135 
1137  void getLocalDiagCopy(Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node> &diag, const Teuchos::ArrayView<const size_t> &offsets) const {}
1138 
1141 
1144 
1147 
1149  //{@
1150 
1152  Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > getMap() const { return Teuchos::null; }
1153 
1157 
1161 
1165 
1169 
1170  void removeEmptyProcessesInPlace(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &newMap) {}
1171 
1172  // @}
1174 
1175 
1177  bool hasMatrix() const { return false; }
1178 
1180  TpetraCrsMatrix(const Teuchos::RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &mtx) {
1182  }
1183 
1185  RCP<const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > getTpetra_CrsMatrix() const { return Teuchos::null; }
1186 
1188  RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > getTpetra_CrsMatrixNonConst() const { return Teuchos::null; } // TODO: remove
1189 
1190 #ifdef HAVE_XPETRA_TPETRA
1191  local_matrix_type getLocalMatrix() const {
1193  TEUCHOS_UNREACHABLE_RETURN(local_matrix_type());
1194  }
1195 
1197  typename local_matrix_type::HostMirror getLocalMatrixHost() const {
1198  TEUCHOS_UNREACHABLE_RETURN(typename local_matrix_type::HostMirror());
1199  }
1202  TEUCHOS_UNREACHABLE_RETURN(local_matrix_type());
1203  }
1204 
1205  void setAllValues(const typename local_matrix_type::row_map_type &ptr,
1206  const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type &ind,
1207  const typename local_matrix_type::values_type &val) {}
1208 #endif
1209 
1211  LocalOrdinal GetStorageBlockSize() const { return 1; }
1212 
1217 
1219 }; // TpetraCrsMatrix class (specialization for GO=long long, NO=EpetraNode)
1220 #endif
1221 
1222 #endif // HAVE_XPETRA_EPETRA
1223 
1224 } // namespace Xpetra
1225 
1226 #define XPETRA_TPETRACRSMATRIX_SHORT
1227 #endif // XPETRA_TPETRACRSMATRIX_DECL_HPP
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import.
void getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView< LocalOrdinal > &Indices, const ArrayView< Scalar > &Values, size_t &NumEntries) const
Extract a list of entries in a specified local row of the matrix. Put into storage allocated by calli...
LocalOrdinal GetStorageBlockSize() const
Returns the block size of the storage mechanism, which is usually 1, except for Tpetra::BlockCrsMatri...
size_t getGlobalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on all nodes.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const
Returns the Map that describes the row distribution in this matrix.
global_size_t getGlobalNumRows() const
Number of global elements in the row map of this matrix.
void replaceDiag(const Vector &diag)
Replace the diagonal entries of the matrix.
RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrixNonConst() const
Get the underlying Tpetra matrix.
TpetraCrsMatrix(const Teuchos::RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying a previously constructed graph.
std::string description() const
A simple one-line description of this object.
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap=Teuchos::null, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap=Teuchos::null, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor for a fused export.
size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const
Returns the current number of entries on this node in the specified local row.
void fillComplete(const RCP< ParameterList > &params=null)
Signal that data entry is complete.
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import.
Xpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::local_matrix_type local_matrix_type
void apply(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp mode, Scalar alpha, Scalar beta, bool sumInterfaceValues, const RCP< Xpetra::Import< LocalOrdinal, GlobalOrdinal, Node > > &regionInterfaceImporter, const Teuchos::ArrayRCP< LocalOrdinal > &regionInterfaceLIDs) const
Computes the matrix-multivector multiplication for region layout matrices.
void setAllToScalar(const Scalar &alpha)
Set all matrix entries equal to scalarThis.
void leftScale(const Vector &x)
Left scale operator with given vector values.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const
Returns the Map associated with the domain of this operator. This will be null until fillComplete() i...
void removeEmptyProcessesInPlace(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newMap)
void replaceDomainMapAndImporter(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newDomainMap, Teuchos::RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > &newImporter)
Replaces the current domainMap and importer with the user-specified objects.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
void replaceLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using local IDs.
global_size_t getGlobalNumCols() const
Number of global columns in the matrix.
LocalOrdinal GetStorageBlockSize() const
Returns the block size of the storage mechanism, which is usually 1, except for Tpetra::BlockCrsMatri...
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using local IDs.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const
Returns the Map that describes the column distribution in this matrix.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const
Returns the Map associated with the range of this operator, which must be compatible with Y...
void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using global IDs.
TpetraCrsMatrix(const Teuchos::RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &mtx)
TpetraCrsMatrix constructor to wrap a Tpetra::CrsMatrix object.
bool isLocallyIndexed() const
If matrix indices are in the local range, this function returns true. Otherwise, this function return...
global_size_t getGlobalNumEntries() const
Returns the global number of entries in this matrix.
RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > mtx_
TpetraImport< LocalOrdinal, GlobalOrdinal, Node > TpetraImportClass
void getLocalDiagOffsets(Teuchos::ArrayRCP< size_t > &offsets) const
Get offsets of the diagonal entries in the matrix.
void getLocalRowView(LocalOrdinal LocalRow, ArrayView< const LocalOrdinal > &indices, ArrayView< const Scalar > &values) const
Extract a const, non-persisting view of local indices in a specified row of the matrix.
void getAllValues(ArrayRCP< Scalar > &values)
Gets the 1D pointer arrays of the graph.
size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const
Returns the current number of entries in the (locally owned) global row.
void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag) const
Get a copy of the diagonal entries owned by this node, with local row idices.
TpetraCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraCrsMatrixClass
void residual(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &R) const
Compute a residual R = B - (*this) * X.
void setAllValues(const ArrayRCP< size_t > &rowptr, const ArrayRCP< LocalOrdinal > &colind, const ArrayRCP< Scalar > &values)
Sets the 1D pointer arrays of the graph.
void setAllValues(const typename local_matrix_type::row_map_type &ptr, const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type &ind, const typename local_matrix_type::values_type &val)
size_t getLocalNumRows() const
Returns the number of matrix rows owned on the calling node.
size_t getLocalNumEntries() const
Returns the local number of entries in this matrix.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const
Returns the Map associated with the range of this operator, which must be compatible with Y...
void rightScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Right scale operator with given vector values.
ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const
Returns the Frobenius norm of the matrix.
void allocateAllValues(size_t numNonZeros, ArrayRCP< size_t > &rowptr, ArrayRCP< LocalOrdinal > &colind, ArrayRCP< Scalar > &values)
Allocates and returns ArrayRCPs of the Crs arrays — This is an Xpetra-only routine.
void replaceLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using local IDs.
void resumeFill(const RCP< ParameterList > &params=null)
void replaceLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using local IDs.
bool isLocallyIndexed() const
If matrix indices are in the local range, this function returns true. Otherwise, this function return...
TpetraCrsMatrix(const Teuchos::RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying a previously constructed graph.
void getLocalRowView(LocalOrdinal LocalRow, ArrayView< const LocalOrdinal > &indices, ArrayView< const Scalar > &values) const
Extract a const, non-persisting view of local indices in a specified row of the matrix.
TpetraExport< LocalOrdinal, GlobalOrdinal, Node > TpetraExportClass
void setObjectLabel(const std::string &objectLabel)
local_matrix_type::HostMirror getLocalMatrixHost() const
Access the local Kokkos::CrsMatrix data.
bool haveGlobalConstants() const
Returns true if globalConstants have been computed; false otherwise.
void setAllValues(const ArrayRCP< size_t > &rowptr, const ArrayRCP< LocalOrdinal > &colind, const ArrayRCP< Scalar > &values)
Sets the 1D pointer arrays of the graph.
void getGlobalRowCopy(GlobalOrdinal GlobalRow, const ArrayView< GlobalOrdinal > &indices, const ArrayView< Scalar > &values, size_t &numEntries) const
Extract a list of entries in a specified global row of this matrix. Put into pre-allocated storage...
void removeEmptyProcessesInPlace(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newMap)
global_size_t getGlobalNumRows() const
Number of global elements in the row map of this matrix.
RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > getCrsGraph() const
Returns the CrsGraph associated with this matrix.
void getAllValues(ArrayRCP< Scalar > &values)
Gets the 1D pointer arrays of the graph.
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Export (using an Importer).
void scale(const Scalar &alpha)
Scale the current values of a matrix, this = alpha*this.
#define XPETRA_TPETRA_ETI_EXCEPTION(cl, obj, go, node)
local_matrix_type::HostMirror getLocalMatrixHost() const
Access the local Kokkos::CrsMatrix data.
local_matrix_type getLocalMatrixDevice() const
Access the local Kokkos::CrsMatrix data.
void fillComplete(const RCP< ParameterList > &params=null)
Signal that data entry is complete.
void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag, const Teuchos::ArrayView< const size_t > &offsets) const
Get a copy of the diagonal entries owned by this node, with local row indices.
void getLocalDiagCopy(Vector &diag) const
Get a copy of the diagonal entries owned by this node, with local row indices.
void leftScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Left scale operator with given vector values.
std::string description() const
A simple one-line description of this object.
void getAllValues(ArrayRCP< const size_t > &rowptr, ArrayRCP< const LocalOrdinal > &colind, ArrayRCP< const Scalar > &values) const
Gets the 1D pointer arrays of the graph.
size_t getLocalNumCols() const
Returns the number of columns connected to the locally owned rows of this matrix. ...
void replaceDomainMapAndImporter(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newDomainMap, Teuchos::RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > &newImporter)
Replaces the current domainMap and importer with the user-specified objects.
void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using local IDs.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const
Returns the Map that describes the row distribution in this matrix.
bool isFillComplete() const
Returns true if the matrix is in compute mode, i.e. if fillComplete() has been called.
bool isLocallyIndexed() const
If matrix indices are in the local range, this function returns true. Otherwise, this function return...
void setAllToScalar(const Scalar &alpha)
Set all matrix entries equal to scalarThis.
Xpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::local_matrix_type local_matrix_type
void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using local IDs.
void getLocalDiagOffsets(Teuchos::ArrayRCP< size_t > &offsets) const
Get offsets of the diagonal entries in the matrix.
std::string description() const
A simple one-line description of this object.
local_matrix_type getLocalMatrixDevice() const
Access the local Kokkos::CrsMatrix data.
bool isFillComplete() const
Returns true if the matrix is in compute mode, i.e. if fillComplete() has been called.
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap=Teuchos::null, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap=Teuchos::null, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor for a fused import.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const
Returns the Map associated with the domain of this operator. This will be null until fillComplete() i...
bool haveGlobalConstants() const
Returns true if globalConstants have been computed; false otherwise.
void residual(const MultiVector &X, const MultiVector &B, MultiVector &R) const
Compute a residual R = B - (*this) * X.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying (possibly different) number of entries in each row.
size_t getLocalNumRows() const
Returns the number of matrix rows owned on the calling node.
void fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > &params=null)
Signal that data entry is complete, specifying domain and range maps.
RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrix() const
Get the underlying Tpetra matrix.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
void getLocalDiagOffsets(Teuchos::ArrayRCP< size_t > &offsets) const
Get offsets of the diagonal entries in the matrix.
bool isGloballyIndexed() const
If matrix indices are in the global range, this function returns true. Otherwise, this function retur...
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Export< LocalOrdinal, GlobalOrdinal, Node > &RowExporter, const Teuchos::RCP< const Export< LocalOrdinal, GlobalOrdinal, Node > > DomainExporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > &params)
Constructor for a fused export.
void setAllValues(const typename local_matrix_type::row_map_type &ptr, const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type &ind, const typename local_matrix_type::values_type &val)
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const
Returns the Map associated with the domain of this operator. This will be null until fillComplete() i...
global_size_t getGlobalNumRows() const
Number of global elements in the row map of this matrix.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const
Returns the Map that describes the column distribution in this matrix.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const
Returns the current number of entries in the (locally owned) global row.
TpetraExport< LocalOrdinal, GlobalOrdinal, Node > TpetraExportClass
RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > getCrsGraph() const
Returns the CrsGraph associated with this matrix.
void residual(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &R) const
Compute a residual R = B - (*this) * X.
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
void getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView< LocalOrdinal > &Indices, const ArrayView< Scalar > &Values, size_t &NumEntries) const
Extract a list of entries in a specified local row of the matrix. Put into storage allocated by calli...
TpetraCrsMatrix(const Teuchos::RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &mtx)
TpetraCrsMatrix constructor to wrap a Tpetra::CrsMatrix object.
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Import (using an Exporter).
void apply(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=ScalarTraits< Scalar >::one(), Scalar beta=ScalarTraits< Scalar >::zero()) const
Computes the sparse matrix-multivector multiplication.
void removeEmptyProcessesInPlace(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newMap)
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap=Teuchos::null, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap=Teuchos::null, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor for a fused import.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, size_t maxNumEntriesPerRow, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying column Map and fixed number of entries for each row.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with some verbosity level to an FancyOStream object.
size_t getLocalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on this node.
TpetraCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraCrsMatrixClass
TpetraCrsMatrix(const local_matrix_type &lclMatrix, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap=Teuchos::null, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap=Teuchos::null, const Teuchos::RCP< Teuchos::ParameterList > &params=null)
Constructor specifying local matrix and 4 maps.
void expertStaticFillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > &importer=Teuchos::null, const RCP< const Export< LocalOrdinal, GlobalOrdinal, Node > > &exporter=Teuchos::null, const RCP< ParameterList > &params=Teuchos::null)
Expert static fill complete.
void allocateAllValues(size_t numNonZeros, ArrayRCP< size_t > &rowptr, ArrayRCP< LocalOrdinal > &colind, ArrayRCP< Scalar > &values)
Allocates and returns ArrayRCPs of the Crs arrays — This is an Xpetra-only routine.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying column Map and number of entries in each row.
size_t getLocalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on this node.
void expertStaticFillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > &importer=Teuchos::null, const RCP< const Export< LocalOrdinal, GlobalOrdinal, Node > > &exporter=Teuchos::null, const RCP< ParameterList > &params=Teuchos::null)
Expert static fill complete.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, size_t maxNumEntriesPerRow, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying column Map and fixed number of entries for each row.
void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView< const GlobalOrdinal > &indices, ArrayView< const Scalar > &values) const
Extract a const, non-persisting view of global indices in a specified row of the matrix.
void resumeFill(const RCP< ParameterList > &params=null)
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with some verbosity level to an FancyOStream object.
void scale(const Scalar &alpha)
Scale the current values of a matrix, this = alpha*this.
void replaceDomainMapAndImporter(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newDomainMap, Teuchos::RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > &newImporter)
Replaces the current domainMap and importer with the user-specified objects.
size_t getGlobalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on all nodes.
void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView< const GlobalOrdinal > &indices, ArrayView< const Scalar > &values) const
Extract a const, non-persisting view of global indices in a specified row of the matrix.
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Export< LocalOrdinal, GlobalOrdinal, Node > &RowExporter, const Teuchos::RCP< const Export< LocalOrdinal, GlobalOrdinal, Node > > DomainExporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > &params)
Constructor for a fused export.
TpetraCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraCrsMatrixClass
size_t global_size_t
Global size_t object.
TpetraVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraVectorClass
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Export (using an Importer).
void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using global IDs.
global_size_t getGlobalNumCols() const
Number of global columns in the matrix.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, size_t maxNumEntriesPerRow, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying fixed number of entries for each row.
virtual ~TpetraCrsMatrix()
Destructor.
TpetraVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraVectorClass
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const
Returns the Map associated with the range of this operator, which must be compatible with Y...
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap=Teuchos::null, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap=Teuchos::null, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor for a fused export.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, size_t maxNumEntriesPerRow, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying fixed number of entries for each row.
bool hasMatrix() const
Does this have an underlying matrix.
void getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView< LocalOrdinal > &Indices, const ArrayView< Scalar > &Values, size_t &NumEntries) const
Extract a list of entries in a specified local row of the matrix. Put into storage allocated by calli...
bool isGloballyIndexed() const
If matrix indices are in the global range, this function returns true. Otherwise, this function retur...
void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView< const GlobalOrdinal > &indices, ArrayView< const Scalar > &values) const
Extract a const, non-persisting view of global indices in a specified row of the matrix.
void setAllValues(const ArrayRCP< size_t > &rowptr, const ArrayRCP< LocalOrdinal > &colind, const ArrayRCP< Scalar > &values)
Sets the 1D pointer arrays of the graph.
size_t getLocalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on this node.
bool isFillComplete() const
Returns true if the matrix is in compute mode, i.e. if fillComplete() has been called.
Tpetra::KokkosCompat::KokkosSerialWrapperNode EpetraNode
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with some verbosity level to an FancyOStream object.
TpetraCrsMatrix(const local_matrix_type &lclMatrix, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap=Teuchos::null, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap=Teuchos::null, const Teuchos::RCP< Teuchos::ParameterList > &params=null)
Constructor specifying local matrix and 4 maps.
void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using global IDs.
void rightScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Right scale operator with given vector values.
RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrix() const
Get the underlying Tpetra matrix.
void getAllValues(ArrayRCP< const size_t > &rowptr, ArrayRCP< const LocalOrdinal > &colind, ArrayRCP< const Scalar > &values) const
Gets the 1D pointer arrays of the graph.
bool hasMatrix() const
Does this have an underlying matrix.
void apply(const MultiVector &X, MultiVector &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=ScalarTraits< Scalar >::one(), Scalar beta=ScalarTraits< Scalar >::zero()) const
Computes the sparse matrix-multivector multiplication.
global_size_t getGlobalNumEntries() const
Returns the global number of entries in this matrix.
size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const
Returns the current number of entries in the (locally owned) global row.
bool haveGlobalConstants() const
Returns true if globalConstants have been computed; false otherwise.
global_size_t getGlobalNumCols() const
Number of global columns in the matrix.
bool supportsRowViews() const
Returns true if getLocalRowView() and getGlobalRowView() are valid for this class.
global_size_t getGlobalNumEntries() const
Returns the global number of entries in this matrix.
void scale(const Scalar &alpha)
Scale the current values of a matrix, this = alpha*this.
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Import< LocalOrdinal, GlobalOrdinal, Node > &RowImporter, const Teuchos::RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > DomainImporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > &params)
Constructor for a fused import.
RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrixNonConst() const
Get the underlying Tpetra matrix.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, size_t maxNumEntriesPerRow, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying fixed number of entries for each row.
void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using global IDs.
size_t getLocalNumRows() const
Returns the number of matrix rows owned on the calling node.
void replaceDiag(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag)
Replace the diagonal entries of the matrix.
void expertStaticFillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > &importer=Teuchos::null, const RCP< const Export< LocalOrdinal, GlobalOrdinal, Node > > &exporter=Teuchos::null, const RCP< ParameterList > &params=Teuchos::null)
Expert static fill complete.
bool hasMatrix() const
Does this have an underlying matrix.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const
Returns the Map that describes the row distribution in this matrix.
TpetraImport< LocalOrdinal, GlobalOrdinal, Node > TpetraImportClass
bool supportsRowViews() const
Returns true if getLocalRowView() and getGlobalRowView() are valid for this class.
CombineMode
Xpetra::Combine Mode enumerable type.
TpetraCrsMatrix(const TpetraCrsMatrix &matrix)
Deep copy constructor.
void apply(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=ScalarTraits< Scalar >::one(), Scalar beta=ScalarTraits< Scalar >::zero()) const
Computes the sparse matrix-multivector multiplication.
size_t getLocalNumCols() const
Returns the number of columns connected to the locally owned rows of this matrix. ...
RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrix() const
Get the underlying Tpetra matrix.
void fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > &params=null)
Signal that data entry is complete, specifying domain and range maps.
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Import< LocalOrdinal, GlobalOrdinal, Node > &RowImporter, const Teuchos::RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > DomainImporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > &params)
Constructor for a fused import.
size_t getLocalNumEntries() const
Returns the local number of entries in this matrix.
void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag) const
Get a copy of the diagonal entries owned by this node, with local row idices.
void getLocalRowView(LocalOrdinal LocalRow, ArrayView< const LocalOrdinal > &indices, ArrayView< const Scalar > &values) const
Extract a const, non-persisting view of local indices in a specified row of the matrix.
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
void getAllValues(ArrayRCP< const size_t > &rowptr, ArrayRCP< const LocalOrdinal > &colind, ArrayRCP< const Scalar > &values) const
Gets the 1D pointer arrays of the graph.
void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using global IDs.
void allocateAllValues(size_t numNonZeros, ArrayRCP< size_t > &rowptr, ArrayRCP< LocalOrdinal > &colind, ArrayRCP< Scalar > &values)
Allocates and returns ArrayRCPs of the Crs arrays — This is an Xpetra-only routine.
KokkosSparse::CrsMatrix< impl_scalar_type, LocalOrdinal, execution_space, void, typename local_graph_type::size_type > local_matrix_type
The specialization of Kokkos::CrsMatrix that represents the part of the sparse matrix on each MPI pro...
void setAllValues(const typename local_matrix_type::row_map_type &ptr, const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type &ind, const typename local_matrix_type::values_type &val)
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const
Returns the Map that describes the column distribution in this matrix.
void leftScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Left scale operator with given vector values.
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import.
TpetraExport< LocalOrdinal, GlobalOrdinal, Node > TpetraExportClass
bool isGloballyIndexed() const
If matrix indices are in the global range, this function returns true. Otherwise, this function retur...
void replaceDiag(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag)
Replace the diagonal entries of the matrix.
bool isFillActive() const
Returns true if the matrix is in edit mode.
RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > getCrsGraph() const
Returns the CrsGraph associated with this matrix.
size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const
Returns the current number of entries on this node in the specified local row.
void fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > &params=null)
Signal that data entry is complete, specifying domain and range maps.
size_t getGlobalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on all nodes.
LocalOrdinal GetStorageBlockSize() const
Returns the block size of the storage mechanism, which is usually 1, except for Tpetra::BlockCrsMatri...
size_t getLocalNumCols() const
Returns the number of columns connected to the locally owned rows of this matrix. ...
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying (possibly different) number of entries in each row.
ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const
Returns the Frobenius norm of the matrix.
ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const
Returns the Frobenius norm of the matrix.
bool isFillActive() const
Returns true if the matrix is in edit mode.
TpetraCrsMatrix(const TpetraCrsMatrix &matrix)
Deep copy constructor.
RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrixNonConst() const
Get the underlying Tpetra matrix.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying column Map and number of entries in each row.
void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using global IDs.
TpetraImport< LocalOrdinal, GlobalOrdinal, Node > TpetraImportClass
void getGlobalRowCopy(GlobalOrdinal GlobalRow, const ArrayView< GlobalOrdinal > &indices, const ArrayView< Scalar > &values, size_t &numEntries) const
Extract a list of entries in a specified global row of this matrix. Put into pre-allocated storage...
void getGlobalRowCopy(GlobalOrdinal GlobalRow, const ArrayView< GlobalOrdinal > &indices, const ArrayView< Scalar > &values, size_t &numEntries) const
Extract a list of entries in a specified global row of this matrix. Put into pre-allocated storage...
void setAllToScalar(const Scalar &alpha)
Set all matrix entries equal to scalarThis.
size_t getLocalNumEntries() const
Returns the local number of entries in this matrix.
bool isFillActive() const
Returns true if the matrix is in edit mode.
bool supportsRowViews() const
Returns true if getLocalRowView() and getGlobalRowView() are valid for this class.
void apply(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp mode, Scalar alpha, Scalar beta, bool sumInterfaceValues, const RCP< Xpetra::Import< LocalOrdinal, GlobalOrdinal, Node > > &regionInterfaceImporter, const Teuchos::ArrayRCP< LocalOrdinal > &regionInterfaceLIDs) const
Computes the matrix-multivector multiplication for region layout matrices.
void rightScale(const Vector &x)
Right scale operator with given vector values.
Xpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::local_matrix_type local_matrix_type
size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const
Returns the current number of entries on this node in the specified local row.
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Import (using an Exporter).
void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag, const Teuchos::ArrayView< const size_t > &offsets) const
Get a copy of the diagonal entries owned by this node, with local row indices.