Xpetra Version of the Day
Loading...
Searching...
No Matches
Xpetra_TpetraBlockCrsMatrix_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_TPETRABLOCKCRSMATRIX_DECL_HPP
11#define XPETRA_TPETRABLOCKCRSMATRIX_DECL_HPP
12
13/* this file is automatically generated - do not edit (see script/tpetra.py) */
14
16
17#include "Tpetra_BlockCrsMatrix.hpp"
18#include "Tpetra_CrsMatrix.hpp"
19
20#include "Xpetra_CrsMatrix.hpp"
25#include "Xpetra_Exceptions.hpp"
26
27namespace Xpetra {
28
29template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = Tpetra::KokkosClassic::DefaultNode::DefaultNodeType>
31 : public CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> //, public TpetraRowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>
32{
33 // The following typedef are used by the XPETRA_DYNAMIC_CAST() macro.
38
39 public:
41
44 size_t maxNumEntriesPerRow,
45 const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null);
46
49 const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc,
50 const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null);
51
55 size_t maxNumEntriesPerRow,
56 const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null);
57
61 const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc,
62 const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null);
63
66 const Teuchos::RCP<Teuchos::ParameterList> &params = Teuchos::null);
67
70 const LocalOrdinal blockSize);
71
74 const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &pointDomainMap,
75 const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &pointRangeMap,
76 const LocalOrdinal blockSize);
77
79 TpetraBlockCrsMatrix(const Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &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 TpetraBlockCrsMatrix(const Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &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 TpetraBlockCrsMatrix(const Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &sourceMatrix,
99
101 TpetraBlockCrsMatrix(const Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &sourceMatrix,
103 const Teuchos::RCP<const Export<LocalOrdinal, GlobalOrdinal, Node> > DomainExporter,
107
109 virtual ~TpetraBlockCrsMatrix();
110
112
114 void insertGlobalValues(GlobalOrdinal globalRow,
116 const ArrayView<const Scalar> &vals);
117
119 void insertLocalValues(LocalOrdinal localRow,
121 const ArrayView<const Scalar> &vals);
122
124 void replaceGlobalValues(GlobalOrdinal globalRow,
126 const ArrayView<const Scalar> &vals);
127
129 void replaceLocalValues(LocalOrdinal localRow,
131 const ArrayView<const Scalar> &vals);
132
134 void setAllToScalar(const Scalar &alpha);
135
137 void scale(const Scalar &alpha);
138
140 //** \warning This is an expert-only routine and should not be called from user code. (not implemented)
141 void allocateAllValues(size_t numNonZeros, ArrayRCP<size_t> &rowptr,
143 ArrayRCP<Scalar> &values);
144
146 void setAllValues(const ArrayRCP<size_t> &rowptr,
147 const ArrayRCP<LocalOrdinal> &colind,
148 const ArrayRCP<Scalar> &values);
149
153 ArrayRCP<const Scalar> &values) const;
154
156 void getAllValues(ArrayRCP<Scalar> &values);
157
159
161 void resumeFill(const RCP<ParameterList> &params = null);
162
164 void fillComplete(const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &domainMap,
165 const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap, const RCP<ParameterList> &params = null);
166
168 void fillComplete(const RCP<ParameterList> &params = null);
169
173
176 const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > &rangeMap,
177 const RCP<const Import<LocalOrdinal, GlobalOrdinal, Node> > &importer = Teuchos::null,
178 const RCP<const Export<LocalOrdinal, GlobalOrdinal, Node> > &exporter = Teuchos::null,
179 const RCP<ParameterList> &params = Teuchos::null);
180
182
185
188
191
194
197
199 size_t getLocalNumRows() const;
200
202 size_t getLocalNumCols() const;
203
206
208 size_t getLocalNumEntries() const;
209
211 size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const;
212
214 size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const;
215
217 size_t getGlobalMaxNumRowEntries() const;
218
220 size_t getLocalMaxNumRowEntries() const;
221
223 bool isLocallyIndexed() const;
224
226 bool isGloballyIndexed() const;
227
229 bool isFillComplete() const;
230
232 bool isFillActive() const;
233
236
238 bool supportsRowViews() const;
239
241 void getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView<LocalOrdinal> &Indices, const ArrayView<Scalar> &Values, size_t &NumEntries) const;
242
244 void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView<const GlobalOrdinal> &indices, ArrayView<const Scalar> &values) const;
245
247 void getGlobalRowCopy(GlobalOrdinal GlobalRow, const ArrayView<GlobalOrdinal> &indices, const ArrayView<Scalar> &values, size_t &numEntries) const;
248
250 void getLocalRowView(LocalOrdinal LocalRow, ArrayView<const LocalOrdinal> &indices, ArrayView<const Scalar> &values) const;
251
253 virtual bool haveGlobalConstants() const;
254
256
259
261 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<Import<LocalOrdinal, GlobalOrdinal, Node> > &regionInterfaceImporter, const Teuchos::ArrayRCP<LocalOrdinal> &regionInterfaceLIDs) const;
262
265
268
270
272 std::string description() const;
273
276
278 void setObjectLabel(const std::string &objectLabel);
279
282
285 const Teuchos::ArrayView<const size_t> &offsets) const;
286
289
291 void getLocalDiagCopy(Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node> &diag, const Kokkos::View<const size_t *, typename Node::device_type, Kokkos::MemoryUnmanaged> &offsets) const;
292
294
297
300
302
305
309
313
317
321
323
325
327 bool hasMatrix() const;
328
330 TpetraBlockCrsMatrix(const Teuchos::RCP<Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &mtx);
331
334
337
338#ifdef HAVE_XPETRA_TPETRA
339 // using local_matrix_type = typename Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::local_matrix_type;
341
343#if KOKKOS_VERSION >= 40799
344 typename local_matrix_type::host_mirror_type getLocalMatrixHost() const;
345#else
346 typename local_matrix_type::HostMirror getLocalMatrixHost() const;
347#endif
348
349 void setAllValues(const typename local_matrix_type::row_map_type &ptr,
350 const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type &ind,
351 const typename local_matrix_type::values_type &val);
352#endif // HAVE_XPETRA_TPETRA
353
355 LocalOrdinal GetStorageBlockSize() const { return mtx_->getBlockSize(); }
356
362 R.update(STS::one(), B, STS::zero());
363 this->apply(X, R, Teuchos::NO_TRANS, -STS::one(), STS::one());
364 }
365
366 private:
368
369}; // TpetraBlockCrsMatrix class
370
371} // namespace Xpetra
372
373#define XPETRA_TPETRABLOCKCRSMATRIX_SHORT
374#endif // XPETRA_TPETRABLOCKCRSMATRIX_DECL_HPP
static const EVerbosityLevel verbLevel_default
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...
virtual void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta)=0
Update multi-vector values with scaled values of A, this = beta*this + alpha*A.
void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using global IDs (not implemented)
void setObjectLabel(const std::string &objectLabel)
local_matrix_type::HostMirror getLocalMatrixHost() const
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const
Returns the Map that describes the row distribution in this matrix.
void setAllToScalar(const Scalar &alpha)
Set all matrix entries equal to scalarThis.
std::string description() const
A simple one-line description of this object.
virtual bool haveGlobalConstants() const
Returns true if globalConstants have been computed; false otherwise.
TpetraImport< LocalOrdinal, GlobalOrdinal, Node > TpetraImportClass
size_t getLocalNumCols() const
Returns the number of columns connected to the locally owned rows of this 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 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.
bool supportsRowViews() const
Returns true if getLocalRowView() and getGlobalRowView() are valid for this class.
bool isLocallyIndexed() const
If matrix indices are in the local range, this function returns true. Otherwise, this function return...
global_size_t getGlobalNumRows() const
Number of global elements in the row map of this matrix.
size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const
Returns the current number of entries on this node in the specified local row.
TpetraBlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraBlockCrsMatrixClass
void leftScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Left scale operator with given vector values.
void getAllValues(ArrayRCP< const size_t > &rowptr, ArrayRCP< const LocalOrdinal > &colind, ArrayRCP< const Scalar > &values) const
Gets the 1D pointer arrays of the graph (not implemented)
size_t getLocalNumRows() const
Returns the number of matrix rows owned on the calling node.
void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag) const
Get a copy of the diagonal entries owned by this node, with local row idices.
TpetraVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraVectorClass
bool hasMatrix() const
Does this have an underlying matrix.
void getLocalDiagOffsets(Teuchos::ArrayRCP< size_t > &offsets) const
Get offsets of the diagonal entries in 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 (not impelmented)
void replaceDiag(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag)
Replace the diagonal entries of the matrix.
global_size_t getGlobalNumEntries() const
Returns the global number of entries in this matrix.
bool isGloballyIndexed() const
If matrix indices are in the global range, this function returns true. Otherwise, this function retur...
TpetraExport< LocalOrdinal, GlobalOrdinal, Node > TpetraExportClass
void rightScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Right scale operator with given vector values.
RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_BlockCrsMatrix() const
Get the underlying Tpetra matrix.
void scale(const Scalar &alpha)
Scale the current values of a matrix, this = alpha*this (not implemented)
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const
Returns the Map that describes the column distribution in this matrix.
void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using local IDs (not implemented)
typename 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=Teuchos::NO_TRANS, Scalar alpha=ScalarTraits< Scalar >::one(), Scalar beta=ScalarTraits< Scalar >::zero()) const
Computes the sparse matrix-multivector multiplication.
RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > getCrsGraph() const
Returns the CrsGraph associated with this matrix.
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...
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....
size_t getLocalNumEntries() const
Returns the local number of entries in this 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.
ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const
Returns the Frobenius norm of the matrix.
LocalOrdinal GetStorageBlockSize() const
Returns the block size of the storage mechanism.
void resumeFill(const RCP< ParameterList > &params=null)
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import.
RCP< Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > mtx_
global_size_t getGlobalNumCols() const
Number of global columns in the matrix.
size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const
Returns the current number of entries in the (locally owned) global row.
bool isFillComplete() const
Returns true if the matrix is in compute mode, i.e. if fillComplete() has been called.
void replaceLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using local 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.
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 doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using global IDs (not implemented)
bool isFillActive() const
Returns true if the matrix is in edit mode.
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 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 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 getLocalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on this node.
RCP< Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_BlockCrsMatrixNonConst() const
Get the underlying Tpetra matrix.
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 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 removeEmptyProcessesInPlace(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newMap)
size_t global_size_t
Global size_t object.
CombineMode
Xpetra::Combine Mode enumerable type.