56 template<
class GlobalOrdinal,
class Node>
59 return *(epetraGraph->getEpetra_CrsGraph());
63 template<
class EpetraGlobalOrdinal,
class Node>
72 template<
class EpetraGlobalOrdinal,
class Node>
73 EpetraCrsGraphT<EpetraGlobalOrdinal, Node>::EpetraCrsGraphT(
const RCP<
const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap,
const RCP<
const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap,
size_t maxNumEntriesPerRow,
ProfileType pftype,
const Teuchos::RCP< Teuchos::ParameterList > &plist)
81 template<
class EpetraGlobalOrdinal,
class Node>
82 void EpetraCrsGraphT<EpetraGlobalOrdinal, Node>::insertGlobalIndices(GlobalOrdinal globalRow,
const ArrayView<const GlobalOrdinal> &indices) {
85 GlobalOrdinal* indices_rawPtr =
const_cast<GlobalOrdinal*
>(indices.getRawPtr());
86 XPETRA_ERR_CHECK(graph_->InsertGlobalIndices(globalRow, indices.size(), indices_rawPtr));
91 template<
class EpetraGlobalOrdinal,
class Node>
92 void EpetraCrsGraphT<EpetraGlobalOrdinal, Node>::insertLocalIndices(
int localRow,
const ArrayView<const int> &indices) {
95 int* indices_rawPtr =
const_cast<int*
>(indices.getRawPtr());
96 XPETRA_ERR_CHECK(graph_->InsertMyIndices(localRow, indices.size(), indices_rawPtr));
100 template<
class EpetraGlobalOrdinal,
class Node>
101 void EpetraCrsGraphT<EpetraGlobalOrdinal, Node>::getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView<const GlobalOrdinal> &Indices)
const {
105 GlobalOrdinal * eIndices;
107 XPETRA_ERR_CHECK(graph_->ExtractGlobalRowView(GlobalRow, numEntries, eIndices));
108 if (numEntries == 0) { eIndices = NULL; }
110 Indices = ArrayView<const GlobalOrdinal>(eIndices, numEntries);
114 template<
class EpetraGlobalOrdinal,
class Node>
115 void EpetraCrsGraphT<EpetraGlobalOrdinal, Node>::getLocalRowView(
int LocalRow, ArrayView<const int> &indices)
const {
122 if (numEntries == 0) { eIndices = NULL; }
124 indices = ArrayView<const int>(eIndices, numEntries);
128 template<
class EpetraGlobalOrdinal,
class Node>
129 void EpetraCrsGraphT<EpetraGlobalOrdinal, Node>::fillComplete(
const RCP<
const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap,
const RCP<
const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap,
const RCP< ParameterList > ¶ms){
132 graph_->FillComplete(toEpetra<EpetraGlobalOrdinal,Node>(domainMap), toEpetra<EpetraGlobalOrdinal,Node>(rangeMap));
133 bool doOptimizeStorage =
true;
134 if (params != null && params->get(
"Optimize Storage",
true) ==
false) doOptimizeStorage =
false;
135 if (doOptimizeStorage) graph_->OptimizeStorage();
139 template<
class EpetraGlobalOrdinal,
class Node>
140 void EpetraCrsGraphT<EpetraGlobalOrdinal, Node>::fillComplete(
const RCP< ParameterList > ¶ms) {
143 graph_->FillComplete();
144 bool doOptimizeStorage =
true;
145 if (params != null && params->get(
"Optimize Storage",
true) ==
false) doOptimizeStorage =
false;
146 if (doOptimizeStorage) graph_->OptimizeStorage();
150 template<
class EpetraGlobalOrdinal,
class Node>
151 std::string EpetraCrsGraphT<EpetraGlobalOrdinal, Node>::description()
const {
XPETRA_MONITOR(
"EpetraCrsGraphT::description");
return "NotImplemented"; }
153 template<
class EpetraGlobalOrdinal,
class Node>
157 out <<
"EpetraCrsGraphT::describe : Warning, verbosity level is ignored by this method." << std::endl;
159 if (rowmap.
Comm().
MyPID() == 0) out <<
"** EpetraCrsGraphT **\n\nrowmap" << std::endl;
166 template<
class GlobalOrdinal,
class Node>
167 RCP<const CrsGraph<int, GlobalOrdinal, Node> >
178 template<
class EpetraGlobalOrdinal,
class Node>
179 void EpetraCrsGraphT<EpetraGlobalOrdinal, Node>::doImport(
const DistObject<GlobalOrdinal, LocalOrdinal, GlobalOrdinal, Node> &source,
180 const Import<LocalOrdinal, GlobalOrdinal, Node> &importer,
CombineMode CM) {
183 XPETRA_DYNAMIC_CAST(
const EpetraCrsGraphT<GlobalOrdinal XPETRA_COMMA Node>, source, tSource,
"Xpetra::EpetraCrsGraphT::doImport only accept Xpetra::EpetraCrsGraphT as input arguments.");
184 XPETRA_DYNAMIC_CAST(
const EpetraImportT<GlobalOrdinal XPETRA_COMMA Node>, importer, tImporter,
"Xpetra::EpetraCrsGraphT::doImport only accept Xpetra::EpetraImportT as input arguments.");
186 RCP<const Epetra_CrsGraph> v = tSource.getEpetra_CrsGraph();
187 int err = graph_->Import(*v, *tImporter.getEpetra_Import(),
toEpetra(CM));
192 template<
class EpetraGlobalOrdinal,
class Node>
193 void EpetraCrsGraphT<EpetraGlobalOrdinal,Node>::doExport(
const DistObject<GlobalOrdinal, LocalOrdinal, GlobalOrdinal, Node> &dest,
194 const Import<LocalOrdinal, GlobalOrdinal, Node>& importer,
CombineMode CM) {
197 XPETRA_DYNAMIC_CAST(
const EpetraCrsGraphT<GlobalOrdinal XPETRA_COMMA Node>, dest, tDest,
"Xpetra::EpetraCrsGraphT::doImport only accept Xpetra::EpetraCrsGraphT as input arguments.");
198 XPETRA_DYNAMIC_CAST(
const EpetraImportT<GlobalOrdinal XPETRA_COMMA Node>, importer, tImporter,
"Xpetra::EpetraCrsGraphT::doImport only accept Xpetra::EpetraImportT as input arguments.");
200 RCP<const Epetra_CrsGraph> v = tDest.getEpetra_CrsGraph();
201 int err = graph_->Export(*v, *tImporter.getEpetra_Import(),
toEpetra(CM));
206 template<
class EpetraGlobalOrdinal,
class Node>
207 void EpetraCrsGraphT<EpetraGlobalOrdinal,Node>::doImport(
const DistObject<GlobalOrdinal, LocalOrdinal, GlobalOrdinal, Node> &source,
208 const Export<LocalOrdinal, GlobalOrdinal, Node>& exporter,
CombineMode CM) {
211 XPETRA_DYNAMIC_CAST(
const EpetraCrsGraphT<GlobalOrdinal XPETRA_COMMA Node>, source, tSource,
"Xpetra::EpetraCrsGraphT::doImport only accept Xpetra::EpetraCrsGraphT as input arguments.");
212 XPETRA_DYNAMIC_CAST(
const EpetraExportT<GlobalOrdinal XPETRA_COMMA Node>, exporter, tExporter,
"Xpetra::EpetraCrsGraphT::doImport only accept Xpetra::EpetraImportT as input arguments.");
214 RCP<const Epetra_CrsGraph> v = tSource.getEpetra_CrsGraph();
215 int err = graph_->Import(*v, *tExporter.getEpetra_Export(),
toEpetra(CM));
221 template<
class EpetraGlobalOrdinal,
class Node>
222 void EpetraCrsGraphT<EpetraGlobalOrdinal,Node>::doExport(
const DistObject<GlobalOrdinal, LocalOrdinal, GlobalOrdinal, Node> &dest,
223 const Export<LocalOrdinal, GlobalOrdinal, Node>& exporter,
CombineMode CM) {
226 XPETRA_DYNAMIC_CAST(
const EpetraCrsGraphT<GlobalOrdinal XPETRA_COMMA Node>, dest, tDest,
"Xpetra::EpetraCrsGraphT::doImport only accept Xpetra::EpetraCrsGraphT as input arguments.");
227 XPETRA_DYNAMIC_CAST(
const EpetraExportT<GlobalOrdinal XPETRA_COMMA Node>, exporter, tExporter,
"Xpetra::EpetraCrsGraphT::doImport only accept Xpetra::EpetraImportT as input arguments.");
229 RCP<const Epetra_CrsGraph> v = tDest.getEpetra_CrsGraph();
230 int err = graph_->Export(*v, *tExporter.getEpetra_Export(),
toEpetra(CM));
235 #ifndef XPETRA_EPETRA_NO_32BIT_GLOBAL_INDICES
236 #ifdef HAVE_XPETRA_TPETRA
237 #include "TpetraCore_config.h"
238 #if ((defined(EPETRA_HAVE_OMP) && !defined(HAVE_TPETRA_INST_OPENMP)) || \
239 (!defined(EPETRA_HAVE_OMP) && !defined(HAVE_TPETRA_INST_SERIAL)))
240 template class EpetraCrsGraphT<int, Xpetra::EpetraNode >;
241 template RCP< const CrsGraph<int, int, Xpetra::EpetraNode > > toXpetra<int, Xpetra::EpetraNode>(
const Epetra_CrsGraph &g);
242 template const Epetra_CrsGraph & toEpetra<int, Xpetra::EpetraNode >(
const RCP< const CrsGraph<int, int, Xpetra::EpetraNode > > &graph);
244 #ifdef HAVE_TPETRA_INST_SERIAL
245 template class EpetraCrsGraphT<int, Kokkos::Compat::KokkosSerialWrapperNode >;
246 template RCP< const CrsGraph<int, int, Kokkos::Compat::KokkosSerialWrapperNode > > toXpetra<int, Kokkos::Compat::KokkosSerialWrapperNode>(
const Epetra_CrsGraph &g);
247 template const Epetra_CrsGraph & toEpetra<int, Kokkos::Compat::KokkosSerialWrapperNode >(
const RCP< const CrsGraph<int, int, Kokkos::Compat::KokkosSerialWrapperNode > > &graph);
249 #ifdef HAVE_TPETRA_INST_PTHREAD
250 template class EpetraCrsGraphT<int, Kokkos::Compat::KokkosThreadsWrapperNode>;
251 template RCP< const CrsGraph<int, int, Kokkos::Compat::KokkosThreadsWrapperNode > > toXpetra<int, Kokkos::Compat::KokkosThreadsWrapperNode>(
const Epetra_CrsGraph &g);
252 template const Epetra_CrsGraph & toEpetra<int, Kokkos::Compat::KokkosThreadsWrapperNode >(
const RCP< const CrsGraph<int, int, Kokkos::Compat::KokkosThreadsWrapperNode > > &graph);
254 #ifdef HAVE_TPETRA_INST_OPENMP
255 template class EpetraCrsGraphT<int, Kokkos::Compat::KokkosOpenMPWrapperNode >;
256 template RCP< const CrsGraph<int, int, Kokkos::Compat::KokkosOpenMPWrapperNode > > toXpetra<int, Kokkos::Compat::KokkosOpenMPWrapperNode>(
const Epetra_CrsGraph &g);
257 template const Epetra_CrsGraph & toEpetra<int, Kokkos::Compat::KokkosOpenMPWrapperNode >(
const RCP< const CrsGraph<int, int, Kokkos::Compat::KokkosOpenMPWrapperNode > > &graph);
259 #ifdef HAVE_TPETRA_INST_CUDA
260 typedef Kokkos::Compat::KokkosCudaWrapperNode default_node_type;
261 template class EpetraCrsGraphT<int, default_node_type >;
262 template RCP< const CrsGraph<int, int, default_node_type > > toXpetra<int, default_node_type>(
const Epetra_CrsGraph &g);
263 template const Epetra_CrsGraph & toEpetra<int, default_node_type >(
const RCP< const CrsGraph<int, int, default_node_type > > &graph);
268 template class EpetraCrsGraphT<int, default_node_type >;
269 template RCP< const CrsGraph<int, int, default_node_type > > toXpetra<int, default_node_type>(
const Epetra_CrsGraph &g);
270 template const Epetra_CrsGraph & toEpetra<int, default_node_type >(
const RCP< const CrsGraph<int, int, default_node_type > > &graph);
271 #endif // HAVE_XPETRA_TPETRA
274 #ifndef XPETRA_EPETRA_NO_64BIT_GLOBAL_INDICES
275 #ifdef HAVE_XPETRA_TPETRA
276 #include "TpetraCore_config.h"
277 #if ((defined(EPETRA_HAVE_OMP) && !defined(HAVE_TPETRA_INST_OPENMP)) || \
278 (!defined(EPETRA_HAVE_OMP) && !defined(HAVE_TPETRA_INST_SERIAL)))
279 template class EpetraCrsGraphT<long long, Xpetra::EpetraNode >;
280 template RCP< const CrsGraph<int, long long, Xpetra::EpetraNode > > toXpetra<long long, Xpetra::EpetraNode>(
const Epetra_CrsGraph &g);
281 template const Epetra_CrsGraph & toEpetra<long long, Xpetra::EpetraNode >(
const RCP< const CrsGraph<int, long long, Xpetra::EpetraNode > > &graph);
283 #ifdef HAVE_TPETRA_INST_SERIAL
284 template class EpetraCrsGraphT<long long, Kokkos::Compat::KokkosSerialWrapperNode >;
285 template RCP< const CrsGraph<int, long long, Kokkos::Compat::KokkosSerialWrapperNode > > toXpetra<long long, Kokkos::Compat::KokkosSerialWrapperNode>(
const Epetra_CrsGraph &g);
286 template const Epetra_CrsGraph & toEpetra<long long, Kokkos::Compat::KokkosSerialWrapperNode >(
const RCP< const CrsGraph<int, long long, Kokkos::Compat::KokkosSerialWrapperNode > > &graph);
288 #ifdef HAVE_TPETRA_INST_PTHREAD
289 template class EpetraCrsGraphT<long long, Kokkos::Compat::KokkosThreadsWrapperNode>;
290 template RCP< const CrsGraph<int, long long, Kokkos::Compat::KokkosThreadsWrapperNode > > toXpetra<long long, Kokkos::Compat::KokkosThreadsWrapperNode>(
const Epetra_CrsGraph &g);
291 template const Epetra_CrsGraph & toEpetra<long long, Kokkos::Compat::KokkosThreadsWrapperNode >(
const RCP< const CrsGraph<int, long long, Kokkos::Compat::KokkosThreadsWrapperNode > > &graph);
293 #ifdef HAVE_TPETRA_INST_OPENMP
294 template class EpetraCrsGraphT<long long, Kokkos::Compat::KokkosOpenMPWrapperNode >;
295 template RCP< const CrsGraph<int, long long, Kokkos::Compat::KokkosOpenMPWrapperNode > > toXpetra<long long, Kokkos::Compat::KokkosOpenMPWrapperNode>(
const Epetra_CrsGraph &g);
296 template const Epetra_CrsGraph & toEpetra<long long, Kokkos::Compat::KokkosOpenMPWrapperNode >(
const RCP< const CrsGraph<int, long long, Kokkos::Compat::KokkosOpenMPWrapperNode > > &graph);
298 #ifdef HAVE_TPETRA_INST_CUDA
299 typedef Kokkos::Compat::KokkosCudaWrapperNode default_node_type;
300 template class EpetraCrsGraphT<long long, default_node_type >;
301 template RCP< const CrsGraph<int, long long, default_node_type > > toXpetra<long long, default_node_type>(
const Epetra_CrsGraph &g);
302 template const Epetra_CrsGraph & toEpetra<long long, default_node_type >(
const RCP< const CrsGraph<int, long long, default_node_type > > &graph);
307 template class EpetraCrsGraphT<long long, default_node_type >;
308 template RCP< const CrsGraph<int, long long, default_node_type > > toXpetra<long long, default_node_type>(
const Epetra_CrsGraph &g);
309 template const Epetra_CrsGraph & toEpetra<long long, default_node_type >(
const RCP< const CrsGraph<int, long long, default_node_type > > &graph);
310 #endif // HAVE_XPETRA_TPETRA
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
const Epetra_CrsGraph & toEpetra(const RCP< const CrsGraph< int, GlobalOrdinal, Node > > &graph)
virtual int MyPID() const =0
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
#define XPETRA_ERR_CHECK(arg)
#define XPETRA_DYNAMIC_CAST(type, obj, newObj, exceptionMsg)
const Epetra_Comm & Comm() const
#define XPETRA_RCP_DYNAMIC_CAST(type, obj, newObj, exceptionMsg)
RCP< const CrsGraph< int, GlobalOrdinal, Node > > toXpetra(const Epetra_CrsGraph &g)
EpetraCrsGraphT(const RCP< const map_type > &rowMap, size_t maxNumEntriesPerRow, ProfileType pftype=DynamicProfile, const RCP< Teuchos::ParameterList > &plist=Teuchos::null)
Constructor specifying fixed number of entries for each row.
CombineMode
Xpetra::Combine Mode enumerable type.
#define XPETRA_MONITOR(funcName)
virtual void Print(std::ostream &os) const