Ifpack Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Ifpack_DenseContainer.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_DENSECONTAINER_H
44 #define IFPACK_DENSECONTAINER_H
45 
46 #include "Ifpack_ConfigDefs.h"
47 #include "Ifpack_Container.h"
48 #include "Epetra_SerialDenseMatrix.h"
49 #include "Epetra_SerialDenseSolver.h"
50 #include "Epetra_IntSerialDenseVector.h"
51 class Epetra_RowMatrix;
52 
54 
118 
119 public:
120 
122 
124  Ifpack_DenseContainer(const int NumRows_in, const int NumVectors_in = 1) :
125  NumRows_(NumRows_in),
126  NumVectors_(NumVectors_in),
130  ComputeFlops_(0.0),
131  ApplyFlops_(0.0),
132  ApplyInverseFlops_(0.0)
133  {}
134 
137  NumRows_(rhs.NumRows()),
138  NumVectors_(rhs.NumVectors()),
141  IsComputed_(rhs.IsComputed())
142  {
143  Matrix_ = rhs.Matrix();
146  LHS_ = rhs.LHS();
147  RHS_ = rhs.RHS();
148  ID_ = rhs.ID();
149  }
150 
153  {}
155 
157 
160  {
161  if (&rhs == this)
162  return(*this);
163 
164  NumRows_ = rhs.NumRows();
165  NumVectors_ = rhs.NumVectors();
166  IsComputed_ = rhs.IsComputed();
168  Matrix_ = rhs.Matrix();
171  LHS_ = rhs.LHS();
172  RHS_ = rhs.RHS();
173  ID_ = rhs.ID();
174 
175  return(*this);
176  }
177 
179 
181 
183  virtual int NumRows() const;
184 
186  virtual int NumVectors() const
187  {
188  return(NumVectors_);
189  }
190 
192  virtual int SetNumVectors(const int NumVectors_in)
193  {
194  if (NumVectors_ == NumVectors_in)
195  return(0);
196 
197  NumVectors_ = NumVectors_in;
200  // zero out vector elements
201  for (int i = 0 ; i < NumRows_ ; ++i)
202  for (int j = 0 ; j < NumVectors_ ; ++j) {
203  LHS_(i,j) = 0.0;
204  RHS_(i,j) = 0.0;
205  }
206  if (NumRows_!=0)
207  {
209  }
210  return(0);
211  }
212 
214  virtual double& LHS(const int i, const int Vector = 0);
215 
217  virtual double& RHS(const int i, const int Vector = 0);
218 
220 
229  virtual int& ID(const int i);
230 
232  virtual int SetMatrixElement(const int row, const int col,
233  const double value);
234 
236  virtual int SetParameters(Teuchos::ParameterList& /* List */)
237  {
238  return(0);
239  }
240 
242  virtual bool IsInitialized() const
243  {
244  return(IsInitialized_);
245  }
246 
248  virtual bool IsComputed() const
249  {
250  return(IsComputed_);
251  }
252 
254  virtual const char* Label() const
255  {
256  return(Label_.c_str());
257  }
258 
260  virtual int SetKeepNonFactoredMatrix(const bool flag)
261  {
262  KeepNonFactoredMatrix_ = flag;
263  return(0);
264  }
265 
267  virtual bool KeepNonFactoredMatrix() const
268  {
269  return(KeepNonFactoredMatrix_);
270  }
271 
273  virtual const Epetra_SerialDenseMatrix& LHS() const
274  {
275  return(LHS_);
276  }
277 
279  virtual const Epetra_SerialDenseMatrix& RHS() const
280  {
281  return(RHS_);
282  }
283 
285  virtual const Epetra_SerialDenseMatrix& Matrix() const
286  {
287  return(Matrix_);
288  }
289 
292  {
293  return(NonFactoredMatrix_);
294  }
295 
297  virtual const Epetra_IntSerialDenseVector& ID() const
298  {
299  return(ID_);
300  }
301 
303 
305  virtual int Initialize();
307 
309  virtual int Compute(const Epetra_RowMatrix& Matrix_in);
310 
312  virtual int Apply();
313 
315  virtual int ApplyInverse();
316 
318 
319  virtual double InitializeFlops() const
320  {
321  return(0.0);
322  }
323 
324  virtual double ComputeFlops() const
325  {
326  return(ComputeFlops_);
327  }
328 
329  virtual double ApplyFlops() const
330  {
331  return(ApplyFlops_);
332  }
333 
334  virtual double ApplyInverseFlops() const
335  {
336  return(ApplyInverseFlops_);
337  }
338 
340  virtual std::ostream& Print(std::ostream& os) const;
341 
342 private:
343 
345  virtual int Extract(const Epetra_RowMatrix& Matrix_in);
346 
348  int NumRows_;
370  std::string Label_;
371 
375  double ApplyFlops_;
378 };
379 
380 #endif
Ifpack_DenseContainer(const int NumRows_in, const int NumVectors_in=1)
Default constructor.
virtual double ComputeFlops() const
Returns the flops in Compute().
Epetra_SerialDenseMatrix NonFactoredMatrix_
Dense matrix, that contains the non-factored matrix.
virtual const Epetra_IntSerialDenseVector & ID() const
Returns the integer dense vector of IDs.
virtual int Compute(const Epetra_RowMatrix &Matrix_in)
Finalizes the linear system matrix and prepares for the application of the inverse.
virtual const char * Label() const
Returns the label of this container.
virtual std::ostream & Print(std::ostream &os) const
Prints basic information on iostream. This function is used by operator&lt;&lt;.
virtual int SetNumVectors(const int NumVectors_in)
Sets the number of vectors for LHS/RHS.
virtual int SetMatrixElement(const int row, const int col, const double value)
Set the matrix element (row,col) to value.
virtual int SetKeepNonFactoredMatrix(const bool flag)
If flag is true, keeps a copy of the non-factored matrix.
virtual const Epetra_SerialDenseMatrix & Matrix() const
Returns the dense matrix or its factors.
int SetVectors(Epetra_SerialDenseMatrix &X, Epetra_SerialDenseMatrix &B)
virtual int NumVectors() const
Returns the number of vectors in LHS/RHS.
Epetra_SerialDenseMatrix RHS_
Dense vector representing the RHS.
virtual double InitializeFlops() const
Returns the flops in Initialize().
virtual double ApplyInverseFlops() const
Returns the flops in ApplyInverse().
Epetra_SerialDenseSolver Solver_
Dense solver (solution will be get using LAPACK).
virtual bool IsInitialized() const
Returns true is the container has been successfully initialized.
virtual ~Ifpack_DenseContainer()
Destructor.
virtual int Extract(const Epetra_RowMatrix &Matrix_in)
Extract the submatrices identified by the ID set int ID().
bool IsInitialized_
If true, the container has been successfully initialized.
Ifpack_DenseContainer(const Ifpack_DenseContainer &rhs)
Copy constructor.
virtual int SetParameters(Teuchos::ParameterList &)
Sets all necessary parameters.
double ApplyFlops_
Flops in Apply().
virtual double & LHS(const int i, const int Vector=0)
Returns the i-th component of the vector Vector of LHS.
virtual int & ID(const int i)
Returns the ID associated to local row i.
virtual int Apply()
Apply the matrix to RHS, results are stored in LHS.
virtual bool IsComputed() const
Returns true is the container has been successfully computed.
int NumRows_
Number of rows in the container.
virtual const Epetra_SerialDenseMatrix & LHS() const
Returns the dense vector containing the LHS.
std::string Label_
Label for this object.
Ifpack_DenseContainer & operator=(const Ifpack_DenseContainer &rhs)
Operator=.
Epetra_IntSerialDenseVector ID_
Sets of local rows.
virtual int Initialize()
Initialize the container.
#define false
int Reshape(int NumRows, int NumCols)
virtual int ApplyInverse()
Apply the inverse of the matrix to RHS, results are stored in LHS.
virtual int NumRows() const
Returns the number of rows of the matrix and LHS/RHS.
Epetra_SerialDenseMatrix Matrix_
Dense matrix.
int NumVectors_
Number of vectors in the container.
bool IsComputed_
If true, the container has been successfully computed.
Ifpack_DenseContainer: a class to define containers for dense matrices.
double ComputeFlops_
Flops in Compute().
Ifpack_Container: a pure virtual class for creating and solving local linear problems.
virtual double & RHS(const int i, const int Vector=0)
Returns the i-th component of the vector Vector of RHS.
double ApplyInverseFlops_
Flops in ApplyInverse().
virtual double ApplyFlops() const
Returns the flops in Apply().
virtual const Epetra_SerialDenseMatrix & NonFactoredMatrix() const
Returns the non-factored dense matrix (only if stored).
#define IFPACK_CHK_ERR(ifpack_err)
bool KeepNonFactoredMatrix_
If true, keeps a copy of the non-factored matrix.
virtual bool KeepNonFactoredMatrix() const
Returns KeepNonFactoredMatrix_.
Epetra_SerialDenseMatrix LHS_
Dense vector representing the LHS.
virtual const Epetra_SerialDenseMatrix & RHS() const
Returns the dense vector containing the RHS.