46 #ifndef XPETRA_TPETRACRSMATRIX_HPP 47 #define XPETRA_TPETRACRSMATRIX_HPP 57 #include "Tpetra_CrsMatrix.hpp" 74 template<class Scalar = CrsMatrix<>::scalar_type,
82 :
public CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>
92 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR 93 #ifdef HAVE_XPETRA_TPETRA 105 :
mtx_ (Teuchos::rcp (new Tpetra::
CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> (
toTpetra(rowMap), maxNumEntriesPerRow,
toTpetra(pftype), params))) { }
109 :
mtx_(Teuchos::rcp(new Tpetra::
CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> (
toTpetra(rowMap), NumEntriesPerRowToAlloc,
toTpetra(pftype), params))) { }
121 :
mtx_(Teuchos::rcp(new Tpetra::
CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >(
toTpetra(graph), params))) { }
130 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null)
132 typedef Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> MyTpetraCrsMatrix;
134 RCP< const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > v = tSourceMatrix.getTpetra_CrsMatrix();
136 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > myDomainMap = domainMap!=Teuchos::null ?
toTpetra(domainMap) : Teuchos::null;
137 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > myRangeMap = rangeMap!=Teuchos::null ?
toTpetra(rangeMap) : Teuchos::null;
138 mtx_=Tpetra::importAndFillCompleteCrsMatrix<MyTpetraCrsMatrix>(tSourceMatrix.getTpetra_CrsMatrix(),
toTpetra(importer),myDomainMap,myRangeMap,params);
139 bool restrictComm=
false;
140 if(!params.is_null()) restrictComm = params->get(
"Restrict Communicator",restrictComm);
141 if(restrictComm &&
mtx_->getRowMap().is_null())
mtx_=Teuchos::null;
150 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null)
152 typedef Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> MyTpetraCrsMatrix;
154 RCP< const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > v = tSourceMatrix.getTpetra_CrsMatrix();
156 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > myDomainMap = domainMap!=Teuchos::null ?
toTpetra(domainMap) : Teuchos::null;
157 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > myRangeMap = rangeMap!=Teuchos::null ?
toTpetra(rangeMap) : Teuchos::null;
158 mtx_=Tpetra::exportAndFillCompleteCrsMatrix<MyTpetraCrsMatrix>(tSourceMatrix.getTpetra_CrsMatrix(),
toTpetra(exporter),myDomainMap,myRangeMap,params);
168 const Teuchos::RCP<Teuchos::ParameterList>& params)
170 typedef Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> MyTpetraCrsMatrix;
172 RCP< const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > v = tSourceMatrix.getTpetra_CrsMatrix();
174 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > myDomainMap = domainMap!=Teuchos::null ?
toTpetra(domainMap) : Teuchos::null;
175 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > myRangeMap = rangeMap!=Teuchos::null ?
toTpetra(rangeMap) : Teuchos::null;
177 mtx_=Tpetra::importAndFillCompleteCrsMatrix<MyTpetraCrsMatrix>(tSourceMatrix.getTpetra_CrsMatrix(),
toTpetra(RowImporter),
toTpetra(*DomainImporter),myDomainMap,myRangeMap,params);
178 bool restrictComm=
false;
179 if(!params.is_null()) restrictComm = params->get(
"Restrict Communicator",restrictComm);
180 if(restrictComm &&
mtx_->getRowMap().is_null())
mtx_=Teuchos::null;
189 const Teuchos::RCP<Teuchos::ParameterList>& params)
191 typedef Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> MyTpetraCrsMatrix;
193 RCP< const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > v = tSourceMatrix.getTpetra_CrsMatrix();
195 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > myDomainMap = domainMap!=Teuchos::null ?
toTpetra(domainMap) : Teuchos::null;
196 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > myRangeMap = rangeMap!=Teuchos::null ?
toTpetra(rangeMap) : Teuchos::null;
198 mtx_=Tpetra::exportAndFillCompleteCrsMatrix<MyTpetraCrsMatrix>(tSourceMatrix.getTpetra_CrsMatrix(),
toTpetra(RowExporter),
toTpetra(*DomainExporter),myDomainMap,myRangeMap,params);
201 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR 202 #ifdef HAVE_XPETRA_TPETRA 225 const local_matrix_type& lclMatrix,
226 const Teuchos::RCP<Teuchos::ParameterList>& params = null)
227 :
mtx_(Teuchos::rcp(new Tpetra::
CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(
toTpetra(rowMap),
toTpetra(colMap), lclMatrix, params))) { }
231 const local_matrix_type& lclMatrix,
232 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& rowMap,
233 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& colMap,
234 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap = Teuchos::null,
235 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& rangeMap = Teuchos::null,
236 const Teuchos::RCP<Teuchos::ParameterList>& params = null)
237 :
mtx_(Teuchos::rcp(new Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(lclMatrix,
toTpetra(rowMap),
toTpetra(colMap),
toTpetra(domainMap),
toTpetra(rangeMap), params))) { }
250 void insertGlobalValues(GlobalOrdinal globalRow,
const ArrayView< const GlobalOrdinal > &cols,
const ArrayView< const Scalar > &vals) {
XPETRA_MONITOR(
"TpetraCrsMatrix::insertGlobalValues");
mtx_->insertGlobalValues(globalRow, cols, vals); }
253 void insertLocalValues(LocalOrdinal localRow,
const ArrayView< const LocalOrdinal > &cols,
const ArrayView< const Scalar > &vals) {
XPETRA_MONITOR(
"TpetraCrsMatrix::insertLocalValues");
mtx_->insertLocalValues(localRow, cols, vals); }
256 void replaceGlobalValues(GlobalOrdinal globalRow,
const ArrayView< const GlobalOrdinal > &cols,
const ArrayView< const Scalar > &vals) {
XPETRA_MONITOR(
"TpetraCrsMatrix::replaceGlobalValues");
mtx_->replaceGlobalValues(globalRow, cols, vals); }
261 const ArrayView<const LocalOrdinal> &cols,
262 const ArrayView<const Scalar> &vals)
265 typedef typename ArrayView<const LocalOrdinal>::size_type size_type;
266 const LocalOrdinal numValid =
267 mtx_->replaceLocalValues (localRow, cols, vals);
268 TEUCHOS_TEST_FOR_EXCEPTION(
269 static_cast<size_type> (numValid) != cols.size (), std::runtime_error,
270 "replaceLocalValues returned " << numValid <<
" != cols.size() = " <<
271 cols.size () <<
".");
282 void allocateAllValues(
size_t numNonZeros,ArrayRCP<size_t> & rowptr, ArrayRCP<LocalOrdinal> & colind, ArrayRCP<Scalar> & values)
283 {
XPETRA_MONITOR(
"TpetraCrsMatrix::allocateAllValues"); rowptr.resize(
getNodeNumRows()+1); colind.resize(numNonZeros); values.resize(numNonZeros);}
286 void setAllValues(
const ArrayRCP<size_t> & rowptr,
const ArrayRCP<LocalOrdinal> & colind,
const ArrayRCP<Scalar> & values)
287 {
XPETRA_MONITOR(
"TpetraCrsMatrix::setAllValues");
mtx_->setAllValues(rowptr,colind,values); }
290 void getAllValues(ArrayRCP<const size_t>& rowptr, ArrayRCP<const LocalOrdinal>& colind, ArrayRCP<const Scalar>& values)
const 291 {
XPETRA_MONITOR(
"TpetraCrsMatrix::getAllValues");
mtx_->getAllValues(rowptr,colind,values); }
294 {
return mtx_->haveGlobalConstants();}
315 RCP<const Tpetra::Import<LocalOrdinal,GlobalOrdinal,Node> > myImport = tImporter.getTpetra_Import();
316 mtx_->replaceDomainMapAndImporter(
toTpetra(newDomainMap),myImport);
324 const RCP<ParameterList> ¶ms=Teuchos::null) {
326 RCP<const Tpetra::Import<LocalOrdinal,GlobalOrdinal,Node> > myImport;
327 RCP<const Tpetra::Export<LocalOrdinal,GlobalOrdinal,Node> > myExport;
329 if(importer!=Teuchos::null) {
331 myImport = tImporter.getTpetra_Import();
333 if(exporter!=Teuchos::null) {
335 myExport = tExporter.getTpetra_Export();
338 mtx_->expertStaticFillComplete(
toTpetra(domainMap),
toTpetra(rangeMap),myImport,myExport,params);
401 void getLocalRowCopy(LocalOrdinal LocalRow,
const ArrayView< LocalOrdinal > &Indices,
const ArrayView< Scalar > &Values,
size_t &NumEntries)
const {
XPETRA_MONITOR(
"TpetraCrsMatrix::getLocalRowCopy");
mtx_->getLocalRowCopy(LocalRow, Indices, Values, NumEntries); }
404 void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView< const GlobalOrdinal > &indices, ArrayView< const Scalar > &values)
const {
XPETRA_MONITOR(
"TpetraCrsMatrix::getGlobalRowView");
mtx_->getGlobalRowView(GlobalRow, indices, values); }
407 void getGlobalRowCopy(GlobalOrdinal GlobalRow,
const ArrayView< GlobalOrdinal > &indices,
const ArrayView< Scalar > &values,
size_t &numEntries)
const {
XPETRA_MONITOR(
"TpetraCrsMatrix::getGlobalRowCopy");
mtx_->getGlobalRowCopy(GlobalRow, indices, values, numEntries); }
410 void getLocalRowView(LocalOrdinal LocalRow, ArrayView< const LocalOrdinal > &indices, ArrayView< const Scalar > &values)
const {
XPETRA_MONITOR(
"TpetraCrsMatrix::getLocalRowView");
mtx_->getLocalRowView(LocalRow, indices, values); }
418 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 {
XPETRA_MONITOR(
"TpetraCrsMatrix::apply");
mtx_->apply(
toTpetra(X),
toTpetra(Y), mode, alpha, beta); }
435 void describe(Teuchos::FancyOStream &out,
const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default)
const {
XPETRA_MONITOR(
"TpetraCrsMatrix::describe");
mtx_->describe(out, verbLevel); }
447 mtx_->getLocalDiagCopy(*tDiag.getTpetra_Vector());
454 mtx_->getLocalDiagOffsets(offsets);
489 RCP< const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > v = tSource.getTpetra_CrsMatrix();
500 RCP< const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > v = tDest.getTpetra_CrsMatrix();
511 RCP< const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > v = tSource.getTpetra_CrsMatrix();
522 RCP< const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > v = tDest.getTpetra_CrsMatrix();
534 template<
class Node2>
535 RCP<CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node2> >
536 clone (
const RCP<Node2> &node2)
const 546 return !
mtx_.is_null ();
550 TpetraCrsMatrix(
const Teuchos::RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &mtx) :
mtx_(mtx) { }
558 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR 559 #ifdef HAVE_XPETRA_TPETRA 560 local_matrix_type getLocalMatrix ()
const {
565 void setAllValues (
const typename local_matrix_type::row_map_type& ptr,
566 const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type& ind,
567 const typename local_matrix_type::values_type& val) {
576 RCP< Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
mtx_;
579 #ifdef HAVE_XPETRA_EPETRA 581 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \ 582 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT)))) 585 template<
class Scalar>
587 :
public CrsMatrix<Scalar,int,int,EpetraNode>
599 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR 600 #ifdef HAVE_XPETRA_TPETRA 642 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null) {
651 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null) {
661 const Teuchos::RCP<Teuchos::ParameterList>& params)
672 const Teuchos::RCP<Teuchos::ParameterList>& params)
677 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR 678 #ifdef HAVE_XPETRA_TPETRA 701 const local_matrix_type& lclMatrix,
702 const Teuchos::RCP<Teuchos::ParameterList>& params = null) {
708 const local_matrix_type& lclMatrix,
709 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& rowMap,
710 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& colMap,
711 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap = Teuchos::null,
712 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& rangeMap = Teuchos::null,
713 const Teuchos::RCP<Teuchos::ParameterList>& params = null) {
739 const ArrayView<const LocalOrdinal> &cols,
740 const ArrayView<const Scalar> &vals)
751 void allocateAllValues(
size_t numNonZeros,ArrayRCP<size_t> & rowptr, ArrayRCP<LocalOrdinal> & colind, ArrayRCP<Scalar> & values) { }
754 void setAllValues(
const ArrayRCP<size_t> & rowptr,
const ArrayRCP<LocalOrdinal> & colind,
const ArrayRCP<Scalar> & values) { }
757 void getAllValues(ArrayRCP<const size_t>& rowptr, ArrayRCP<const LocalOrdinal>& colind, ArrayRCP<const Scalar>& values)
const { }
785 const RCP<ParameterList> ¶ms=Teuchos::null) { }
793 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getRowMap()
const {
return Teuchos::null; }
796 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getColMap()
const {
return Teuchos::null; }
799 RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > >
getCrsGraph()
const {
return Teuchos::null; }
841 typename ScalarTraits< Scalar >::magnitudeType
getFrobeniusNorm()
const {
return ScalarTraits< Scalar >::magnitude(ScalarTraits< Scalar >::zero()); }
847 void getLocalRowCopy(
LocalOrdinal LocalRow,
const ArrayView< LocalOrdinal > &Indices,
const ArrayView< Scalar > &Values,
size_t &NumEntries)
const { }
864 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 { }
867 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getDomainMap()
const {
return Teuchos::null; }
870 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getRangeMap()
const {
return Teuchos::null; }
881 void describe(Teuchos::FancyOStream &out,
const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default)
const { }
907 Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getMap()
const {
return Teuchos::null; }
929 template<
class Node2>
930 RCP<CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node2> >
931 clone (
const RCP<Node2> &node2)
const {
return Teuchos::null; }
940 TpetraCrsMatrix(
const Teuchos::RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &mtx) {
945 RCP<const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
getTpetra_CrsMatrix()
const {
return Teuchos::null; }
950 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR 951 #ifdef HAVE_XPETRA_TPETRA 952 local_matrix_type getLocalMatrix ()
const {
957 void setAllValues (
const typename local_matrix_type::row_map_type& ptr,
958 const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type& ind,
959 const typename local_matrix_type::values_type& val) {
969 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_LONG_LONG))) || \ 970 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_LONG_LONG)))) 973 template<
class Scalar>
975 :
public CrsMatrix<Scalar,int,long long,EpetraNode>
987 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR 988 #ifdef HAVE_XPETRA_TPETRA 1030 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null) {
1039 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null) {
1049 const Teuchos::RCP<Teuchos::ParameterList>& params)
1060 const Teuchos::RCP<Teuchos::ParameterList>& params)
1065 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR 1066 #ifdef HAVE_XPETRA_TPETRA 1089 const local_matrix_type& lclMatrix,
1090 const Teuchos::RCP<Teuchos::ParameterList>& params = null) {
1096 const local_matrix_type& lclMatrix,
1097 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& rowMap,
1098 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& colMap,
1099 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap = Teuchos::null,
1100 const Teuchos::RCP<
const Map<LocalOrdinal,GlobalOrdinal,Node> >& rangeMap = Teuchos::null,
1101 const Teuchos::RCP<Teuchos::ParameterList>& params = null) {
1102 XPETRA_TPETRA_ETI_EXCEPTION(
typeid(TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,EpetraNode>).name() ,
typeid(TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,EpetraNode>).name(),
"long long",
typeid(
EpetraNode).name() );
1127 const ArrayView<const LocalOrdinal> &cols,
1128 const ArrayView<const Scalar> &vals)
1139 void allocateAllValues(
size_t numNonZeros,ArrayRCP<size_t> & rowptr, ArrayRCP<LocalOrdinal> & colind, ArrayRCP<Scalar> & values) { }
1142 void setAllValues(
const ArrayRCP<size_t> & rowptr,
const ArrayRCP<LocalOrdinal> & colind,
const ArrayRCP<Scalar> & values) { }
1145 void getAllValues(ArrayRCP<const size_t>& rowptr, ArrayRCP<const LocalOrdinal>& colind, ArrayRCP<const Scalar>& values)
const { }
1172 const RCP<ParameterList> ¶ms=Teuchos::null) { }
1180 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getRowMap()
const {
return Teuchos::null; }
1183 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getColMap()
const {
return Teuchos::null; }
1186 RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > >
getCrsGraph()
const {
return Teuchos::null; }
1228 typename ScalarTraits< Scalar >::magnitudeType
getFrobeniusNorm()
const {
return ScalarTraits< Scalar >::magnitude(ScalarTraits< Scalar >::zero()); }
1234 void getLocalRowCopy(
LocalOrdinal LocalRow,
const ArrayView< LocalOrdinal > &Indices,
const ArrayView< Scalar > &Values,
size_t &NumEntries)
const { }
1251 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 { }
1254 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getDomainMap()
const {
return Teuchos::null; }
1257 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getRangeMap()
const {
return Teuchos::null; }
1268 void describe(Teuchos::FancyOStream &out,
const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default)
const { }
1294 Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > >
getMap()
const {
return Teuchos::null; }
1316 template<
class Node2>
1317 RCP<CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node2> >
1318 clone (
const RCP<Node2> &node2)
const {
return Teuchos::null; }
1327 TpetraCrsMatrix(
const Teuchos::RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &mtx) {
1332 RCP<const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
getTpetra_CrsMatrix()
const {
return Teuchos::null; }
1337 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR 1338 #ifdef HAVE_XPETRA_TPETRA 1339 local_matrix_type getLocalMatrix ()
const {
1344 void setAllValues (
const typename local_matrix_type::row_map_type& ptr,
1345 const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type& ind,
1346 const typename local_matrix_type::values_type& val) {
1356 #endif // HAVE_XPETRA_EPETRA 1360 #define XPETRA_TPETRACRSMATRIX_SHORT 1361 #endif // XPETRA_TPETRACRSMATRIX_HPP void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import.
RCP< CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node2 > > clone(const RCP< Node2 > &node2) const
TpetraCrsMatrix(const Teuchos::RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying a previously constructed graph.
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, 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)
Constructor for a fused export.
TpetraCrsMatrix(const TpetraCrsMatrix &matrix)
Deep copy constructor.
bool isFillComplete() const
Returns true if the matrix is in compute mode, i.e. if fillComplete() has been called.
void fillComplete(const RCP< ParameterList > ¶ms=null)
Signal that data entry is complete.
void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag, const Teuchos::ArrayView< const size_t > &offsets) const
Get a copy of the diagonal entries owned by this node, with local row indices.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, size_t maxNumEntriesPerRow, ProfileType pftype=DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying column Map and fixed number of entries for each row.
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import.
void setAllToScalar(const Scalar &alpha)
Set all matrix entries equal to scalarThis.
void getAllValues(ArrayRCP< const size_t > &rowptr, ArrayRCP< const LocalOrdinal > &colind, ArrayRCP< const Scalar > &values) const
Gets the 1D pointer arrays of the graph.
void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag, const Teuchos::ArrayView< const size_t > &offsets) const
Get a copy of the diagonal entries owned by this node, with local row indices.
TpetraExport< LocalOrdinal, GlobalOrdinal, Node > TpetraExportClass
TpetraVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraVectorClass
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const
Returns the Map that describes the column distribution in this matrix.
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.
size_t getNodeNumCols() const
Returns the number of columns connected to the locally owned rows of this matrix. ...
void removeEmptyProcessesInPlace(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newMap)
bool isFillComplete() const
Returns true if the matrix is in compute mode, i.e. if fillComplete() has been called.
TpetraImport< LocalOrdinal, GlobalOrdinal, Node > TpetraImportClass
Teuchos::RCP< NodeType > getNode()
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.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, size_t maxNumEntriesPerRow, ProfileType pftype=DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying fixed number of entries for each row.
void getLocalDiagOffsets(Teuchos::ArrayRCP< size_t > &offsets) const
Get offsets of the diagonal entries in the matrix.
void replaceLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using local IDs.
size_t getNodeNumRows() const
Returns the number of matrix rows owned on the calling node.
TpetraCrsMatrix(const Teuchos::RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying a previously constructed graph.
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
size_t getNodeNumRows() const
Returns the number of matrix rows owned on the calling node.
void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using local IDs.
TpetraExport< LocalOrdinal, GlobalOrdinal, Node > TpetraExportClass
void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using global IDs.
TpetraCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraCrsMatrixClass
TpetraCrsMatrix(const Teuchos::RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &mtx)
TpetraCrsMatrix constructor to wrap a Tpetra::CrsMatrix object.
global_size_t getGlobalNumEntries() const
Returns the global number of entries in this matrix.
GlobalOrdinal global_ordinal_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 Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrix() const
Get the underlying Tpetra matrix.
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.
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...
bool isLocallyIndexed() const
If matrix indices are in the local range, this function returns true. Otherwise, this function return...
RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > getCrsGraph() const
Returns the CrsGraph associated with this matrix.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, ProfileType pftype=DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying column Map and number of entries in each row.
RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > mtx_
RCP< CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node2 > > clone(const RCP< Node2 > &node2) const
void setAllValues(const ArrayRCP< size_t > &rowptr, const ArrayRCP< LocalOrdinal > &colind, const ArrayRCP< Scalar > &values)
Sets the 1D pointer arrays of the graph.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const
Returns the Map that describes the row distribution in this matrix.
void fillComplete(const RCP< ParameterList > ¶ms=null)
Signal that data entry is complete.
void rightScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Right scale operator with given vector values.
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 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.
size_t getGlobalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on all nodes.
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 replaceLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using local IDs.
void resumeFill(const RCP< ParameterList > ¶ms=null)
void replaceLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using local IDs.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
TpetraCrsMatrix(const Teuchos::RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying a previously constructed graph.
size_t getNodeNumEntries() const
Returns the local number of entries in this matrix.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, size_t maxNumEntriesPerRow, ProfileType pftype=DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying fixed number of entries for each row.
size_t getGlobalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on all nodes.
LocalOrdinal local_ordinal_type
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...
void setAllValues(const ArrayRCP< size_t > &rowptr, const ArrayRCP< LocalOrdinal > &colind, const ArrayRCP< Scalar > &values)
Sets the 1D pointer arrays of the graph.
RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrixNonConst() const
Get the underlying Tpetra matrix.
Kokkos::Compat::KokkosSerialWrapperNode EpetraNode
void removeEmptyProcessesInPlace(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newMap)
virtual ~TpetraCrsMatrix()
Destructor.
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Export (using an Importer).
void scale(const Scalar &alpha)
Scale the current values of a matrix, this = alpha*this.
global_size_t getGlobalNumRows() const
Number of global elements in the row map of this matrix.
#define XPETRA_TPETRA_ETI_EXCEPTION(cl, obj, go, node)
void fillComplete(const RCP< ParameterList > ¶ms=null)
Signal that data entry is complete.
size_t getNodeMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on this node.
bool haveGlobalConstants() const
Returns true if globalConstants have been computed; false otherwise.
size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const
Returns the current number of entries on this node in the specified local row.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const
Returns the Map that describes the column distribution in this matrix.
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...
bool supportsRowViews() const
Returns true if getLocalRowView() and getGlobalRowView() are valid for this class.
std::string description() const
A simple one-line description of this object.
void leftScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Left scale operator with given vector values.
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.
void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using local IDs.
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...
bool isFillComplete() const
Returns true if the matrix is in compute mode, i.e. if fillComplete() has been called.
size_t getNodeNumEntries() const
Returns the local number of entries in this matrix.
void setAllToScalar(const Scalar &alpha)
Set all matrix entries equal to scalarThis.
global_size_t getGlobalNumCols() const
Number of global columns in the matrix.
RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > getCrsGraph() const
Returns the CrsGraph associated with this matrix.
void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using local IDs.
void getAllValues(ArrayRCP< const size_t > &rowptr, ArrayRCP< const LocalOrdinal > &colind, ArrayRCP< const Scalar > &values) const
Gets the 1D pointer arrays of the graph.
TpetraVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraVectorClass
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, ProfileType pftype=DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying (possibly different) number of entries in each row.
bool haveGlobalConstants() const
Returns true if globalConstants have been computed; false otherwise.
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, 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)
Constructor for a fused import.
std::string description() const
A simple one-line description of this object.
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Export (using an Importer).
TpetraImport< LocalOrdinal, GlobalOrdinal, Node > TpetraImportClass
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.
bool supportsRowViews() const
Returns true if getLocalRowView() and getGlobalRowView() are valid for this class.
TpetraCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraCrsMatrixClass
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.
bool haveGlobalConstants() const
Returns true if globalConstants have been computed; false otherwise.
void fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > ¶ms=null)
Signal that data entry is complete, specifying domain and range maps.
virtual ~TpetraCrsMatrix()
Destructor.
bool isGloballyIndexed() const
If matrix indices are in the global range, this function returns true. Otherwise, this function retur...
bool isLocallyIndexed() const
If matrix indices are in the local range, this function returns true. Otherwise, this function return...
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Export< LocalOrdinal, GlobalOrdinal, Node > &RowExporter, const Teuchos::RCP< const Export< LocalOrdinal, GlobalOrdinal, Node > > DomainExporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > ¶ms)
Constructor for a fused export.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, ProfileType pftype=DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying (possibly different) number of entries in each row.
bool isGloballyIndexed() const
If matrix indices are in the global range, this function returns true. Otherwise, this function retur...
global_size_t getGlobalNumCols() const
Number of global columns in the matrix.
RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrixNonConst() const
Get the underlying Tpetra matrix.
bool isGloballyIndexed() const
If matrix indices are in the global range, this function returns true. Otherwise, this function retur...
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 getNodeNumCols() const
Returns the number of columns connected to the locally owned rows of this matrix. ...
RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrix() const
Get the underlying Tpetra matrix.
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
TpetraCrsMatrix(const Teuchos::RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &mtx)
TpetraCrsMatrix constructor to wrap a Tpetra::CrsMatrix object.
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Import (using an Exporter).
TpetraCrsMatrix(const Teuchos::RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &mtx)
TpetraCrsMatrix constructor to wrap a Tpetra::CrsMatrix object.
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.
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Import< LocalOrdinal, GlobalOrdinal, Node > &RowImporter, const Teuchos::RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > DomainImporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > ¶ms)
Constructor for a fused import.
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 removeEmptyProcessesInPlace(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newMap)
bool hasMatrix() const
Does this have an underlying matrix.
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, 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)
Constructor for a fused import.
size_t getNodeNumRows() const
Returns the number of matrix rows owned on the calling node.
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.
global_size_t getGlobalNumCols() const
Number of global columns in the matrix.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, ProfileType pftype=DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying column Map and number of entries in each row.
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 > ¶ms=Teuchos::null)
Expert static fill complete.
bool hasMatrix() const
Does this have an underlying 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 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 > ¶ms=Teuchos::null)
Expert static fill complete.
size_t getNodeNumCols() const
Returns the number of columns connected to the locally owned rows of this matrix. ...
void getAllValues(ArrayRCP< const size_t > &rowptr, ArrayRCP< const LocalOrdinal > &colind, ArrayRCP< const Scalar > &values) const
Gets the 1D pointer arrays of the graph.
TpetraVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraVectorClass
void resumeFill(const RCP< ParameterList > ¶ms=null)
void rightScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Right scale operator with given vector values.
bool isFillActive() const
Returns true if the matrix is in edit mode.
void scale(const Scalar &alpha)
Scale the current values of a matrix, this = alpha*this.
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.
#define XPETRA_DYNAMIC_CAST(type, obj, newObj, exceptionMsg)
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.
void getLocalDiagOffsets(Teuchos::ArrayRCP< size_t > &offsets) const
Get offsets of the diagonal entries in the matrix.
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Export< LocalOrdinal, GlobalOrdinal, Node > &RowExporter, const Teuchos::RCP< const Export< LocalOrdinal, GlobalOrdinal, Node > > DomainExporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > ¶ms)
Constructor for a fused export.
size_t global_size_t
Global size_t object.
RCP< CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node2 > > clone(const RCP< Node2 > &node2) const
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Export (using an Importer).
void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using global IDs.
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...
TpetraImport< LocalOrdinal, GlobalOrdinal, Node > TpetraImportClass
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, 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)
Constructor for a fused export.
RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrix() const
Get the underlying Tpetra matrix.
size_t getNodeMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on this node.
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Export< LocalOrdinal, GlobalOrdinal, Node > &RowExporter, const Teuchos::RCP< const Export< LocalOrdinal, GlobalOrdinal, Node > > DomainExporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > ¶ms)
Constructor for a fused export.
ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const
Returns the Frobenius norm of the matrix.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
global_size_t getGlobalNumRows() const
Number of global elements in the row map of this matrix.
void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag) const
Get a copy of the diagonal entries owned by this node, with local row idices.
RCP< const Tpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > toTpetra(const RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph)
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 > > getDomainMap() const
Returns the Map associated with the domain of this operator. This will be null until fillComplete() i...
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, ProfileType pftype=DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying (possibly different) number of entries in each row.
virtual ~TpetraCrsMatrix()
Destructor.
void setAllValues(const ArrayRCP< size_t > &rowptr, const ArrayRCP< LocalOrdinal > &colind, const ArrayRCP< Scalar > &values)
Sets the 1D pointer arrays of the graph.
RCP< const CrsGraph< int, GlobalOrdinal, Node > > toXpetra(const Epetra_CrsGraph &g)
void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using global IDs.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, size_t maxNumEntriesPerRow, ProfileType pftype=DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying column Map and fixed number of entries for each row.
bool hasMatrix() const
Does this have an underlying matrix.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const
Returns the Map that describes the row distribution in this matrix.
void rightScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Right scale operator with given vector values.
void resumeFill(const RCP< ParameterList > ¶ms=null)
bool isFillActive() const
Returns true if the matrix is in edit mode.
size_t getNodeMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on this node.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, size_t maxNumEntriesPerRow, ProfileType pftype=DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying column Map and fixed number of entries for each row.
size_t getNodeNumEntries() const
Returns the local number of entries in this matrix.
global_size_t getGlobalNumEntries() const
Returns the global number of entries in this matrix.
void scale(const Scalar &alpha)
Scale the current values of a matrix, this = alpha*this.
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Import< LocalOrdinal, GlobalOrdinal, Node > &RowImporter, const Teuchos::RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > DomainImporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > ¶ms)
Constructor for a fused import.
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 insertGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using global IDs.
void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag, const Teuchos::ArrayView< const size_t > &offsets) const
Get a copy of the diagonal entries owned by this node, with local row indices.
void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag) const
Get a copy of the diagonal entries owned by this node, with local row idices.
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 > ¶ms=Teuchos::null)
Expert static fill complete.
size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const
Returns the current number of entries on this node in the specified local row.
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.
void getLocalDiagOffsets(Teuchos::ArrayRCP< size_t > &offsets) const
Get offsets of the diagonal entries in the matrix.
size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const
Returns the current number of entries on this node in the specified local row.
ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const
Returns the Frobenius norm of the matrix.
CombineMode
Xpetra::Combine Mode enumerable type.
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 getGlobalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on all nodes.
bool isLocallyIndexed() const
If matrix indices are in the local range, this function returns true. Otherwise, this function return...
TpetraCrsMatrix(const TpetraCrsMatrix &matrix)
Deep copy constructor.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const
Returns the Map that describes the column distribution in this matrix.
#define XPETRA_MONITOR(funcName)
void leftScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Left scale operator with given vector values.
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Import (using an Exporter).
void fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > ¶ms=null)
Signal that data entry is complete, specifying domain and range maps.
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Import< LocalOrdinal, GlobalOrdinal, Node > &RowImporter, const Teuchos::RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > DomainImporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > ¶ms)
Constructor for a fused import.
global_size_t getGlobalNumRows() const
Number of global elements in the row map of this matrix.
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, 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)
Constructor for a fused export.
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.
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.
RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrixNonConst() const
Get the underlying Tpetra matrix.
global_size_t getGlobalNumEntries() const
Returns the global number of entries in this matrix.
void leftScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Left scale operator with given vector values.
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import.
std::string description() const
A simple one-line description of this object.
ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const
Returns the Frobenius norm of the matrix.
TpetraCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraCrsMatrixClass
void fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > ¶ms=null)
Signal that data entry is complete, specifying domain and range maps.
bool supportsRowViews() const
Returns true if getLocalRowView() and getGlobalRowView() are valid for this class.
TpetraCrsMatrix(const TpetraCrsMatrix &matrix)
Deep copy constructor.
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 insertGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using global IDs.
RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > getCrsGraph() const
Returns the CrsGraph associated with this matrix.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, size_t maxNumEntriesPerRow, ProfileType pftype=DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying fixed number of entries for each row.
void setAllToScalar(const Scalar &alpha)
Set all matrix entries equal to scalarThis.
TpetraCrsMatrix(const Teuchos::RCP< const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, 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)
Constructor for a fused import.
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...
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const
Returns the Map that describes the row distribution in this matrix.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, ProfileType pftype=DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Constructor specifying column Map and number of entries in each row.
TpetraExport< LocalOrdinal, GlobalOrdinal, Node > TpetraExportClass
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Import (using an Exporter).
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag) const
Get a copy of the diagonal entries owned by this node, with local row idices.
bool isFillActive() const
Returns true if the matrix is in edit mode.