49 #ifndef XPETRA_MATRIX_HPP 50 #define XPETRA_MATRIX_HPP 66 #include <Teuchos_SerialDenseMatrix.hpp> 67 #include <Teuchos_Hashtable.hpp> 93 template <class Scalar = Operator<>::scalar_type,
97 class Matrix :
public Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > {
110 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR 111 #ifdef HAVE_XPETRA_TPETRA 112 typedef typename CrsMatrix::local_matrix_type local_matrix_type;
130 RCP<MatrixView> view = rcp(
new MatrixView(rowMap, colMap));
135 void CreateView(
const viewLabel_t viewLabel,
const RCP<const Matrix>& A,
bool transposeA =
false,
const RCP<const Matrix>& B = Teuchos::null,
bool transposeB =
false) {
136 RCP<const Map> domainMap = Teuchos::null;
137 RCP<const Map> rangeMap = Teuchos::null;
139 const size_t blkSize = 1;
140 std::vector<size_t> stridingInfo(1, blkSize);
141 LocalOrdinal stridedBlockId = -1;
144 if (A->IsView(viewLabel)) {
145 rangeMap = transposeA ? A->getColMap(viewLabel) : A->getRowMap(viewLabel);
146 domainMap = transposeA ? A->getRowMap(viewLabel) : A->getColMap(viewLabel);
149 rangeMap = transposeA ? A->getDomainMap() : A->getRangeMap();
150 domainMap = transposeA ? A->getRangeMap() : A->getDomainMap();
152 if (viewLabel ==
"stridedMaps") {
158 if (B != Teuchos::null ) {
161 if (B->IsView(viewLabel)) {
162 domainMap = transposeB ? B->getRowMap(viewLabel) : B->getColMap(viewLabel);
165 domainMap = transposeB ? B->getRangeMap() : B->getDomainMap();
167 if (viewLabel ==
"stridedMaps")
181 int last = out.getOutputToRootOnly();
182 Teuchos::OSTab tab(out);
183 out.setOutputToRootOnly(0);
184 Teuchos::Array<viewLabel_t> viewLabels;
185 Teuchos::Array<RCP<MatrixView> > viewList;
187 out <<
"views associated with this operator" << std::endl;
188 for (
int i=0; i<viewLabels.size(); ++i)
189 out << viewLabels[i] << std::endl;
190 out.setOutputToRootOnly(last);
234 virtual void insertGlobalValues(GlobalOrdinal globalRow,
const ArrayView<const GlobalOrdinal> &cols,
const ArrayView<const Scalar> &vals) = 0;
244 virtual void insertLocalValues(LocalOrdinal localRow,
const ArrayView<const LocalOrdinal> &cols,
const ArrayView<const Scalar> &vals) = 0;
253 const ArrayView<const GlobalOrdinal> &cols,
254 const ArrayView<const Scalar> &vals) = 0;
261 const ArrayView<const LocalOrdinal> &cols,
262 const ArrayView<const Scalar> &vals) = 0;
268 virtual void scale(
const Scalar &alpha)= 0;
283 virtual void resumeFill(
const RCP< ParameterList > ¶ms=null) = 0;
296 virtual void fillComplete(
const RCP<const Map> &domainMap,
const RCP<const Map> &rangeMap,
const RCP<ParameterList> ¶ms = null) =0;
312 virtual void fillComplete(
const RCP<ParameterList> ¶ms=null) =0;
404 const ArrayView<LocalOrdinal> &Indices,
405 const ArrayView<Scalar> &Values,
419 virtual void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView<const GlobalOrdinal> &indices, ArrayView<const Scalar> &values)
const =0;
431 virtual void getLocalRowView(LocalOrdinal LocalRow, ArrayView<const LocalOrdinal> &indices, ArrayView<const Scalar> &values)
const =0;
439 virtual typename ScalarTraits<Scalar>::magnitudeType
getFrobeniusNorm()
const = 0;
471 virtual const Teuchos::RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > >
getMap()
const = 0;
514 virtual void describe(Teuchos::FancyOStream &out,
const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default)
const =0;
520 virtual RCP<const CrsGraph>
getCrsGraph()
const =0;
534 std::vector<size_t> stridingInfo;
535 stridingInfo.push_back(Teuchos::as<size_t>(blksize));
536 LocalOrdinal stridedBlockId = -1;
552 CreateView(
"stridedMaps", stridedRangeMap, stridedDomainMap);
558 if(
IsView(
"stridedMaps")==
true) {
561 TEUCHOS_TEST_FOR_EXCEPTION(rangeMap == Teuchos::null,
Exceptions::BadCast,
"Xpetra::Matrix::GetFixedBlockSize(): rangeMap is not of type StridedMap");
562 TEUCHOS_TEST_FOR_EXCEPTION(domainMap == Teuchos::null,
Exceptions::BadCast,
"Xpetra::Matrix::GetFixedBlockSize(): domainMap is not of type StridedMap");
563 TEUCHOS_TEST_FOR_EXCEPTION(domainMap->getFixedBlockSize() != rangeMap->getFixedBlockSize(),
Exceptions::RuntimeError,
"Xpetra::Matrix::GetFixedBlockSize(): block size of rangeMap and domainMap are different.");
564 return Teuchos::as<LocalOrdinal>(domainMap->getFixedBlockSize());
583 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR 584 #ifdef HAVE_XPETRA_TPETRA 585 virtual local_matrix_type getLocalMatrix ()
const = 0;
589 #warning "Xpetra Kokkos interface for CrsMatrix is enabled (HAVE_XPETRA_KOKKOS_REFACTOR) but Tpetra is disabled. The Kokkos interface needs Tpetra to be enabled, too." 605 #define XPETRA_MATRIX_SHORT 606 #endif //XPETRA_MATRIX_DECL_HPP virtual void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)=0
Insert matrix entries, using local IDs.
virtual size_t getNodeMaxNumRowEntries() const =0
Returns the maximum number of entries across all rows/columns on this node.
virtual size_t getNodeNumRows() const =0
Returns the number of matrix rows owned on the calling node.
virtual size_t getNodeNumEntries() const =0
Returns the local number of entries in this matrix.
virtual global_size_t getGlobalNumRows() const =0
Returns the number of global rows in this matrix.
virtual const RCP< const Map > & getRowMap() const
Returns the Map that describes the row distribution in this matrix.
virtual void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView< const GlobalOrdinal > &indices, ArrayView< const Scalar > &values) const =0
Extract a const, non-persisting view of global indices in a specified row of the matrix.
void RemoveView(const viewLabel_t viewLabel)
virtual global_size_t getGlobalNumCols() const =0
Returns the number of global columns in the matrix.
LocalOrdinal local_ordinal_type
virtual ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const =0
Get Frobenius norm of the matrix.
virtual void resumeFill(const RCP< ParameterList > ¶ms=null)=0
void PrintViews(Teuchos::FancyOStream &out) const
Print all of the views associated with the Matrix.
Exception throws to report errors in the internal logical of the program.
virtual size_t getGlobalMaxNumRowEntries() const =0
Returns the maximum number of entries across all rows/columns on all nodes.
Teuchos::Hashtable< viewLabel_t, RCP< MatrixView > > operatorViewTable_
virtual bool isGloballyIndexed() const =0
If matrix indices are in the global range, this function returns true. Otherwise, this function retur...
virtual void doExport(const Matrix &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)=0
Export.
Node node_type
The Kokkos Node type.
virtual const RCP< const Map > & getColMap() const
Returns the Map that describes the column distribution in this matrix. This might be null until fillC...
Exception indicating invalid cast attempted.
GlobalOrdinal global_ordinal_type
The global index type.
virtual const RCP< const Map > & getColMap(viewLabel_t viewLabel) const
Returns the Map that describes the column distribution in this matrix.
virtual const Teuchos::RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const =0
Multiplies this matrix by a MultiVector.
virtual void SetMaxEigenvalueEstimate(Scalar const &sigma)
virtual void fillComplete(const RCP< const Map > &domainMap, const RCP< const Map > &rangeMap, const RCP< ParameterList > ¶ms=null)=0
Signal that data entry is complete, specifying domain and range maps.
Xpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > CrsMatrix
virtual ~Matrix()
Destructor.
virtual void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag) const =0
Get a copy of the diagonal entries owned by this node, with local row idices.
virtual void getLocalRowView(LocalOrdinal LocalRow, ArrayView< const LocalOrdinal > &indices, ArrayView< const Scalar > &values) const =0
Extract a const, non-persisting view of local indices in a specified row of the matrix.
virtual void getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView< LocalOrdinal > &Indices, const ArrayView< Scalar > &Values, size_t &NumEntries) const =0
Extract a list of entries in a specified local row of the matrix. Put into storage allocated by calli...
viewLabel_t defaultViewLabel_
virtual void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)=0
Insert matrix entries, using global IDs.
void SetFixedBlockSize(LocalOrdinal blksize, GlobalOrdinal offset=0)
virtual global_size_t getGlobalNumDiags() const =0
Returns the number of global diagonal entries, based on global row/column index comparisons.
virtual void rightScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)=0
Right scale matrix using the given vector entries.
LocalOrdinal local_ordinal_type
The local index type.
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const =0
Print the object with some verbosity level to an FancyOStream object.
Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > Map
virtual bool isFillComplete() const =0
Returns true if fillComplete() has been called and the matrix is in compute mode. ...
const viewLabel_t & GetDefaultViewLabel() const
static RCP< StridedMap > Build(UnderlyingLib lib, global_size_t numGlobalElements, GlobalOrdinal indexBase, std::vector< size_t > &stridingInfo, const Teuchos::RCP< const Teuchos::Comm< int > > &comm, LocalOrdinal stridedBlockId=-1, GlobalOrdinal offset=0, LocalGlobal lg=Xpetra::GloballyDistributed, const Teuchos::RCP< Node > &node=defaultArgNode())
Map constructor with Xpetra-defined contiguous uniform distribution.
void CreateView(viewLabel_t viewLabel, const RCP< const Map > &rowMap, const RCP< const Map > &colMap)
virtual void replaceLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)=0
Replace matrix entries, using local IDs.
viewLabel_t currentViewLabel_
const viewLabel_t SwitchToDefaultView()
virtual Scalar GetMaxEigenvalueEstimate() const
virtual void leftScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)=0
Left scale matrix using the given vector entries.
virtual size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const =0
Returns the current number of entries on this node in the specified local row.
Xpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > CrsGraph
size_t global_size_t
Global size_t object.
LocalOrdinal GetFixedBlockSize() const
virtual void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)=0
Replace matrix entries, using global IDs.
virtual void doImport(const Matrix &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)=0
Import.
virtual RCP< const CrsGraph > getCrsGraph() const =0
Returns the CrsGraph associated with this matrix.
virtual size_t getNodeNumDiags() const =0
Returns the number of local diagonal entries, based on global row/column index comparisons.
virtual std::string description() const =0
Return a simple one-line description of this object.
virtual bool isLocallyIndexed() const =0
If matrix indices are in the local range, this function returns true. Otherwise, this function return...
const viewLabel_t & GetCurrentViewLabel() const
virtual global_size_t getGlobalNumEntries() const =0
Returns the global number of entries in this matrix.
virtual Teuchos::RCP< const Map > getRangeMap() const =0
The Map associated with the range of this operator, which must be compatible with Y...
bool IsView(const viewLabel_t viewLabel) const
CombineMode
Xpetra::Combine Mode enumerable type.
virtual const RCP< const Map > & getRowMap(viewLabel_t viewLabel) const
Returns the Map that describes the row distribution in this matrix.
void CreateView(const viewLabel_t viewLabel, const RCP< const Matrix > &A, bool transposeA=false, const RCP< const Matrix > &B=Teuchos::null, bool transposeB=false)
virtual void scale(const Scalar &alpha)=0
Scale the current values of a matrix, this = alpha*this.
Xpetra::MatrixView< Scalar, LocalOrdinal, GlobalOrdinal, Node > MatrixView
virtual void setAllToScalar(const Scalar &alpha)=0
Set all matrix entries equal to scalar.
virtual Teuchos::RCP< const Map > getDomainMap() const =0
The Map associated with the domain of this operator, which must be compatible with X...
Xpetra::CrsMatrixFactory< Scalar, LocalOrdinal, GlobalOrdinal, Node > CrsMatrixFactory
Xpetra-specific matrix class.
Class that stores a strided map.
const viewLabel_t SwitchToView(const viewLabel_t viewLabel)
GlobalOrdinal global_ordinal_type