PlayaEpetraMatrixMatrixSum.cpp
00001
00002
00003
00004
00005 #include "PlayaEpetraMatrix.hpp"
00006 #include "PlayaEpetraMatrixMatrixSum.hpp"
00007 #include "PlayaEpetraVector.hpp"
00008 #include "PlayaExceptions.hpp"
00009 #include "EpetraExt_MatrixMatrix.h"
00010
00011 #ifndef HAVE_TEUCHOS_EXPLICIT_INSTANTIATION
00012 #include "PlayaVectorImpl.hpp"
00013 #include "PlayaLinearOperatorImpl.hpp"
00014 #endif
00015
00016
00017
00018 namespace Playa
00019 {
00020 using namespace Teuchos;
00021
00022
00023 LinearOperator<double> epetraMatrixMatrixSum(
00024 const LinearOperator<double>& A,
00025 const LinearOperator<double>& B)
00026 {
00027
00028
00029 RCP<const Epetra_CrsMatrix> A_crs = EpetraMatrix::getConcretePtr(A);
00030
00031
00032
00033 RCP<const Epetra_CrsMatrix> B_crs = EpetraMatrix::getConcretePtr(B);
00034
00035 bool transA = false;
00036 bool transB = false;
00037
00038 TEUCHOS_TEST_FOR_EXCEPTION(A.range() != B.range(), RuntimeError,
00039 "incompatible operand ranges in epetraMatrixMatrixSum()"
00040 << std::endl << "A.range()=" << A.range()
00041 << std::endl << "B.range()=" << B.range()
00042 );
00043
00044
00045 TEUCHOS_TEST_FOR_EXCEPTION(A.domain() != B.domain(), RuntimeError,
00046 "incompatible operand domains in epetraMatrixMatrixSum()"
00047 << std::endl << "A.domain()=" << A.domain()
00048 << std::endl << "B.domain()=" << B.domain()
00049 );
00050
00051
00052
00053 const Epetra_Map* rowmap
00054 = transA ? &(A_crs->DomainMap()) : &(A_crs->RowMap());
00055
00056
00057 RCP<Epetra_CrsMatrix> C = rcp(new Epetra_CrsMatrix(Copy, *rowmap, 1));
00058 Epetra_CrsMatrix* CPtr = C.get();
00059
00060
00061 int ierr
00062 = EpetraExt::MatrixMatrix::Add(
00063 *A_crs, transA, 1.0,
00064 *B_crs, transB, 1.0, CPtr);
00065 TEUCHOS_TEST_FOR_EXCEPTION(ierr != 0, RuntimeError,
00066 "EpetraExt Matrix-matrix add failed with error code ierr=" << ierr);
00067
00068
00069 C->FillComplete();
00070
00071
00072 RCP<LinearOperatorBase<double> > rtn
00073 = rcp(new EpetraMatrix(C, B.domain(), A.range()));
00074 return rtn;
00075
00076 }
00077
00078 }