46 #ifndef XPETRA_TPETRACRSMATRIX_DECL_HPP
47 #define XPETRA_TPETRACRSMATRIX_DECL_HPP
57 #include "Tpetra_CrsMatrix.hpp"
58 #include "Tpetra_replaceDiagonalCrsMatrix.hpp"
69 template<
class Scalar,
74 :
public CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>
84 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
85 #ifdef HAVE_XPETRA_TPETRA
115 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null);
122 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null);
130 const Teuchos::RCP<Teuchos::ParameterList>& params);
138 const Teuchos::RCP<Teuchos::ParameterList>& params);
140 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
141 #ifdef HAVE_XPETRA_TPETRA
164 const local_matrix_type& lclMatrix,
165 const Teuchos::RCP<Teuchos::ParameterList>& params = null);
169 const local_matrix_type& lclMatrix,
174 const Teuchos::RCP<Teuchos::ParameterList>& params = null);
187 void insertGlobalValues(GlobalOrdinal globalRow,
const ArrayView< const GlobalOrdinal > &cols,
const ArrayView< const Scalar > &vals);
190 void insertLocalValues(LocalOrdinal localRow,
const ArrayView< const LocalOrdinal > &cols,
const ArrayView< const Scalar > &vals);
193 void replaceGlobalValues(GlobalOrdinal globalRow,
const ArrayView< const GlobalOrdinal > &cols,
const ArrayView< const Scalar > &vals);
198 const ArrayView<const LocalOrdinal> &cols,
199 const ArrayView<const Scalar> &vals);
205 void scale(
const Scalar &alpha);
209 void allocateAllValues(
size_t numNonZeros,ArrayRCP<size_t> & rowptr, ArrayRCP<LocalOrdinal> & colind, ArrayRCP<Scalar> & values)
213 void setAllValues(
const ArrayRCP<size_t> & rowptr,
const ArrayRCP<LocalOrdinal> & colind,
const ArrayRCP<Scalar> & values)
217 void getAllValues(ArrayRCP<const size_t>& rowptr, ArrayRCP<const LocalOrdinal>& colind, ArrayRCP<const Scalar>& values)
const
229 void resumeFill(
const RCP< ParameterList > ¶ms=null);
235 void fillComplete(
const RCP< ParameterList > ¶ms=null);
246 const RCP<ParameterList> ¶ms=Teuchos::null);
254 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getRowMap()
const;
257 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getColMap()
const;
260 RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > >
getCrsGraph()
const;
311 void getLocalRowCopy(LocalOrdinal LocalRow,
const ArrayView< LocalOrdinal > &Indices,
const ArrayView< Scalar > &Values,
size_t &NumEntries)
const;
314 void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView< const GlobalOrdinal > &indices, ArrayView< const Scalar > &values)
const;
317 void getGlobalRowCopy(GlobalOrdinal GlobalRow,
const ArrayView< GlobalOrdinal > &indices,
const ArrayView< Scalar > &values,
size_t &numEntries)
const;
320 void getLocalRowView(LocalOrdinal LocalRow, ArrayView< const LocalOrdinal > &indices, ArrayView< const Scalar > &values)
const;
328 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;
331 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getDomainMap()
const;
334 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getRangeMap()
const;
345 void describe(Teuchos::FancyOStream &out,
const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default)
const;
381 Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getMap()
const;
409 TpetraCrsMatrix(
const Teuchos::RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &mtx);
412 RCP<const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
getTpetra_CrsMatrix()
const;
417 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
418 #ifdef HAVE_XPETRA_TPETRA
419 local_matrix_type getLocalMatrix ()
const {
424 void setAllValues (
const typename local_matrix_type::row_map_type& ptr,
425 const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type& ind,
426 const typename local_matrix_type::values_type& val) {
433 void residual(
const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & X,
434 const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & B,
435 MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & R)
const;
440 RCP< Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
mtx_;
443 #ifdef HAVE_XPETRA_EPETRA
445 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
446 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
449 template<
class Scalar>
451 :
public CrsMatrix<Scalar,int,int,EpetraNode>
463 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
464 #ifdef HAVE_XPETRA_TPETRA
506 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null) {
515 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null) {
525 const Teuchos::RCP<Teuchos::ParameterList>& params)
536 const Teuchos::RCP<Teuchos::ParameterList>& params)
541 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
542 #ifdef HAVE_XPETRA_TPETRA
565 const local_matrix_type& lclMatrix,
566 const Teuchos::RCP<Teuchos::ParameterList>& params = null) {
572 const local_matrix_type& lclMatrix,
573 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& rowMap,
574 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& colMap,
575 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap = Teuchos::null,
576 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& rangeMap = Teuchos::null,
577 const Teuchos::RCP<Teuchos::ParameterList>& params = null) {
603 const ArrayView<const LocalOrdinal> &cols,
604 const ArrayView<const Scalar> &vals)
615 void allocateAllValues(
size_t numNonZeros,ArrayRCP<size_t> & rowptr, ArrayRCP<LocalOrdinal> & colind, ArrayRCP<Scalar> & values) { }
618 void setAllValues(
const ArrayRCP<size_t> & rowptr,
const ArrayRCP<LocalOrdinal> & colind,
const ArrayRCP<Scalar> & values) { }
621 void getAllValues(ArrayRCP<const size_t>& rowptr, ArrayRCP<const LocalOrdinal>& colind, ArrayRCP<const Scalar>& values)
const { }
648 const RCP<ParameterList> ¶ms=Teuchos::null) { }
656 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getRowMap()
const {
return Teuchos::null; }
659 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getColMap()
const {
return Teuchos::null; }
662 RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > >
getCrsGraph()
const {
return Teuchos::null; }
707 typename ScalarTraits< Scalar >::magnitudeType
getFrobeniusNorm()
const {
return ScalarTraits< Scalar >::magnitude(ScalarTraits< Scalar >::zero()); }
713 void getLocalRowCopy(
LocalOrdinal LocalRow,
const ArrayView< LocalOrdinal > &Indices,
const ArrayView< Scalar > &Values,
size_t &NumEntries)
const { }
730 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 { }
733 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getDomainMap()
const {
return Teuchos::null; }
736 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getRangeMap()
const {
return Teuchos::null; }
747 void describe(Teuchos::FancyOStream &out,
const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default)
const { }
782 Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getMap()
const {
return Teuchos::null; }
810 TpetraCrsMatrix(
const Teuchos::RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &mtx) {
815 RCP<const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
getTpetra_CrsMatrix()
const {
return Teuchos::null; }
820 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
821 #ifdef HAVE_XPETRA_TPETRA
822 local_matrix_type getLocalMatrix ()
const {
824 TEUCHOS_UNREACHABLE_RETURN(local_matrix_type());
827 void setAllValues (
const typename local_matrix_type::row_map_type& ptr,
828 const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type& ind,
829 const typename local_matrix_type::values_type& val) { }
842 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_LONG_LONG))) || \
843 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_LONG_LONG))))
846 template<
class Scalar>
848 :
public CrsMatrix<Scalar,int,long long,EpetraNode>
860 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
861 #ifdef HAVE_XPETRA_TPETRA
903 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null) {
912 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null) {
922 const Teuchos::RCP<Teuchos::ParameterList>& params)
933 const Teuchos::RCP<Teuchos::ParameterList>& params)
938 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
939 #ifdef HAVE_XPETRA_TPETRA
962 const local_matrix_type& lclMatrix,
963 const Teuchos::RCP<Teuchos::ParameterList>& params = null) {
969 const local_matrix_type& lclMatrix,
970 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& rowMap,
971 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& colMap,
972 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap = Teuchos::null,
973 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& rangeMap = Teuchos::null,
974 const Teuchos::RCP<Teuchos::ParameterList>& params = null) {
975 XPETRA_TPETRA_ETI_EXCEPTION(
typeid(TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,EpetraNode>).name() ,
typeid(TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,EpetraNode>).name(),
"long long",
typeid(
EpetraNode).name() );
1000 const ArrayView<const LocalOrdinal> &cols,
1001 const ArrayView<const Scalar> &vals)
1012 void allocateAllValues(
size_t numNonZeros,ArrayRCP<size_t> & rowptr, ArrayRCP<LocalOrdinal> & colind, ArrayRCP<Scalar> & values) { }
1015 void setAllValues(
const ArrayRCP<size_t> & rowptr,
const ArrayRCP<LocalOrdinal> & colind,
const ArrayRCP<Scalar> & values) { }
1018 void getAllValues(ArrayRCP<const size_t>& rowptr, ArrayRCP<const LocalOrdinal>& colind, ArrayRCP<const Scalar>& values)
const { }
1045 const RCP<ParameterList> ¶ms=Teuchos::null) { }
1053 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getRowMap()
const {
return Teuchos::null; }
1056 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getColMap()
const {
return Teuchos::null; }
1059 RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > >
getCrsGraph()
const {
return Teuchos::null; }
1104 typename ScalarTraits< Scalar >::magnitudeType
getFrobeniusNorm()
const {
return ScalarTraits< Scalar >::magnitude(ScalarTraits< Scalar >::zero()); }
1110 void getLocalRowCopy(
LocalOrdinal LocalRow,
const ArrayView< LocalOrdinal > &Indices,
const ArrayView< Scalar > &Values,
size_t &NumEntries)
const { }
1127 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 { }
1130 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getDomainMap()
const {
return Teuchos::null; }
1133 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getRangeMap()
const {
return Teuchos::null; }
1144 void describe(Teuchos::FancyOStream &out,
const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default)
const { }
1178 Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getMap()
const {
return Teuchos::null; }
1206 TpetraCrsMatrix(
const Teuchos::RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &mtx) {
1211 RCP<const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
getTpetra_CrsMatrix()
const {
return Teuchos::null; }
1216 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
1217 #ifdef HAVE_XPETRA_TPETRA
1218 local_matrix_type getLocalMatrix ()
const {
1220 TEUCHOS_UNREACHABLE_RETURN(local_matrix_type());
1223 void setAllValues (
const typename local_matrix_type::row_map_type& ptr,
1224 const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type& ind,
1225 const typename local_matrix_type::values_type& val) { }
1238 #endif // HAVE_XPETRA_EPETRA
1242 #define XPETRA_TPETRACRSMATRIX_SHORT
1243 #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...
size_t getGlobalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on all nodes.
global_size_t getGlobalNumRows() const
Number of global elements in the row map of this 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 > ¶ms=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 > ¶ms=Teuchos::null)
Constructor for a fused export.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
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 > ¶ms=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.
void setAllToScalar(const Scalar &alpha)
Set all matrix entries equal to scalarThis.
TpetraExport< LocalOrdinal, GlobalOrdinal, Node > TpetraExportClass
TpetraVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraVectorClass
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)
TpetraImport< LocalOrdinal, GlobalOrdinal, Node > TpetraImportClass
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 getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag) const
Get a copy of the diagonal entries owned by this node, with local row indices.
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.
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 insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using local IDs.
TpetraExport< LocalOrdinal, GlobalOrdinal, Node > TpetraExportClass
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 replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using global IDs.
TpetraCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraCrsMatrixClass
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.
size_t getNodeNumRows() const
Returns the number of matrix rows owned on the calling node.
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.
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.
RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > mtx_
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 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 rightScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Right scale operator with given vector values.
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 > ¶ms=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...
size_t getNodeNumEntries() const
Returns the local number of entries in this matrix.
void setObjectLabel(const std::string &objectLabel)
void setObjectLabel(const std::string &objectLabel)
TpetraCrsMatrix(const Teuchos::RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=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.
size_t getNodeNumEntries() const
Returns the local number of entries in this matrix.
void setObjectLabel(const std::string &objectLabel)
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)
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const
Returns the Map that describes the row distribution in this matrix.
virtual ~TpetraCrsMatrix()
Destructor.
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 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)
void fillComplete(const RCP< ParameterList > ¶ms=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.
ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const
Returns the Frobenius norm of the matrix.
size_t getNodeNumCols() const
Returns the number of columns connected to the locally owned rows of this matrix. ...
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.
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.
size_t getNodeMaxNumRowEntries() 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.
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.
void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using local IDs.
TpetraVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraVectorClass
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const
Returns the Map that describes the row distribution in this matrix.
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.
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 > ¶ms=Teuchos::null)
Constructor for a fused import.
TpetraImport< LocalOrdinal, GlobalOrdinal, Node > TpetraImportClass
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...
TpetraCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraCrsMatrixClass
bool haveGlobalConstants() const
Returns true if globalConstants have been computed; false otherwise.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying (possibly different) number of entries in each row.
void fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > ¶ms=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.
RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > getCrsGraph() const
Returns the CrsGraph associated with this 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 > ¶ms)
Constructor for a fused export.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const
Returns the Map that describes the column distribution in this matrix.
global_size_t getGlobalNumRows() const
Number of global elements in the row map of this matrix.
size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const
Returns the current number of entries in the (locally owned) global row.
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 > ¶ms=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 > ¶ms=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.
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 > ¶ms=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.
RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > getCrsGraph() const
Returns the CrsGraph associated with this 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 > ¶ms=Teuchos::null)
Constructor specifying column Map and number of entries in each row.
size_t getNodeNumCols() const
Returns the number of columns connected to the locally owned rows of this 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 > ¶ms=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 > ¶ms=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.
TpetraVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraVectorClass
void resumeFill(const RCP< ParameterList > ¶ms=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 rightScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Right scale operator with given vector values.
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 > ¶ms)
Constructor for a fused export.
size_t global_size_t
Global size_t object.
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 > ¶ms=Teuchos::null)
Constructor specifying fixed number of entries for each row.
TpetraImport< LocalOrdinal, GlobalOrdinal, Node > TpetraImportClass
virtual ~TpetraCrsMatrix()
Destructor.
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 > ¶ms=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 > ¶ms=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...
size_t getNodeMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on this node.
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.
virtual ~TpetraCrsMatrix()
Destructor.
void setAllValues(const ArrayRCP< size_t > &rowptr, const ArrayRCP< LocalOrdinal > &colind, const ArrayRCP< Scalar > &values)
Sets the 1D pointer arrays of the graph.
bool isFillComplete() const
Returns true if the matrix is in compute mode, i.e. if fillComplete() has been called.
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.
ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const
Returns the Frobenius norm of the matrix.
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.
void resumeFill(const RCP< ParameterList > ¶ms=null)
bool hasMatrix() const
Does this have an underlying matrix.
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.
size_t getNodeNumEntries() const
Returns the local 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 > ¶ms)
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 > ¶ms=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.
void replaceDiag(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag)
Replace the diagonal entries of the 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 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 > ¶ms=Teuchos::null)
Expert static fill complete.
bool hasMatrix() const
Does this have an underlying matrix.
void replaceDiag(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag)
Replace the diagonal entries of the matrix.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const
Returns the Map that describes the row distribution in this matrix.
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.
void leftScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Left scale operator with given vector values.
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 > ¶ms=null)
Signal that data entry is complete, specifying domain and range maps.
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 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 > ¶ms)
Constructor for a fused import.
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.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const
Returns the Map that describes the column distribution in this 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.
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.
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.
size_t getNodeNumRows() const
Returns the number of matrix rows owned on the calling node.
size_t getNodeNumRows() const
Returns the number of matrix rows owned on the calling node.
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.
TpetraCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraCrsMatrixClass
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 > ¶ms=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.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying (possibly different) number of entries in each row.
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 > ¶ms=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.
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...
size_t getNodeMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on this node.
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.
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.
TpetraExport< LocalOrdinal, GlobalOrdinal, Node > TpetraExportClass
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.
size_t getNodeNumCols() const
Returns the number of columns connected to the locally owned rows of this matrix. ...