00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046 #ifndef XPETRA_TPETRACRSGRAPH_HPP
00047 #define XPETRA_TPETRACRSGRAPH_HPP
00048
00049
00050
00051 #include "Xpetra_TpetraConfigDefs.hpp"
00052
00053 #include "Xpetra_CrsGraph.hpp"
00054
00055 #include "Tpetra_CrsGraph.hpp"
00056
00057 #include "Xpetra_Utils.hpp"
00058
00059 #include "Xpetra_TpetraMap.hpp"
00060
00061 #include "Xpetra_TpetraImport.hpp"
00062
00063 #include "Xpetra_TpetraExport.hpp"
00064
00065 namespace Xpetra {
00066
00067
00068 template <class LocalOrdinal, class GlobalOrdinal, class Node>
00069 RCP<const CrsGraph<LocalOrdinal, GlobalOrdinal, Node> >
00070 toXpetra (RCP<const Tpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > graph);
00071
00072 template <class LocalOrdinal, class GlobalOrdinal, class Node>
00073 RCP< const Tpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > >
00074 toTpetra (const RCP<const CrsGraph< LocalOrdinal, GlobalOrdinal, Node> >& graph);
00075
00076 template <class LocalOrdinal = CrsGraph<>::local_ordinal_type,
00077 class GlobalOrdinal = typename CrsGraph<LocalOrdinal>::global_ordinal_type,
00078 class Node = typename CrsGraph<LocalOrdinal, GlobalOrdinal>::node_type>
00079 class TpetraCrsGraph
00080 : public CrsGraph<LocalOrdinal,GlobalOrdinal,Node>
00081 {
00082
00083
00084 typedef TpetraCrsGraph<LocalOrdinal,GlobalOrdinal,Node> TpetraCrsGraphClass;
00085 typedef Map<LocalOrdinal,GlobalOrdinal,Node> map_type;
00086
00087 public:
00088
00090
00091
00093 TpetraCrsGraph(const RCP< const map_type > &rowMap, size_t maxNumEntriesPerRow, ProfileType pftype=DynamicProfile, const RCP< ParameterList > ¶ms=null)
00094 : graph_(Teuchos::rcp(new Tpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node >(toTpetra(rowMap), maxNumEntriesPerRow, toTpetra(pftype), params))) { }
00095
00097 TpetraCrsGraph(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, ProfileType pftype=DynamicProfile, const RCP< ParameterList > ¶ms=null)
00098 : graph_(Teuchos::rcp(new Tpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node >(toTpetra(rowMap), NumEntriesPerRowToAlloc, toTpetra(pftype), params))) { }
00099
00101 TpetraCrsGraph(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, size_t maxNumEntriesPerRow, ProfileType pftype=DynamicProfile, const RCP< ParameterList > ¶ms=null)
00102 : graph_(Teuchos::rcp(new Tpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node >(toTpetra(rowMap), toTpetra(colMap), maxNumEntriesPerRow, toTpetra(pftype), params))) { }
00103
00105 TpetraCrsGraph(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, ProfileType pftype=DynamicProfile, const RCP< ParameterList > ¶ms=null)
00106 : graph_(Teuchos::rcp(new Tpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node >(toTpetra(rowMap), toTpetra(colMap), NumEntriesPerRowToAlloc, toTpetra(pftype), params))) { }
00107
00109 virtual ~TpetraCrsGraph() { }
00110
00112
00114
00115
00117 void insertGlobalIndices(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &indices) { XPETRA_MONITOR("TpetraCrsGraph::insertGlobalIndices"); graph_->insertGlobalIndices(globalRow, indices); }
00118
00120 void insertLocalIndices(const LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &indices) { XPETRA_MONITOR("TpetraCrsGraph::insertLocalIndices"); graph_->insertLocalIndices(localRow, indices); }
00121
00123 void removeLocalIndices(LocalOrdinal localRow) { XPETRA_MONITOR("TpetraCrsGraph::removeLocalIndices"); graph_->removeLocalIndices(localRow); }
00124
00126
00128
00129
00131 void fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > ¶ms=null) { XPETRA_MONITOR("TpetraCrsGraph::fillComplete"); graph_->fillComplete(toTpetra(domainMap), toTpetra(rangeMap), params); }
00132
00134 void fillComplete(const RCP< ParameterList > ¶ms=null) { XPETRA_MONITOR("TpetraCrsGraph::fillComplete"); graph_->fillComplete(params); }
00135
00137
00139
00140
00142 RCP< const Comm< int > > getComm() const { XPETRA_MONITOR("TpetraCrsGraph::getComm"); return graph_->getComm(); }
00143
00145 RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const { XPETRA_MONITOR("TpetraCrsGraph::getRowMap"); return toXpetra(graph_->getRowMap()); }
00146
00148 RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const { XPETRA_MONITOR("TpetraCrsGraph::getColMap"); return toXpetra(graph_->getColMap()); }
00149
00151 RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const { XPETRA_MONITOR("TpetraCrsGraph::getDomainMap"); return toXpetra(graph_->getDomainMap()); }
00152
00154 RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const { XPETRA_MONITOR("TpetraCrsGraph::getRangeMap"); return toXpetra(graph_->getRangeMap()); }
00155
00157 RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > getImporter() const { XPETRA_MONITOR("TpetraCrsGraph::getImporter"); return toXpetra(graph_->getImporter()); }
00158
00160 RCP< const Export< LocalOrdinal, GlobalOrdinal, Node > > getExporter() const { XPETRA_MONITOR("TpetraCrsGraph::getExporter"); return toXpetra(graph_->getExporter()); }
00161
00163 global_size_t getGlobalNumRows() const { XPETRA_MONITOR("TpetraCrsGraph::getGlobalNumRows"); return graph_->getGlobalNumRows(); }
00164
00166 global_size_t getGlobalNumCols() const { XPETRA_MONITOR("TpetraCrsGraph::getGlobalNumCols"); return graph_->getGlobalNumCols(); }
00167
00169 size_t getNodeNumRows() const { XPETRA_MONITOR("TpetraCrsGraph::getNodeNumRows"); return graph_->getNodeNumRows(); }
00170
00172 size_t getNodeNumCols() const { XPETRA_MONITOR("TpetraCrsGraph::getNodeNumCols"); return graph_->getNodeNumCols(); }
00173
00175 GlobalOrdinal getIndexBase() const { XPETRA_MONITOR("TpetraCrsGraph::getIndexBase"); return graph_->getIndexBase(); }
00176
00178 global_size_t getGlobalNumEntries() const { XPETRA_MONITOR("TpetraCrsGraph::getGlobalNumEntries"); return graph_->getGlobalNumEntries(); }
00179
00181 size_t getNodeNumEntries() const { XPETRA_MONITOR("TpetraCrsGraph::getNodeNumEntries"); return graph_->getNodeNumEntries(); }
00182
00184 size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const { XPETRA_MONITOR("TpetraCrsGraph::getNumEntriesInGlobalRow"); return graph_->getNumEntriesInGlobalRow(globalRow); }
00185
00187 size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const { XPETRA_MONITOR("TpetraCrsGraph::getNumEntriesInLocalRow"); return graph_->getNumEntriesInLocalRow(localRow); }
00188
00190 size_t getNumAllocatedEntriesInGlobalRow(GlobalOrdinal globalRow) const { XPETRA_MONITOR("TpetraCrsGraph::getNumAllocatedEntriesInGlobalRow"); return graph_->getNumAllocatedEntriesInGlobalRow(globalRow); }
00191
00193 size_t getNumAllocatedEntriesInLocalRow(LocalOrdinal localRow) const { XPETRA_MONITOR("TpetraCrsGraph::getNumAllocatedEntriesInLocalRow"); return graph_->getNumAllocatedEntriesInLocalRow(localRow); }
00194
00196 global_size_t getGlobalNumDiags() const { XPETRA_MONITOR("TpetraCrsGraph::getGlobalNumDiags"); return graph_->getGlobalNumDiags(); }
00197
00199 size_t getNodeNumDiags() const { XPETRA_MONITOR("TpetraCrsGraph::getNodeNumDiags"); return graph_->getNodeNumDiags(); }
00200
00202 size_t getGlobalMaxNumRowEntries() const { XPETRA_MONITOR("TpetraCrsGraph::getGlobalMaxNumRowEntries"); return graph_->getGlobalMaxNumRowEntries(); }
00203
00205 size_t getNodeMaxNumRowEntries() const { XPETRA_MONITOR("TpetraCrsGraph::getNodeMaxNumRowEntries"); return graph_->getNodeMaxNumRowEntries(); }
00206
00208 bool hasColMap() const { XPETRA_MONITOR("TpetraCrsGraph::hasColMap"); return graph_->hasColMap(); }
00209
00211 bool isLowerTriangular() const { XPETRA_MONITOR("TpetraCrsGraph::isLowerTriangular"); return graph_->isLowerTriangular(); }
00212
00214 bool isUpperTriangular() const { XPETRA_MONITOR("TpetraCrsGraph::isUpperTriangular"); return graph_->isUpperTriangular(); }
00215
00217 bool isLocallyIndexed() const { XPETRA_MONITOR("TpetraCrsGraph::isLocallyIndexed"); return graph_->isLocallyIndexed(); }
00218
00220 bool isGloballyIndexed() const { XPETRA_MONITOR("TpetraCrsGraph::isGloballyIndexed"); return graph_->isGloballyIndexed(); }
00221
00223 bool isFillComplete() const { XPETRA_MONITOR("TpetraCrsGraph::isFillComplete"); return graph_->isFillComplete(); }
00224
00226 bool isStorageOptimized() const { XPETRA_MONITOR("TpetraCrsGraph::isStorageOptimized"); return graph_->isStorageOptimized(); }
00227
00229 void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView< const GlobalOrdinal > &Indices) const { XPETRA_MONITOR("TpetraCrsGraph::getGlobalRowView"); graph_->getGlobalRowView(GlobalRow, Indices); }
00230
00232 void getLocalRowView(LocalOrdinal LocalRow, ArrayView< const LocalOrdinal > &indices) const { XPETRA_MONITOR("TpetraCrsGraph::getLocalRowView"); graph_->getLocalRowView(LocalRow, indices); }
00233
00235
00237
00238
00240 std::string description() const { XPETRA_MONITOR("TpetraCrsGraph::description"); return graph_->description(); }
00241
00243 void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const { XPETRA_MONITOR("TpetraCrsGraph::describe"); graph_->describe(out, verbLevel); }
00244
00246
00248
00249
00251 ArrayRCP< const size_t > getNodeRowPtrs() const { XPETRA_MONITOR("TpetraCrsGraph::getNodeRowPtrs"); return graph_->getNodeRowPtrs(); }
00252
00254
00256
00257
00259 Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const { XPETRA_MONITOR("TpetraCrsGraph::getMap"); return rcp( new TpetraMap< LocalOrdinal, GlobalOrdinal, Node >(graph_->getMap()) ); }
00260
00262 void doImport(const DistObject<GlobalOrdinal, LocalOrdinal, GlobalOrdinal, Node> &source,
00263 const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM) {
00264 XPETRA_MONITOR("TpetraCrsGraph::doImport");
00265
00266 XPETRA_DYNAMIC_CAST(const TpetraCrsGraphClass, source, tSource, "Xpetra::TpetraCrsGraph::doImport only accept Xpetra::TpetraCrsGraph as input arguments.");
00267 RCP< const Tpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > v = tSource.getTpetra_CrsGraph();
00268
00269
00270 graph_->doImport(*v, toTpetra(importer), toTpetra(CM));
00271 }
00272
00274 void doExport(const DistObject<GlobalOrdinal, LocalOrdinal, GlobalOrdinal, Node> &dest,
00275 const Import< LocalOrdinal, GlobalOrdinal, Node >& importer, CombineMode CM) {
00276 XPETRA_MONITOR("TpetraCrsGraph::doExport");
00277
00278 XPETRA_DYNAMIC_CAST(const TpetraCrsGraphClass, dest, tDest, "Xpetra::TpetraCrsGraph::doImport only accept Xpetra::TpetraCrsGraph as input arguments.");
00279 RCP< const Tpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > v = tDest.getTpetra_CrsGraph();
00280 graph_->doExport(*v, toTpetra(importer), toTpetra(CM));
00281
00282 }
00283
00285 void doImport(const DistObject<GlobalOrdinal, LocalOrdinal, GlobalOrdinal, Node> &source,
00286 const Export< LocalOrdinal, GlobalOrdinal, Node >& exporter, CombineMode CM) {
00287 XPETRA_MONITOR("TpetraCrsGraph::doImport");
00288
00289 XPETRA_DYNAMIC_CAST(const TpetraCrsGraphClass, source, tSource, "Xpetra::TpetraCrsGraph::doImport only accept Xpetra::TpetraCrsGraph as input arguments.");
00290 RCP< const Tpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > v = tSource.getTpetra_CrsGraph();
00291
00292 graph_->doImport(*v, toTpetra(exporter), toTpetra(CM));
00293
00294 }
00295
00297 void doExport(const DistObject<GlobalOrdinal, LocalOrdinal, GlobalOrdinal, Node> &dest,
00298 const Export< LocalOrdinal, GlobalOrdinal, Node >& exporter, CombineMode CM) {
00299 XPETRA_MONITOR("TpetraCrsGraph::doExport");
00300
00301 XPETRA_DYNAMIC_CAST(const TpetraCrsGraphClass, dest, tDest, "Xpetra::TpetraCrsGraph::doImport only accept Xpetra::TpetraCrsGraph as input arguments.");
00302 RCP< const Tpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > v = tDest.getTpetra_CrsGraph();
00303
00304 graph_->doExport(*v, toTpetra(exporter), toTpetra(CM));
00305
00306 }
00307
00308
00309
00311
00312
00314 TpetraCrsGraph(const Teuchos::RCP<Tpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > &graph) : graph_(graph) { }
00315
00317 RCP< const Tpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > getTpetra_CrsGraph() const { return graph_; }
00318
00320
00321 private:
00322 RCP< Tpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > graph_;
00323 };
00324
00325
00326 template <class LocalOrdinal, class GlobalOrdinal, class Node>
00327 RCP<const CrsGraph<LocalOrdinal, GlobalOrdinal, Node> >
00328 toXpetra (RCP<const Tpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > graph)
00329 {
00330
00331
00332 if (graph.is_null ()) {
00333 return Teuchos::null;
00334 }
00335 RCP<Tpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > tGraph =
00336 Teuchos::rcp_const_cast<Tpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > (graph);
00337 return rcp (new Xpetra::TpetraCrsGraph<LocalOrdinal, GlobalOrdinal, Node> (tGraph));
00338 }
00339
00340 template <class LocalOrdinal, class GlobalOrdinal, class Node>
00341 RCP< const Tpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > >
00342 toTpetra (const RCP<const CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > &graph)
00343 {
00344 typedef TpetraCrsGraph<LocalOrdinal, GlobalOrdinal, Node> TpetraCrsGraphClass;
00345 XPETRA_RCP_DYNAMIC_CAST(const TpetraCrsGraphClass, graph, tpetraCrsGraph, "toTpetra");
00346 return tpetraCrsGraph->getTpetra_CrsGraph ();
00347 }
00348
00349 }
00350
00351 #define XPETRA_TPETRACRSGRAPH_SHORT
00352 #endif // XPETRA_TPETRACRSGRAPH_HPP