Tpetra Matrix/Vector Services Version of the Day
Functions
Tpetra::MatrixMatrix Namespace Reference

Collection of matrix-matrix operations. More...

Functions

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class SpMatOps >
void Multiply (const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps > &A, bool transposeA, const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps > &B, bool transposeB, CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps > &C, bool call_FillComplete_on_result=true)
template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class SpMatOps >
void Add (const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps > &A, bool transposeA, Scalar scalarA, CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps > &B, Scalar scalarB)
template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class SpMatOps >
void Add (const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps > &A, bool transposeA, Scalar scalarA, const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps > &B, bool transposeB, Scalar scalarB, RCP< CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps > > C)

Detailed Description

Collection of matrix-matrix operations.


Function Documentation

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class SpMatOps >
void Tpetra::MatrixMatrix::Multiply ( const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps > &  A,
bool  transposeA,
const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps > &  B,
bool  transposeB,
CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps > &  C,
bool  call_FillComplete_on_result = true 
)

Given CrsMatrix objects A, B and C, form the product C = A*B. In a parallel setting, A and B need not have matching distributions, but C needs to have the same row-map as A. At this time C=AT*B and C=A*BT are known to not work. However, C=A*B and C=AT*BT are known to work, Kurtis Nusbaum 03/24/2011

Parameters:
AInput, must already have had 'FillComplete()' called.
transposeAInput, whether to use transpose of matrix A.
BInput, must already have had 'FillComplete()' called.
transposeBInput, whether to use transpose of matrix B.
CResult. On entry to this method, it doesn't matter whether FillComplete() has already been called on C or not. If it has, then C's graph must already contain all nonzero locations that will be produced when forming the product A*B. On exit, C.FillComplete() will have been called, unless the last argument to this function is specified to be false.
call_FillComplete_on_resultOptional argument, defaults to true. Power users may specify this argument to be false if they *DON'T* want this function to call C.FillComplete. (It is often useful to allow this function to call C.FillComplete, in cases where one or both of the input matrices are rectangular and it is not trivial to know which maps to use for the domain- and range-maps.)

Definition at line 61 of file TpetraExt_MatrixMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class SpMatOps >
void Tpetra::MatrixMatrix::Add ( const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps > &  A,
bool  transposeA,
Scalar  scalarA,
CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps > &  B,
Scalar  scalarB 
)

Given CrsMatrix objects A and B, form the sum B = a*A + b*B Currently not functional.

Parameters:
AInput, must already have had 'FillComplete()' called.
transposeAInput, whether to use transpose of matrix A.
scalarAInput, scalar multiplier for matrix A.
BResult. On entry to this method, it doesn't matter whether FillComplete() has already been called on B or not. If it has, then B's graph must already contain all nonzero locations that will be produced when forming the sum.
scalarBInput, scalar multiplier for matrix B.

Definition at line 210 of file TpetraExt_MatrixMatrix_def.hpp.

template<class Scalar , class LocalOrdinal , class GlobalOrdinal , class Node , class SpMatOps >
void Tpetra::MatrixMatrix::Add ( const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps > &  A,
bool  transposeA,
Scalar  scalarA,
const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps > &  B,
bool  transposeB,
Scalar  scalarB,
RCP< CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps > >  C 
)

Given CrsMatrix objects A and B, form the sum C = a*A + b*B

Parameters:
AInput, must already have had 'FillComplete()' called.
transposeAInput, whether to use transpose of matrix A.
scalarAInput, scalar multiplier for matrix A.
BInput, must already have had 'FillComplete()' called.
transposeBInput, whether to use transpose of matrix B.
scalarBInput, scalar multiplier for matrix B.
CResult. On entry to this method, C can be NULL or a pointer to an unfilled or filled matrix. If C is NULL then a new object is allocated and must be deleted by the user. If C is not NULL and FillComplete has already been called then the sparsity pattern is assumed to be fixed and compatible with the sparsity of A+B. If FillComplete has not been called then the sum is completed and the function returns without calling FillComplete on C.

Definition at line 274 of file TpetraExt_MatrixMatrix_def.hpp.

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines