42 #ifndef TPETRA_MATRIXMATRIX_DECL_HPP 43 #define TPETRA_MATRIXMATRIX_DECL_HPP 46 #include <Teuchos_RCP.hpp> 47 #include <Teuchos_Array.hpp> 48 #include "Tpetra_ConfigDefs.hpp" 49 #include "Tpetra_CrsMatrix.hpp" 50 #include "Tpetra_Vector.hpp" 51 #include "TpetraExt_MMHelpers.hpp" 61 namespace MatrixMatrix {
93 template <
class Scalar,
98 const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& A,
100 const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& B,
102 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& C,
103 bool call_FillComplete_on_result=
true,
104 const std::string & label = std::string());
119 template <
class Scalar,
124 const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& A,
127 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& B,
160 template <
class Scalar,
164 Teuchos::RCP<CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
165 add (
const Scalar& alpha,
166 const bool transposeA,
167 const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& A,
169 const bool transposeB,
170 const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& B,
171 const Teuchos::RCP<
const Map<LocalOrdinal, GlobalOrdinal, Node> >& domainMap=Teuchos::null,
172 const Teuchos::RCP<
const Map<LocalOrdinal, GlobalOrdinal, Node> >& rangeMap=Teuchos::null,
173 const Teuchos::RCP<Teuchos::ParameterList>& params=Teuchos::null);
204 template <
class Scalar,
209 const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& A,
212 const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& B,
215 RCP<CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > C);
239 template <
class Scalar,
244 const Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node> & Dinv,
245 const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& A,
246 const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& B,
247 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& C,
248 bool call_FillComplete_on_result=
true,
249 const std::string & label = std::string());
255 template<
class Scalar,
259 void mult_AT_B_newmatrix(
260 const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& A,
261 const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& B,
262 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& C,
263 const std::string & label = std::string());
266 template<
class Scalar,
271 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Aview,
272 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Bview,
273 CrsWrapper<Scalar, LocalOrdinal, GlobalOrdinal, Node>& C,
274 const std::string & label = std::string());
276 template<
class Scalar,
280 void mult_A_B_newmatrix(
281 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Aview,
282 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Bview,
283 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& C,
284 const std::string & label = std::string());
287 template<
class Scalar,
291 void jacobi_A_B_newmatrix(
293 const Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node> & Dinv,
294 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Aview,
295 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Bview,
296 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& C,
297 const std::string & label = std::string());
300 template<
class Scalar,
304 void import_and_extract_views(
305 const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& M,
306 RCP<
const Map<LocalOrdinal, GlobalOrdinal, Node> > targetMap,
307 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Mview,
308 RCP<
const Import<LocalOrdinal,GlobalOrdinal, Node> > prototypeImporter = Teuchos::null,
309 bool userAssertsThereAreNoRemotes=
false,
310 const std::string & label = std::string());
312 template<
class Scalar,
316 void setMaxNumEntriesPerRow(
317 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Mview);
323 #endif // TPETRA_MATRIXMATRIX_DECL_HPP Namespace Tpetra contains the class and methods constituting the Tpetra library.
void Add(const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, bool transposeA, Scalar scalarA, CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, Scalar scalarB)
void Jacobi(Scalar omega, const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Dinv, const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &C, bool call_FillComplete_on_result=true, const std::string &label=std::string())
void Multiply(const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, bool transposeA, const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, bool transposeB, CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &C, bool call_FillComplete_on_result=true, const std::string &label=std::string())
Sparse matrix-matrix multiply.
Teuchos::RCP< CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > add(const Scalar &alpha, const bool transposeA, const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta, const bool transposeB, const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, 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)
Compute the sparse matrix sum C = scalarA * Op(A) + scalarB * Op(B), where Op(X) is either X or its t...