10 #ifndef IFPACK2_CREATEOVERLAPGRAPH_HPP
11 #define IFPACK2_CREATEOVERLAPGRAPH_HPP
13 #include "Ifpack2_ConfigDefs.hpp"
14 #include "Tpetra_CrsGraph.hpp"
15 #include "Tpetra_CrsMatrix.hpp"
16 #include "Tpetra_Import.hpp"
17 #include "Teuchos_RefCountPtr.hpp"
37 template <
class GraphType>
40 const int overlapLevel) {
43 typedef Tpetra::Map<
typename GraphType::local_ordinal_type,
44 typename GraphType::global_ordinal_type,
45 typename GraphType::node_type>
47 typedef Tpetra::Import<
typename GraphType::local_ordinal_type,
48 typename GraphType::global_ordinal_type,
49 typename GraphType::node_type>
52 overlapLevel < 0, std::invalid_argument,
53 "Ifpack2::createOverlapGraph: overlapLevel must be >= 0, "
54 "but you specified overlapLevel = "
55 << overlapLevel <<
".");
57 const int numProcs = inputGraph->getMap()->getComm()->getSize();
58 if (overlapLevel == 0 || numProcs < 2) {
62 RCP<const map_type> overlapRowMap = inputGraph->getRowMap();
63 RCP<const map_type> domainMap = inputGraph->getDomainMap();
64 RCP<const map_type> rangeMap = inputGraph->getRangeMap();
66 RCP<GraphType> overlapGraph;
67 RCP<const GraphType> oldGraph;
68 RCP<const map_type> oldRowMap;
69 for (
int level = 0; level < overlapLevel; ++level) {
70 oldGraph = overlapGraph;
71 oldRowMap = overlapRowMap;
73 RCP<const import_type> overlapImporter;
75 overlapImporter = inputGraph->getImporter();
77 overlapImporter = oldGraph->getImporter();
80 overlapRowMap = overlapImporter->getTargetMap();
81 if (level < overlapLevel - 1) {
82 overlapGraph =
rcp(
new GraphType(overlapRowMap, 0));
86 overlapGraph =
rcp(
new GraphType(overlapRowMap, overlapRowMap, 0));
89 overlapGraph->doImport(*inputGraph, *overlapImporter, Tpetra::INSERT);
90 overlapGraph->fillComplete(domainMap, rangeMap);
105 template <
class MatrixType>
108 const int overlapLevel) {
111 typedef typename MatrixType::map_type map_type;
112 typedef Tpetra::Import<
typename MatrixType::local_ordinal_type,
113 typename MatrixType::global_ordinal_type,
114 typename MatrixType::node_type>
118 overlapLevel < 0, std::invalid_argument,
119 "Ifpack2::createOverlapMatrix: overlapLevel must be >= 0, "
120 "but you specified overlapLevel = "
121 << overlapLevel <<
".");
123 const int numProcs = inputMatrix->getMap()->getComm()->getSize();
124 if (overlapLevel == 0 || numProcs < 2) {
128 RCP<const map_type> overlapRowMap = inputMatrix->getRowMap();
129 RCP<const map_type> domainMap = inputMatrix->getDomainMap();
130 RCP<const map_type> rangeMap = inputMatrix->getRangeMap();
132 RCP<MatrixType> overlapMatrix;
133 RCP<const MatrixType> oldMatrix;
134 RCP<const map_type> oldRowMap;
135 for (
int level = 0; level < overlapLevel; ++level) {
136 oldMatrix = overlapMatrix;
137 oldRowMap = overlapRowMap;
139 RCP<const import_type> overlapImporter;
141 overlapImporter = inputMatrix->getGraph()->getImporter();
143 overlapImporter = oldMatrix->getGraph()->getImporter();
146 overlapRowMap = overlapImporter->getTargetMap();
147 if (level < overlapLevel - 1) {
148 overlapMatrix =
rcp(
new MatrixType(overlapRowMap, 0));
152 overlapMatrix =
rcp(
new MatrixType(overlapRowMap, overlapRowMap, 0));
155 overlapMatrix->doImport(*inputMatrix, *overlapImporter, Tpetra::INSERT);
156 overlapMatrix->fillComplete(domainMap, rangeMap);
159 return overlapMatrix;
164 #endif // IFPACK2_CREATEOVERLAPGRAPH_HPP
Teuchos::RCP< const GraphType > createOverlapGraph(const Teuchos::RCP< const GraphType > &inputGraph, const int overlapLevel)
Construct an overlapped graph for use with Ifpack2 preconditioners.
Definition: Ifpack2_CreateOverlapGraph.hpp:39
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Teuchos::RCP< const MatrixType > createOverlapMatrix(const Teuchos::RCP< const MatrixType > &inputMatrix, const int overlapLevel)
Construct an overlapped matrix for use with Ifpack2 preconditioners.
Definition: Ifpack2_CreateOverlapGraph.hpp:107