Ifpack Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
Ifpack_IlukGraph.h
Go to the documentation of this file.
1/*@HEADER
2// ***********************************************************************
3//
4// Ifpack: Object-Oriented Algebraic Preconditioner Package
5// Copyright (2002) Sandia Corporation
6//
7// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8// license for use of this work by or on behalf of the U.S. Government.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38//
39// ***********************************************************************
40//@HEADER
41*/
42
43#ifndef _IFPACK_ILUK_GRAPH_H_
44#define _IFPACK_ILUK_GRAPH_H_
45
46#include "Ifpack_ConfigDefs.h"
47#include "Epetra_Object.h"
48#include "Epetra_CrsGraph.h"
49#include "Epetra_Import.h"
50
51#include "Teuchos_RefCountPtr.hpp"
52
53namespace Teuchos {
54 class ParameterList;
55}
56
58
78
79 // Give ostream << function some access to private and protected data/functions.
80
81 friend std::ostream& operator << (std::ostream& os, const Ifpack_IlukGraph& A);
82
83 public:
84
86
99 Ifpack_IlukGraph(const Epetra_CrsGraph & Graph_in, int LevelFill_in, int LevelOverlap_in);
100
102 Ifpack_IlukGraph(const Ifpack_IlukGraph & Graph_in);
103
105 virtual ~Ifpack_IlukGraph();
106
108 /* This method is only available if the Teuchos package is enabled.
109 This method recogizes two parameter names: Level_fill and Level_overlap.
110 Both are case insensitive, and in both cases the ParameterEntry must
111 have type int.
112 */
113 int SetParameters(const Teuchos::ParameterList& parameterlist,
114 bool cerr_warning_if_unused=false);
115
117 /*
118 \return Integer error code, set to 0 if successful.
119
120 */
121 virtual int ConstructFilledGraph();
122
124 /*
125 \return Integer error code, set to 0 if successful.
126
127 */
128 virtual int ConstructOverlapGraph();
129
131 virtual int LevelFill() const {return(LevelFill_);};
132
134 virtual int LevelOverlap() const {return(LevelOverlap_);};
135
136#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
138 int NumGlobalBlockRows() const {
140 return (int)(NumGlobalBlockRows_);
141 else
142 throw "Ifpack_IlukGraph::NumGlobalBlockRows: GlobalIndices not int.";
143 }
144
146 int NumGlobalBlockCols() const {
148 return (int)(NumGlobalBlockCols_);
149 else
150 throw "Ifpack_IlukGraph::NumGlobalBlockCols: GlobalIndices not int.";
151 }
152
154 int NumGlobalRows() const {
156 return (int)(NumGlobalRows_);
157 else
158 throw "Ifpack_IlukGraph::NumGlobalRows: GlobalIndices not int.";
159 }
160
162 int NumGlobalCols() const {
164 return (int)(NumGlobalCols_);
165 else
166 throw "Ifpack_IlukGraph::NumGlobalCols: GlobalIndices not int.";
167 }
169 int NumGlobalNonzeros() const {
171 return (int)(NumGlobalNonzeros_);
172 else
173 throw "Ifpack_IlukGraph::NumGlobalNonzeros: GlobalIndices not int.";
174 }
175
177 virtual int NumGlobalBlockDiagonals() const {
179 return (int)(NumGlobalBlockDiagonals_);
180 else
181 throw "Ifpack_IlukGraph::NumGlobalBlockDiagonals: GlobalIndices not int.";
182 }
183#endif
184
186 long long NumGlobalBlockRows64() const {return(NumGlobalBlockRows_);};
187
189 long long NumGlobalBlockCols64() const {return(NumGlobalBlockCols_);};
190
192 long long NumGlobalRows64() const {return(NumGlobalRows_);};
193
195 long long NumGlobalCols64() const {return(NumGlobalCols_);};
197 long long NumGlobalNonzeros64() const {return(NumGlobalNonzeros_);};
198
200 virtual long long NumGlobalBlockDiagonals64() const {return(NumGlobalBlockDiagonals_);};
201
203 int NumMyBlockRows() const {return(NumMyBlockRows_);};
204
206 int NumMyBlockCols() const {return(NumMyBlockCols_);};
207
208
210 int NumMyRows() const {return(NumMyRows_);};
211
213 int NumMyCols() const {return(NumMyCols_);};
214
216 int NumMyNonzeros() const {return(NumMyNonzeros_);};
217
219 virtual int NumMyBlockDiagonals() const {return(NumMyBlockDiagonals_);};
220
222#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
223 int IndexBase() const {
225 return (int) IndexBase64();
226 throw "Ifpack_IlukGraph::IndexBase: GlobalIndices not int.";
227 }
228#endif
229 long long IndexBase64() const {return(IndexBase_);};
230
232 virtual Epetra_CrsGraph & L_Graph() {return(*L_Graph_);};
233
235 virtual Epetra_CrsGraph & U_Graph() {return(*U_Graph_);};
236
238 virtual Epetra_CrsGraph & L_Graph() const {return(*L_Graph_);};
239
241 virtual Epetra_CrsGraph & U_Graph() const {return(*U_Graph_);};
242
244 virtual Epetra_Import * OverlapImporter() const {return(&*OverlapImporter_);};
245
247 virtual Epetra_CrsGraph * OverlapGraph() const {return(&*OverlapGraph_);};
248
250 virtual const Epetra_BlockMap & DomainMap() const {return(DomainMap_);};
251
253 virtual const Epetra_BlockMap & RangeMap() const{return(RangeMap_);};
254
256 virtual const Epetra_Comm & Comm() const{return(Comm_);};
257
258 private:
259
260
265 Teuchos::RefCountPtr<Epetra_CrsGraph> OverlapGraph_;
266 Teuchos::RefCountPtr<Epetra_BlockMap> OverlapRowMap_;
267 Teuchos::RefCountPtr<Epetra_Import> OverlapImporter_;
270 Teuchos::RefCountPtr<Epetra_CrsGraph> L_Graph_;
271 Teuchos::RefCountPtr<Epetra_CrsGraph> U_Graph_;
272 long long IndexBase_;
273 long long NumGlobalRows_;
274 long long NumGlobalCols_;
287
288
289 };
290
292std::ostream& operator << (std::ostream& os, const Ifpack_IlukGraph& A);
293
294#endif /* _IFPACK_ILUK_GRAPH_H_ */
std::ostream & operator<<(std::ostream &os, const Ifpack_IlukGraph &A)
<< operator will work for Ifpack_IlukGraph.
bool GlobalIndicesInt() const
const Epetra_BlockMap & RowMap() const
Ifpack_IlukGraph: A class for constructing level filled graphs for use with ILU(k) class precondition...
virtual int LevelFill() const
Returns the level of fill used to construct this graph.
int NumGlobalNonzeros() const
Returns the number of nonzero entries in the global graph.
long long NumGlobalCols_
int SetParameters(const Teuchos::ParameterList &parameterlist, bool cerr_warning_if_unused=false)
Set parameters using Teuchos::ParameterList object.
long long NumGlobalNonzeros64() const
Returns the number of nonzero entries in the global graph.
int NumGlobalBlockRows() const
Returns the number of global matrix rows.
virtual ~Ifpack_IlukGraph()
Ifpack_IlukGraph Destructor.
Teuchos::RefCountPtr< Epetra_BlockMap > OverlapRowMap_
virtual Epetra_CrsGraph & L_Graph() const
Returns the graph of lower triangle of the ILU(k) graph as a Epetra_CrsGraph.
int NumGlobalRows() const
Returns the number of global matrix rows.
long long NumGlobalNonzeros_
virtual Epetra_CrsGraph & U_Graph() const
Returns the graph of upper triangle of the ILU(k) graph as a Epetra_CrsGraph.
long long NumGlobalRows_
virtual int NumMyBlockDiagonals() const
Returns the number of diagonal entries found in the local input graph.
int IndexBase() const
Returns the index base for row and column indices for this graph.
const Epetra_BlockMap & DomainMap_
virtual const Epetra_BlockMap & DomainMap() const
Returns the Epetra_BlockMap object associated with the domain of this matrix operator.
virtual int ConstructOverlapGraph()
Does the actual construction of the overlap matrix graph.
long long NumGlobalEntries_
const Epetra_Comm & Comm_
const Epetra_CrsGraph & Graph_
virtual Epetra_CrsGraph & U_Graph()
Returns the graph of upper triangle of the ILU(k) graph as a Epetra_CrsGraph.
long long NumGlobalBlockRows64() const
Returns the number of global matrix rows.
int NumMyBlockRows() const
Returns the number of local matrix rows.
Teuchos::RefCountPtr< Epetra_Import > OverlapImporter_
Teuchos::RefCountPtr< Epetra_CrsGraph > OverlapGraph_
long long IndexBase64() const
int NumMyBlockCols() const
Returns the number of local matrix columns.
Teuchos::RefCountPtr< Epetra_CrsGraph > L_Graph_
int NumGlobalCols() const
Returns the number of global matrix columns.
int NumMyCols() const
Returns the number of local matrix columns.
long long NumGlobalBlockCols64() const
Returns the number of global matrix columns.
virtual long long NumGlobalBlockDiagonals64() const
Returns the number of diagonal entries found in the global input graph.
virtual int NumGlobalBlockDiagonals() const
Returns the number of diagonal entries found in the global input graph.
Teuchos::RefCountPtr< Epetra_CrsGraph > U_Graph_
virtual const Epetra_Comm & Comm() const
Returns the Epetra_BlockMap object associated with the range of this matrix operator.
const Epetra_BlockMap & RangeMap_
long long NumGlobalRows64() const
Returns the number of global matrix rows.
int NumMyRows() const
Returns the number of local matrix rows.
virtual Epetra_CrsGraph * OverlapGraph() const
Returns the the overlapped graph.
long long NumGlobalBlockRows_
virtual const Epetra_BlockMap & RangeMap() const
Returns the Epetra_BlockMap object associated with the range of this matrix operator.
int NumMyNonzeros() const
Returns the number of nonzero entries in the local graph.
long long NumGlobalCols64() const
Returns the number of global matrix columns.
friend std::ostream & operator<<(std::ostream &os, const Ifpack_IlukGraph &A)
<< operator will work for Ifpack_IlukGraph.
virtual int LevelOverlap() const
Returns the level of overlap used to construct this graph.
virtual Epetra_CrsGraph & L_Graph()
Returns the graph of lower triangle of the ILU(k) graph as a Epetra_CrsGraph.
virtual int ConstructFilledGraph()
Does the actual construction of the graph.
long long NumGlobalBlockDiagonals_
long long NumGlobalBlockCols_
int NumGlobalBlockCols() const
Returns the number of global matrix columns.
virtual Epetra_Import * OverlapImporter() const
Returns the importer used to create the overlapped graph.