49 #ifndef IFPACK2_ILUK_GRAPH_HPP
50 #define IFPACK2_ILUK_GRAPH_HPP
55 #include <Ifpack2_ConfigDefs.hpp>
57 #include <Teuchos_CommHelpers.hpp>
58 #include <Tpetra_CrsGraph.hpp>
59 #include <Tpetra_Import.hpp>
60 #include <Ifpack2_CreateOverlapGraph.hpp>
61 #include <Ifpack2_Parameters.hpp>
96 template<
class GraphType>
99 typedef typename GraphType::local_ordinal_type local_ordinal_type;
100 typedef typename GraphType::global_ordinal_type global_ordinal_type;
101 typedef typename GraphType::node_type node_type;
104 typedef Tpetra::RowGraph<local_ordinal_type,
108 typedef Tpetra::CrsGraph<local_ordinal_type,
124 const int levelOverlap);
167 return OverlapGraph_;
174 typedef typename GraphType::map_type map_type;
197 void constructOverlapGraph();
205 size_t NumMyDiagonals_;
206 size_t NumGlobalDiagonals_;
210 template<
class GraphType>
214 const int levelOverlap)
216 LevelFill_ (levelFill),
217 LevelOverlap_ (levelOverlap),
219 NumGlobalDiagonals_ (0)
223 template<
class GraphType>
228 template<
class GraphType>
232 getParameter (parameterlist,
"iluk level-of-fill", LevelFill_);
233 getParameter (parameterlist,
"iluk level-of-overlap", LevelOverlap_);
237 template<
class GraphType>
242 if (OverlapGraph_ == Teuchos::null) {
243 OverlapGraph_ = createOverlapGraph<GraphType> (Graph_, LevelOverlap_);
248 template<
class GraphType>
258 size_t NumIn, NumL, NumU;
261 constructOverlapGraph();
264 const int MaxNumIndices = OverlapGraph_->getNodeMaxNumRowEntries ();
268 const int NumMyRows = OverlapGraph_->getRowMap ()->getNodeNumElements ();
271 const int MaxNumEntriesPerRow = (LevelFill_ == 0)
273 : MaxNumIndices + 5*LevelFill_;
275 OverlapGraph_->getRowMap (), MaxNumEntriesPerRow));
277 OverlapGraph_->getRowMap (), MaxNumEntriesPerRow));
279 Array<local_ordinal_type> L (MaxNumIndices);
280 Array<local_ordinal_type> U (MaxNumIndices);
286 for (
int i = 0; i< NumMyRows; ++i) {
287 ArrayView<const local_ordinal_type> my_indices;
288 OverlapGraph_->getLocalRowView (i, my_indices);
295 NumIn = my_indices.size();
297 for (
size_t j = 0; j < NumIn; ++j) {
298 const local_ordinal_type k = my_indices[j];
322 ArrayView<const local_ordinal_type> L_view = L.view (0, NumL);
323 L_Graph_->insertLocalIndices (i, L_view);
326 ArrayView<const local_ordinal_type> U_view = U.view (0, NumU);
327 U_Graph_->insertLocalIndices (i, U_view);
331 if (LevelFill_ > 0) {
333 RCP<const map_type> L_DomainMap = OverlapGraph_->getRowMap ();
334 RCP<const map_type> L_RangeMap = Graph_->getRangeMap ();
335 RCP<const map_type> U_DomainMap = Graph_->getDomainMap ();
336 RCP<const map_type> U_RangeMap = OverlapGraph_->getRowMap ();
337 RCP<Teuchos::ParameterList> params = Teuchos::parameterList ();
338 params->set (
"Optimize Storage",
false);
339 L_Graph_->fillComplete (L_DomainMap, L_RangeMap, params);
340 U_Graph_->fillComplete (U_DomainMap, U_RangeMap, params);
341 L_Graph_->resumeFill ();
342 U_Graph_->resumeFill ();
348 int MaxRC = NumMyRows;
349 std::vector<std::vector<int> > Levels(MaxRC);
350 std::vector<int> LinkList(MaxRC);
351 std::vector<int> CurrentLevel(MaxRC);
352 Array<local_ordinal_type> CurrentRow (MaxRC + 1);
353 std::vector<int> LevelsRowU(MaxRC);
355 for (
int i = 0; i < NumMyRows; ++i) {
360 size_t LenL = L_Graph_->getNumEntriesInLocalRow(i);
361 size_t LenU = U_Graph_->getNumEntriesInLocalRow(i);
362 size_t Len = LenL + LenU + 1;
364 CurrentRow.resize(Len);
366 L_Graph_->getLocalRowCopy(i, CurrentRow(), LenL);
367 CurrentRow[LenL] = i;
369 ArrayView<local_ordinal_type> URowView = CurrentRow.view (LenL+1, LenU);
371 U_Graph_->getLocalRowCopy (i, URowView, LenU);
376 for (
size_t j=0; j<Len-1; j++) {
377 LinkList[CurrentRow[j]] = CurrentRow[j+1];
378 CurrentLevel[CurrentRow[j]] = 0;
381 LinkList[CurrentRow[Len-1]] = NumMyRows;
382 CurrentLevel[CurrentRow[Len-1]] = 0;
386 First = CurrentRow[0];
389 int PrevInList = Next;
390 int NextInList = LinkList[Next];
393 ArrayView<const local_ordinal_type> IndicesU;
394 U_Graph_->getLocalRowView (RowU, IndicesU);
396 int LengthRowU = IndicesU.size ();
402 for (ii = 0; ii < LengthRowU; ) {
403 int CurInList = IndicesU[ii];
404 if (CurInList < NextInList) {
406 int NewLevel = CurrentLevel[RowU] + Levels[RowU][ii+1] + 1;
407 if (NewLevel <= LevelFill_) {
408 LinkList[PrevInList] = CurInList;
409 LinkList[CurInList] = NextInList;
410 PrevInList = CurInList;
411 CurrentLevel[CurInList] = NewLevel;
415 else if (CurInList == NextInList) {
416 PrevInList = NextInList;
417 NextInList = LinkList[PrevInList];
418 int NewLevel = CurrentLevel[RowU] + Levels[RowU][ii+1] + 1;
419 CurrentLevel[CurInList] = std::min (CurrentLevel[CurInList], NewLevel);
423 PrevInList = NextInList;
424 NextInList = LinkList[PrevInList];
427 Next = LinkList[Next];
432 CurrentRow.resize (0);
439 CurrentRow.push_back (Next);
440 Next = LinkList[Next];
446 L_Graph_->removeLocalIndices (i);
447 if (CurrentRow.size() > 0) {
448 L_Graph_->insertLocalIndices (i, CurrentRow ());
454 Next != i, std::runtime_error,
455 "Ifpack2::IlukGraph::initialize: FATAL: U has zero diagonal")
457 LevelsRowU[0] = CurrentLevel[Next];
458 Next = LinkList[Next];
462 CurrentRow.resize (0);
465 while (Next < NumMyRows) {
466 LevelsRowU[LenU+1] = CurrentLevel[Next];
467 CurrentRow.push_back (Next);
469 Next = LinkList[Next];
476 U_Graph_->removeLocalIndices (i);
478 U_Graph_->insertLocalIndices (i, CurrentRow ());
482 Levels[i] = std::vector<int> (LenU+1);
483 for (
size_t jj=0; jj<LenU+1; jj++) {
484 Levels[i][jj] = LevelsRowU[jj];
490 RCP<const map_type> L_DomainMap = OverlapGraph_->getRowMap ();
491 RCP<const map_type> L_RangeMap = Graph_->getRangeMap ();
492 RCP<const map_type> U_DomainMap = Graph_->getDomainMap ();
493 RCP<const map_type> U_RangeMap = OverlapGraph_->getRowMap ();
494 L_Graph_->fillComplete (L_DomainMap, L_RangeMap);
495 U_Graph_->fillComplete (U_DomainMap, U_RangeMap);
497 reduceAll<int, size_t> (* (L_DomainMap->getComm ()), REDUCE_SUM, 1,
498 &NumMyDiagonals_, &NumGlobalDiagonals_);
void initialize()
Set up the graph structure of the L and U factors.
Definition: Ifpack2_IlukGraph.hpp:249
void setParameters(const Teuchos::ParameterList ¶meterlist)
Set parameters.
Definition: Ifpack2_IlukGraph.hpp:230
Tpetra::CrsGraph< local_ordinal_type, global_ordinal_type, node_type > crs_graph_type
Tpetra::CrsGraph specialization used by this class.
Definition: Ifpack2_IlukGraph.hpp:110
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
void getParameter(const Teuchos::ParameterList ¶ms, const std::string &name, T &value)
Set a value from a ParameterList if a parameter with the specified name exists.
Definition: Ifpack2_Parameters.hpp:59
Tpetra::RowGraph< local_ordinal_type, global_ordinal_type, node_type > row_graph_type
Tpetra::RowGraph specialization used by this class.
Definition: Ifpack2_IlukGraph.hpp:106
virtual ~IlukGraph()
IlukGraph Destructor.
Definition: Ifpack2_IlukGraph.hpp:224
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Teuchos::RCP< crs_graph_type > getL_Graph() const
Returns the graph of lower triangle of the ILU(k) graph as a Tpetra::CrsGraph.
Definition: Ifpack2_IlukGraph.hpp:156
IlukGraph(const Teuchos::RCP< const GraphType > &G, const int levelFill, const int levelOverlap)
Constructor.
Definition: Ifpack2_IlukGraph.hpp:212
TEUCHOS_DEPRECATED void reduceAll(const Comm< Ordinal > &comm, const EReductionType reductType, const Packet &send, Packet *globalReduct)
int getLevelFill() const
The level of fill used to construct this graph.
Definition: Ifpack2_IlukGraph.hpp:150
int getLevelOverlap() const
The level of overlap used to construct this graph.
Definition: Ifpack2_IlukGraph.hpp:153
Construct a level filled graph for use in computing an ILU(k) incomplete factorization.
Definition: Ifpack2_IlukGraph.hpp:97
size_t getNumGlobalDiagonals() const
Returns the global number of diagonals in the ILU(k) graph.
Definition: Ifpack2_IlukGraph.hpp:171
Teuchos::RCP< const crs_graph_type > getOverlapGraph() const
Returns the the overlapped graph.
Definition: Ifpack2_IlukGraph.hpp:166
Teuchos::RCP< crs_graph_type > getU_Graph() const
Returns the graph of upper triangle of the ILU(k) graph as a Tpetra::CrsGraph.
Definition: Ifpack2_IlukGraph.hpp:161