42 #ifndef TPETRA_BLOCKCRSMATRIX_DECL_HPP
43 #define TPETRA_BLOCKCRSMATRIX_DECL_HPP
48 #include "Tpetra_CrsGraph.hpp"
49 #include "Tpetra_RowMatrix.hpp"
50 #include "Tpetra_BlockMultiVector_decl.hpp"
132 template<
class Scalar,
143 using STS = Teuchos::ScalarTraits<Scalar>;
184 typedef ::Tpetra::MultiVector<Scalar, LO, GO, node_type>
mv_type;
192 Kokkos::MemoryTraits<Kokkos::Unmanaged> >
198 Kokkos::MemoryTraits<Kokkos::Unmanaged> >
204 Kokkos::MemoryTraits<Kokkos::Unmanaged> >
210 Kokkos::MemoryTraits<Kokkos::Unmanaged> >
257 Teuchos::RCP<const map_type>
getRowMap ()
const;
260 Teuchos::RCP<const map_type>
getColMap ()
const;
282 Teuchos::ETransp mode = Teuchos::NO_TRANS,
283 Scalar alpha = Teuchos::ScalarTraits<Scalar>::one (),
284 Scalar beta = Teuchos::ScalarTraits<Scalar>::zero ())
const;
328 describe (Teuchos::FancyOStream& out,
329 const Teuchos::EVerbosityLevel verbLevel)
const;
339 virtual Teuchos::RCP<const ::Tpetra::RowGraph<LO,GO,Node> >
getGraph ()
const;
348 applyBlock (
const BlockMultiVector<Scalar, LO, GO, Node>& X,
349 BlockMultiVector<Scalar, LO, GO, Node>& Y,
350 Teuchos::ETransp mode = Teuchos::NO_TRANS,
351 const Scalar alpha = Teuchos::ScalarTraits<Scalar>::one (),
352 const Scalar beta = Teuchos::ScalarTraits<Scalar>::zero ());
359 const ::Tpetra::MultiVector<Scalar,LO,GO,Node> &B,
360 const ::Tpetra::MultiVector<Scalar,LO,GO,Node> &D,
361 const Scalar& dampingFactor,
364 const bool zeroInitialGuess)
const;
371 const ::Tpetra::MultiVector<Scalar,LO,GO,Node>& B,
372 const ::Tpetra::MultiVector<Scalar,LO,GO,Node>& D,
373 const Teuchos::ArrayView<LO>& rowIndices,
374 const Scalar& dampingFactor,
377 const bool zeroInitialGuess)
const;
400 BlockMultiVector<Scalar, LO, GO, Node>& Solution,
402 Kokkos::MemoryUnmanaged>& D_inv,
436 const LO numColInds)
const;
468 const LO numColInds)
const;
508 Teuchos::ArrayView<const LO> &indices,
509 Teuchos::ArrayView<const Scalar> &values)
const;
514 const Teuchos::ArrayView<LO> &Indices,
515 const Teuchos::ArrayView<Scalar> &Values,
516 size_t &NumEntries)
const;
519 getLocalBlock (
const LO localRowInd,
const LO localColInd)
const;
548 const LO numColInds)
const;
557 const ptrdiff_t offsets[],
559 const LO numOffsets)
const;
568 const ptrdiff_t offsets[],
570 const LO numOffsets)
const;
615 return (*errs_).is_null () ? std::string (
"") : (*errs_)->str ();
651 Kokkos::MemoryUnmanaged>& offsets)
const;
653 #ifdef TPETRA_ENABLE_DEPRECATED_CODE
654 void TPETRA_DEPRECATED
661 #endif // TPETRA_ENABLE_DEPRECATED_CODE
678 Kokkos::MemoryUnmanaged>& diag,
680 Kokkos::MemoryUnmanaged>& offsets)
const;
697 Kokkos::MemoryUnmanaged>& diag,
698 const Teuchos::ArrayView<const size_t>& offsets)
const;
700 #ifdef TPETRA_ENABLE_DEPRECATED_CODE
701 void TPETRA_DEPRECATED
716 const Teuchos::ArrayView<const size_t>& offsets)
const;
717 #endif // TPETRA_ENABLE_DEPRECATED_CODE
725 const LO numColInds)
const;
730 const ptrdiff_t offsets[],
732 const LO numOffsets)
const;
748 virtual bool checkSizes (const ::Tpetra::SrcDistObject& source);
751 #ifdef TPETRA_ENABLE_DEPRECATED_CODE
753 #else // TPETRA_ENABLE_DEPRECATED_CODE
755 #endif // TPETRA_ENABLE_DEPRECATED_CODE
757 const size_t numSameIDs,
764 #ifdef TPETRA_ENABLE_DEPRECATED_CODE
766 #else // TPETRA_ENABLE_DEPRECATED_CODE
768 #endif // TPETRA_ENABLE_DEPRECATED_CODE
774 Kokkos::DualView<
size_t*,
776 size_t& constantNumPackets,
780 #ifdef TPETRA_ENABLE_DEPRECATED_CODE
782 #else // TPETRA_ENABLE_DEPRECATED_CODE
784 #endif // TPETRA_ENABLE_DEPRECATED_CODE
789 Kokkos::DualView<
size_t*,
791 const size_t constantNumPackets,
799 Teuchos::RCP<crs_graph_type> graphRCP_;
839 typename crs_graph_type::local_graph_type::row_map_type::HostMirror ptrHost_;
846 typename crs_graph_type::local_graph_type::entries_type::HostMirror indHost_;
853 typename Kokkos::DualView<impl_scalar_type*, device_type> val_;
876 Teuchos::RCP<Teuchos::RCP<BMV> > X_colMap_;
880 Teuchos::RCP<Teuchos::RCP<BMV> > Y_rowMap_;
889 Teuchos::RCP<Teuchos::RCP<typename crs_graph_type::import_type> > pointImporter_;
905 Teuchos::RCP<bool> localError_;
914 Teuchos::RCP<Teuchos::RCP<std::ostringstream> > errs_;
917 std::ostream& markLocalErrorAndGetStream ();
922 template<
class Device>
924 #if defined(KOKKOS_ENABLE_CUDA)
930 static constexpr
bool value =
931 std::is_same<typename Device::execution_space, Kokkos::Cuda>::value;
933 static constexpr
bool value =
false;
934 #endif // defined(KOKKOS_ENABLE_CUDA)
950 val_.modify_device();
954 template<
class MemorySpace>
957 if (is_cuda<MemorySpace>::value) {
968 return val_.need_sync_host();
974 return val_.need_sync_device();
978 template<
class MemorySpace>
981 if (is_cuda<MemorySpace>::value) {
1002 template<
class MemorySpace>
1005 if (is_cuda<MemorySpace>::value) {
1014 typename Kokkos::DualView<impl_scalar_type*, device_type>::t_host getValuesHost ()
const {
1015 return val_.view_host();
1019 typename Kokkos::DualView<impl_scalar_type*, device_type>::t_dev getValuesDevice ()
const {
1020 return val_.view_device();
1040 template<
class MemorySpace>
1041 typename std::conditional<is_cuda<MemorySpace>::value,
1042 typename Kokkos::DualView<impl_scalar_type*, device_type>::t_dev,
1043 typename Kokkos::DualView<impl_scalar_type*, device_type>::t_host>::type
1047 return Kokkos::Impl::if_c<
1048 is_cuda<MemorySpace>::value,
1049 typename Kokkos::DualView<impl_scalar_type*, device_type>::t_dev,
1050 typename Kokkos::DualView<impl_scalar_type*, device_type>::t_host
1051 >::select (this->getValuesDevice (), this->getValuesHost ());
1070 const Teuchos::ETransp mode,
1140 findRelOffsetOfColumnIndex (
const LO localRowIndex,
1141 const LO colIndexToFind,
1142 const LO hint = 0)
const;
1146 LO offsetPerBlock ()
const;
1149 getConstLocalBlockFromInput (
const impl_scalar_type* val,
const size_t pointOffset)
const;
1152 getNonConstLocalBlockFromInput (
impl_scalar_type* val,
const size_t pointOffset)
const;
1155 getConstLocalBlockFromAbsOffset (
const size_t absBlockOffset)
const;
1158 getNonConstLocalBlockFromAbsOffset (
const size_t absBlockOffset)
const;
1164 getConstLocalBlockFromRelOffset (
const LO lclMeshRow,
1165 const size_t relMeshOffset)
const;
1169 virtual Teuchos::RCP<const Teuchos::Comm<int> >
getComm()
const;
1171 #ifdef TPETRA_ENABLE_DEPRECATED_CODE
1172 virtual TPETRA_DEPRECATED Teuchos::RCP<Node> getNode()
const;
1174 #endif // TPETRA_ENABLE_DEPRECATED_CODE
1235 #ifdef TPETRA_ENABLE_DEPRECATED_CODE
1252 virtual size_t TPETRA_DEPRECATED getNodeNumDiags()
const;
1264 virtual bool TPETRA_DEPRECATED isLowerTriangular ()
const;
1276 virtual bool TPETRA_DEPRECATED isUpperTriangular ()
const;
1277 #endif // TPETRA_ENABLE_DEPRECATED_CODE
1305 const Teuchos::ArrayView<GO> &Indices,
1306 const Teuchos::ArrayView<Scalar> &Values,
1307 size_t& NumEntries)
const;
1335 Teuchos::ArrayView<const GO>& indices,
1336 Teuchos::ArrayView<const Scalar>& values)
const;
1360 virtual void leftScale (const ::Tpetra::Vector<Scalar, LO, GO, Node>& x);
1367 virtual void rightScale (const ::Tpetra::Vector<Scalar, LO, GO, Node>& x);
1377 virtual typename ::Tpetra::RowMatrix<Scalar, LO, GO, Node>::mag_type
1384 #endif // TPETRA_BLOCKCRSMATRIX_DECL_HPP
virtual size_t getNumEntriesInGlobalRow(GO globalRow) const
The current number of entries on the calling process in the specified global row. ...
::Tpetra::MultiVector< Scalar, LO, GO, node_type > mv_type
The implementation of MultiVector that this class uses.
typename BMV::impl_scalar_type impl_scalar_type
The implementation type of entries in the matrix.
void getLocalDiagCopy(const Kokkos::View< impl_scalar_type ***, device_type, Kokkos::MemoryUnmanaged > &diag, const Kokkos::View< const size_t *, device_type, Kokkos::MemoryUnmanaged > &offsets) const
Variant of getLocalDiagCopy() that uses precomputed offsets and puts diagonal blocks in a 3-D Kokkos:...
LO replaceLocalValues(const LO localRowInd, const LO colInds[], const Scalar vals[], const LO numColInds) const
Replace values at the given (mesh, i.e., block) column indices, in the given (mesh, i.e., block) row.
void gaussSeidelCopy(MultiVector< Scalar, LO, GO, Node > &X, const ::Tpetra::MultiVector< Scalar, LO, GO, Node > &B, const ::Tpetra::MultiVector< Scalar, LO, GO, Node > &D, const Scalar &dampingFactor, const ESweepDirection direction, const int numSweeps, const bool zeroInitialGuess) const
Version of gaussSeidel(), with fewer requirements on X.
LO replaceLocalValuesByOffsets(const LO localRowInd, const ptrdiff_t offsets[], const Scalar vals[], const LO numOffsets) const
Like replaceLocalValues, but avoids computing row offsets.
Sparse matrix whose entries are small dense square blocks, all of the same dimensions.
size_t getNodeNumRows() const
get the local number of block rows
void modify_device()
Mark the matrix's valueas as modified in device space.
typename DistObject< Scalar, LO, GO, Node >::buffer_device_type buffer_device_type
Kokkos::Device specialization for communication buffers.
virtual bool isGloballyIndexed() const
Whether matrix indices are globally indexed.
size_t getNumEntriesInLocalRow(const LO localRowInd) const
Return the number of entries in the given row on the calling process.
LO local_ordinal_type
The type of local indices.
virtual void getGlobalRowView(GO GlobalRow, Teuchos::ArrayView< const GO > &indices, Teuchos::ArrayView< const Scalar > &values) const
Get a constant, nonpersisting, globally indexed view of the given row of the matrix.
virtual typename::Tpetra::RowMatrix< Scalar, LO, GO, Node >::mag_type getFrobeniusNorm() const
The Frobenius norm of the matrix.
virtual global_size_t getGlobalNumCols() const
The global number of columns of this matrix.
::Tpetra::CrsGraph< LO, GO, node_type > crs_graph_type
The implementation of CrsGraph that this class uses.
void setAllToScalar(const Scalar &alpha)
Set all matrix entries equal to alpha.
Node node_type
The Node type.
virtual bool isLocallyIndexed() const
Whether matrix indices are locally indexed.
bool localError() const
Whether this object had an error on the calling process.
Declaration of the Tpetra::CrsMatrix class.
void sync_device()
Sync the matrix's values to device space.
Scalar scalar_type
The type of entries in the matrix (that is, of each entry in each block).
virtual void leftScale(const ::Tpetra::Vector< Scalar, LO, GO, Node > &x)
Scale the RowMatrix on the left with the given Vector x.
virtual bool hasColMap() const
Whether this matrix has a well-defined column Map.
Kokkos::View< impl_scalar_type *, Kokkos::LayoutRight, device_type, Kokkos::MemoryTraits< Kokkos::Unmanaged > > little_vec_type
The type used to access nonconst vector blocks.
virtual size_t getNodeNumEntries() const
The local number of stored (structurally nonzero) entries.
::Tpetra::Map< LO, GO, node_type > map_type
The implementation of Map that this class uses.
virtual size_t getNodeNumCols() const
The number of columns needed to apply the forward operator on this node.
virtual ~BlockCrsMatrix()
Destructor (declared virtual for memory safety).
MultiVector for multiple degrees of freedom per mesh point.
GO global_ordinal_type
The type of global indices.
size_t global_size_t
Global size_t object.
BlockCrsMatrix()
Default constructor: Makes an empty block matrix.
virtual bool isFillComplete() const
Whether fillComplete() has been called.
void apply(const mv_type &X, mv_type &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=Teuchos::ScalarTraits< Scalar >::one(), Scalar beta=Teuchos::ScalarTraits< Scalar >::zero()) const
For this matrix A, compute Y := beta * Y + alpha * Op(A) * X.
LO absMaxLocalValues(const LO localRowInd, const LO colInds[], const Scalar vals[], const LO numColInds) const
Like sumIntoLocalValues, but for the ABSMAX combine mode.
virtual size_t getGlobalMaxNumRowEntries() const
The maximum number of entries in any row over all processes in the matrix's communicator.
device_type::memory_space memory_space
The Kokkos memory space that this class uses.
size_t getNodeMaxNumRowEntries() const
Maximum number of entries in any row of the matrix, on this process.
std::string errorMessages() const
The current stream of error messages.
LO sumIntoLocalValuesByOffsets(const LO localRowInd, const ptrdiff_t offsets[], const Scalar vals[], const LO numOffsets) const
Like sumIntoLocalValues, but avoids computing row offsets.
virtual void rightScale(const ::Tpetra::Vector< Scalar, LO, GO, Node > &x)
Scale the RowMatrix on the right with the given Vector x.
ESweepDirection
Sweep direction for Gauss-Seidel or Successive Over-Relaxation (SOR).
Teuchos::RCP< const map_type > getRangeMap() const
Get the (point) range Map of this matrix.
bool need_sync_device() const
Whether the matrix's values need sync'ing to device space.
Teuchos::RCP< const map_type > getRowMap() const
get the (mesh) map for the rows of this block matrix.
std::conditional< is_cuda< MemorySpace >::value, typename Kokkos::DualView< impl_scalar_type *, device_type >::t_dev, typename Kokkos::DualView< impl_scalar_type *, device_type >::t_host >::type getValues()
Get the host or device View of the matrix's values (val_).
void modify_host()
Mark the matrix's valueas as modified in host space.
virtual global_size_t getGlobalNumEntries() const
The global number of stored (structurally nonzero) entries.
Sets up and executes a communication plan for a Tpetra DistObject.
void applyBlock(const BlockMultiVector< Scalar, LO, GO, Node > &X, BlockMultiVector< Scalar, LO, GO, Node > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, const Scalar alpha=Teuchos::ScalarTraits< Scalar >::one(), const Scalar beta=Teuchos::ScalarTraits< Scalar >::zero())
Version of apply() that takes BlockMultiVector input and output.
CombineMode
Rule for combining data in an Import or Export.
global_size_t getGlobalNumRows() const
get the global number of block rows
CrsGraphType::global_ordinal_type getGlobalNumDiags(const CrsGraphType &G)
Number of populated diagonal entries in the given sparse graph, over all processes in the graph's (MP...
virtual void getGlobalRowCopy(GO GlobalRow, const Teuchos::ArrayView< GO > &Indices, const Teuchos::ArrayView< Scalar > &Values, size_t &NumEntries) const
Get a copy of the given global row's entries.
LO getLocalRowOffsets(const LO localRowInd, ptrdiff_t offsets[], const LO colInds[], const LO numColInds) const
Get relative offsets corresponding to the given rows, given by local row index.
void sync_host()
Sync the matrix's values to host space.
Abstract base class for objects that can be the source of an Import or Export operation.
device_type::execution_space execution_space
The Kokkos execution space that this class uses.
Kokkos::View< const impl_scalar_type *, Kokkos::LayoutRight, device_type, Kokkos::MemoryTraits< Kokkos::Unmanaged > > const_little_vec_type
The type used to access const vector blocks.
Kokkos::View< impl_scalar_type **, Kokkos::LayoutRight, device_type, Kokkos::MemoryTraits< Kokkos::Unmanaged > > little_block_type
The type used to access nonconst matrix blocks.
virtual bool supportsRowViews() const
Whether this object implements getLocalRowView() and getGlobalRowView().
LO absMaxLocalValuesByOffsets(const LO localRowInd, const ptrdiff_t offsets[], const Scalar vals[], const LO numOffsets) const
Like sumIntoLocalValuesByOffsets, but for the ABSMAX combine mode.
virtual Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
The communicator over which this matrix is distributed.
LO getLocalRowView(const LO localRowInd, const LO *&colInds, Scalar *&vals, LO &numInds) const
Get a view of the (mesh, i.e., block) row, using local (mesh, i.e., block) indices.
void getLocalDiagOffsets(const Kokkos::View< size_t *, device_type, Kokkos::MemoryUnmanaged > &offsets) const
Get offsets of the diagonal entries in the matrix.
char packet_type
Implementation detail; tells.
void sync()
Sync the matrix's values to the given memory space.
void localGaussSeidel(const BlockMultiVector< Scalar, LO, GO, Node > &Residual, BlockMultiVector< Scalar, LO, GO, Node > &Solution, const Kokkos::View< impl_scalar_type ***, device_type, Kokkos::MemoryUnmanaged > &D_inv, const Scalar &omega, const ESweepDirection direction) const
Local Gauss-Seidel solve, given a factorized diagonal.
A read-only, row-oriented interface to a sparse matrix.
void modify()
Mark the matrix's values as modified in the given memory space.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
Print a description of this object to the given output stream.
A distributed dense vector.
LO getBlockSize() const
The number of degrees of freedom per mesh point.
bool need_sync() const
Whether the matrix's values need sync'ing to the given memory space.
Teuchos::RCP< const map_type > getDomainMap() const
Get the (point) domain Map of this matrix.
bool hasTransposeApply() const
Whether it is valid to apply the transpose or conjugate transpose of this matrix. ...
void getLocalRowCopy(LO LocalRow, const Teuchos::ArrayView< LO > &Indices, const Teuchos::ArrayView< Scalar > &Values, size_t &NumEntries) const
Not implemented.
Kokkos::View< const impl_scalar_type **, Kokkos::LayoutRight, device_type, Kokkos::MemoryTraits< Kokkos::Unmanaged > > const_little_block_type
The type used to access const matrix blocks.
Node::device_type device_type
The Kokkos::Device specialization that this class uses.
Teuchos::RCP< const map_type > getColMap() const
get the (mesh) map for the columns of this block matrix.
virtual Teuchos::RCP< const ::Tpetra::RowGraph< LO, GO, Node > > getGraph() const
Get the (mesh) graph.
Base class for distributed Tpetra objects that support data redistribution.
typename mv_type::impl_scalar_type impl_scalar_type
The implementation type of entries in the object.
void reorderedGaussSeidelCopy(::Tpetra::MultiVector< Scalar, LO, GO, Node > &X, const ::Tpetra::MultiVector< Scalar, LO, GO, Node > &B, const ::Tpetra::MultiVector< Scalar, LO, GO, Node > &D, const Teuchos::ArrayView< LO > &rowIndices, const Scalar &dampingFactor, const ESweepDirection direction, const int numSweeps, const bool zeroInitialGuess) const
Version of reorderedGaussSeidel(), with fewer requirements on X.
LO sumIntoLocalValues(const LO localRowInd, const LO colInds[], const Scalar vals[], const LO numColInds) const
Sum into values at the given (mesh, i.e., block) column indices, in the given (mesh, i.e., block) row.
std::string description() const
One-line description of this object.
bool need_sync_host() const
Whether the matrix's values need sync'ing to host space.
virtual GO getIndexBase() const
The index base for global indices in this matrix.