44#ifndef TPETRA_EPETRAROWMATRIX_HPP
45#define TPETRA_EPETRAROWMATRIX_HPP
47#include "TpetraCore_config.h"
49#if defined(HAVE_TPETRA_EPETRA)
51#include <Epetra_Comm.h>
52#include <Epetra_BasicRowMatrix.h>
53#include <Tpetra_CrsMatrix.hpp>
54#include <Teuchos_TestForException.hpp>
67std::shared_ptr<Epetra_Comm>
68makeEpetraCommFromTeuchosComm (
const Teuchos::Comm<int>& teuchosComm);
75template<
class EpetraGlobalOrdinalType,
class TpetraMapType>
77tpetraToEpetraMapTmpl (
const TpetraMapType& tpetraMap)
79 using Teuchos::ArrayView;
80 typedef typename TpetraMapType::global_ordinal_type TGO;
81 typedef typename TpetraMapType::local_ordinal_type LO;
82 typedef EpetraGlobalOrdinalType EGO;
84 const TGO gblNumInds =
static_cast<TGO
> (tpetraMap.getGlobalNumElements ());
85 const LO lclNumInds =
static_cast<LO
> (tpetraMap.getLocalNumElements ());
86 ArrayView<const TGO> global_index_list = tpetraMap.getLocalElementList ();
88 std::vector<EGO> global_index_list_epetra;
89 const EGO* global_index_list_epetra_ptr = NULL;
90 if (std::is_same<TGO, EGO>::value) {
91 global_index_list_epetra_ptr =
92 reinterpret_cast<const EGO*
> (global_index_list.getRawPtr ());
95 global_index_list_epetra.resize (lclNumInds);
96 for (LO k = 0; k < lclNumInds; ++k) {
98 global_index_list_epetra[k] =
static_cast<EGO
> (global_index_list[k]);
100 global_index_list_epetra_ptr = global_index_list_epetra.data ();
102 const EGO indexBase = tpetraMap.getIndexBase ();
103 std::shared_ptr<Epetra_Comm> epetraComm =
104 Tpetra::Details::makeEpetraCommFromTeuchosComm (* (tpetraMap.getComm ()));
112 return Epetra_Map (
static_cast<EGO
> (gblNumInds),
113 static_cast<int> (lclNumInds),
114 global_index_list_epetra_ptr, indexBase, *epetraComm);
122template<
class TpetraMatrixType>
123class EpetraRowMatrix :
public Epetra_BasicRowMatrix {
125 EpetraRowMatrix(
const Teuchos::RCP<TpetraMatrixType> &mat,
const Epetra_Comm &comm);
126 virtual ~EpetraRowMatrix() {};
128 int ExtractMyRowCopy(
int MyRow,
int Length,
int & NumEntries,
double *Values,
int * Indices)
const;
131 int ExtractMyEntryView(
int CurEntry,
double * & Value,
int & RowIndex,
int & ColIndex);
134 int ExtractMyEntryView(
int CurEntry,
double const * & Value,
int & RowIndex,
int & ColIndex)
const;
136 int NumMyRowEntries(
int MyRow,
int & NumEntries)
const;
139 Teuchos::RCP<TpetraMatrixType> tpetra_matrix_;
142template<
class TpetraMatrixType>
143EpetraRowMatrix<TpetraMatrixType>::EpetraRowMatrix(
144 const Teuchos::RCP<TpetraMatrixType> &mat,
const Epetra_Comm &comm
146 : Epetra_BasicRowMatrix(comm),
150 typedef typename TpetraMatrixType::map_type tpetra_map_type;
151#if ! defined(EPETRA_NO_64BIT_GLOBAL_INDICES)
153 typedef long long EGO;
154#elif ! defined(EPETRA_NO_32BIT_GLOBAL_INDICES)
158# error "Epetra was not configured correctly. Neither 64-bit nor 32-bit indices are enabled."
160 const char tfecfFuncName[] =
"EpetraRowMatrix: ";
162 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC
163 (mat.is_null (), std::invalid_argument,
164 "The input Tpetra matrix is null.");
165 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC
166 (mat->getRowMap ().is_null (), std::invalid_argument,
167 "The input Tpetra matrix's row Map is null.");
168 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC
169 (mat->getColMap ().is_null (), std::invalid_argument,
170 "The input Tpetra matrix's column Map is null.");
172 RCP<const tpetra_map_type> tpetraRowMap = mat->getRowMap ();
173 Epetra_Map epetraRowMap =
174 tpetraToEpetraMapTmpl<EGO, tpetra_map_type> (*tpetraRowMap);
175 RCP<const tpetra_map_type> tpetraColMap = mat->getColMap ();
176 Epetra_Map epetraColMap =
177 tpetraToEpetraMapTmpl<EGO, tpetra_map_type> (*tpetraColMap);
178 this->SetMaps (epetraRowMap, epetraColMap);
181template<
class TpetraMatrixType>
182int EpetraRowMatrix<TpetraMatrixType>::ExtractMyRowCopy(
int MyRow,
int Length,
int & NumEntries,
double *Values,
int * Indices)
const
184 using inds_view =
typename TpetraMatrixType::nonconst_local_inds_host_view_type;
185 using vals_view =
typename TpetraMatrixType::nonconst_values_host_view_type;
186 static_assert (std::is_same<typename TpetraMatrixType::scalar_type, double>::value,
187 "This code assumes that Tpetra::CrsMatrix's scalar_type is int.");
188 static_assert (std::is_same<typename TpetraMatrixType::local_ordinal_type, int>::value,
189 "This code assumes that Tpetra::CrsMatrix's local_ordinal_type is int.");
190 inds_view IndicesView(Indices, Length);
191 vals_view ValuesView(Values, Length);
192 size_t num_entries = NumEntries;
193 tpetra_matrix_->getLocalRowCopy(MyRow, IndicesView, ValuesView, num_entries);
194 NumEntries = num_entries;
198template<
class TpetraMatrixType>
199int EpetraRowMatrix<TpetraMatrixType>::ExtractMyEntryView(
int CurEntry,
double * & Value,
int & RowIndex,
int & ColIndex)
205template<
class TpetraMatrixType>
206int EpetraRowMatrix<TpetraMatrixType>::ExtractMyEntryView(
int CurEntry,
double const * & Value,
int & RowIndex,
int & ColIndex)
const
212template<
class TpetraMatrixType>
213int EpetraRowMatrix<TpetraMatrixType>::NumMyRowEntries(
int MyRow,
int & NumEntries)
const
215 NumEntries = tpetra_matrix_->getNumEntriesInLocalRow(MyRow);
Implementation details of Tpetra.
Namespace Tpetra contains the class and methods constituting the Tpetra library.