Epetra_SerialDenseMatrix.h

Go to the documentation of this file.
00001 
00002 //@HEADER
00003 /*
00004 ************************************************************************
00005 
00006               Epetra: Linear Algebra Services Package 
00007                 Copyright (2001) Sandia Corporation
00008 
00009 Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00010 license for use of this work by or on behalf of the U.S. Government.
00011 
00012 This library is free software; you can redistribute it and/or modify
00013 it under the terms of the GNU Lesser General Public License as
00014 published by the Free Software Foundation; either version 2.1 of the
00015 License, or (at your option) any later version.
00016  
00017 This library is distributed in the hope that it will be useful, but
00018 WITHOUT ANY WARRANTY; without even the implied warranty of
00019 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00020 Lesser General Public License for more details.
00021  
00022 You should have received a copy of the GNU Lesser General Public
00023 License along with this library; if not, write to the Free Software
00024 Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00025 USA
00026 Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00027 
00028 ************************************************************************
00029 */
00030 //@HEADER
00031 
00032 #ifndef EPETRA_SERIALDENSEMATRIX_H
00033 #define EPETRA_SERIALDENSEMATRIX_H
00034 
00035 #include "Epetra_Object.h" 
00036 #include "Epetra_CompObject.h"
00037 #include "Epetra_BLAS.h"
00038 #include "Epetra_SerialDenseOperator.h"
00039 class Epetra_SerialSymDenseMatrix;
00040 class Epetra_VbrMatrix;
00041 
00043 
00093 //=========================================================================
00094 class Epetra_SerialDenseMatrix : public Epetra_CompObject, public Epetra_Object, public Epetra_SerialDenseOperator, public Epetra_BLAS {
00095 
00096   public:
00097   
00099 
00100 
00101 
00106   Epetra_SerialDenseMatrix(bool set_object_label=true);
00107 
00109 
00120   Epetra_SerialDenseMatrix(int NumRows, int NumCols, bool set_object_label=true);
00121   
00123 
00138   Epetra_SerialDenseMatrix(Epetra_DataAccess CV, double* A_in, int LDA_in, int NumRows, int NumCols,
00139                            bool set_object_label=true);
00140   
00142   
00143   Epetra_SerialDenseMatrix(const Epetra_SerialDenseMatrix& Source);
00144 
00146   virtual ~Epetra_SerialDenseMatrix ();
00148 
00150 
00151 
00152 
00164   int Shape(int NumRows, int NumCols);
00165   
00167 
00180   int Reshape(int NumRows, int NumCols);
00182 
00184 
00185 
00187 
00206   int Multiply(char TransA, char TransB, double ScalarAB, 
00207          const Epetra_SerialDenseMatrix& A, 
00208          const Epetra_SerialDenseMatrix& B,
00209          double ScalarThis);
00210 
00212   /* This method is intended to imitate the semantics of the matrix-vector
00213     multiplication provided by Epetra's sparse matrices. The 'vector' arguments
00214     are actually matrices; this method will return an error if the
00215     dimensions of 'x' are not compatible. 'y' will be reshaped if necessary.
00216   */
00217   int Multiply(bool transA,
00218                const Epetra_SerialDenseMatrix& x,
00219                Epetra_SerialDenseMatrix& y);
00220 
00222 
00243   int Multiply(char SideA, double ScalarAB, 
00244          const Epetra_SerialSymDenseMatrix& A, 
00245          const Epetra_SerialDenseMatrix& B,
00246          double ScalarThis);
00247 
00249 
00257   int Scale(double ScalarA);
00258 
00260 
00263   virtual double NormOne() const;
00264 
00266   virtual double NormInf() const;
00267 
00269 
00271 
00272 
00274 
00280     Epetra_SerialDenseMatrix & operator = (const Epetra_SerialDenseMatrix& Source);
00281 
00283 
00286     bool operator==(const Epetra_SerialDenseMatrix& rhs) const;
00287 
00289 
00291     bool operator!=(const Epetra_SerialDenseMatrix& rhs) const
00292     { return !(*this == rhs); }
00293 
00295 
00302     Epetra_SerialDenseMatrix & operator += (const Epetra_SerialDenseMatrix& Source);
00303 
00305 
00314       double& operator () (int RowIndex, int ColIndex);
00315 
00317 
00326     const double& operator () (int RowIndex, int ColIndex) const;
00327 
00329 
00339     double* operator [] (int ColIndex);
00340 
00342 
00352     const double* operator [] (int ColIndex) const;
00353     
00355 
00361   int Random();
00362 
00364   int M() const {return(M_);};
00365 
00367   int N() const {return(N_);};
00368 
00370   double* A() const {return(A_);};
00371 
00373   double* A() {return(A_);};
00374 
00376   int LDA() const {return(LDA_);};
00377 
00379   Epetra_DataAccess CV() const {return(CV_);};
00381   
00383 
00384 
00385   virtual void Print(ostream& os) const;
00387 
00389 
00390 
00392 
00395   virtual double OneNorm() const {return(NormOne());};
00396 
00398   virtual double InfNorm() const {return(NormInf());};
00400 
00402 
00403  
00405 
00414     virtual int SetUseTranspose(bool UseTranspose_in) { UseTranspose_ = UseTranspose_in; return (0); }
00415  
00417 
00425   virtual int Apply(const Epetra_SerialDenseMatrix& X, Epetra_SerialDenseMatrix& Y);
00426  
00428 
00437     virtual int ApplyInverse(const Epetra_SerialDenseMatrix & X, Epetra_SerialDenseMatrix & Y)
00438     {
00439       (void)X;//prevents unused variable compiler warning
00440       (void)Y;
00441       return (-1);
00442     }
00443  
00445     virtual const char * Label() const { return Epetra_Object::Label(); }
00446  
00448     virtual bool UseTranspose() const { return UseTranspose_; }
00449  
00451     virtual bool HasNormInf() const { return true; }
00452 
00454     virtual int RowDim() const { return M(); } 
00455  
00457     virtual int ColDim() const { return N(); } 
00459 
00460  protected:
00461 
00462   void CopyMat(double* Source, int Source_LDA, int NumRows, int NumCols,
00463                double* Target, int Target_LDA, bool add=false);
00464   void CleanupData();
00465 
00466   int M_;
00467   int N_;
00468   bool A_Copied_;
00469   Epetra_DataAccess CV_;
00470 
00471   //For performance reasons, it's better if Epetra_VbrMatrix can access the
00472   //LDA_ and A_ members of this class directly without going through an
00473   //accessor method. Rather than making them public members, we'll make
00474   //Epetra_VbrMatrix a friend class.
00475 
00476   friend class Epetra_VbrMatrix;
00477 
00478   int LDA_;
00479   double* A_;
00480 
00481   bool UseTranspose_;
00482 };
00483 
00484 // inlined definitions of op() and op[]
00485 //=========================================================================
00486 inline double& Epetra_SerialDenseMatrix::operator () (int RowIndex, int ColIndex) {
00487 #ifdef HAVE_EPETRA_ARRAY_BOUNDS_CHECK
00488   if (RowIndex >= M_ || RowIndex < 0)
00489     throw ReportError("Row index = " +toString(RowIndex) + 
00490                       " Out of Range 0 - " + toString(M_-1),-1);
00491   if (ColIndex >= N_ || ColIndex < 0) 
00492     throw ReportError("Column index = " +toString(ColIndex) + 
00493                       " Out of Range 0 - " + toString(N_-1),-2);
00494 #endif
00495    return(A_[ColIndex*LDA_ + RowIndex]);
00496 }
00497 //=========================================================================
00498 inline const double& Epetra_SerialDenseMatrix::operator () (int RowIndex, int ColIndex) const {
00499 #ifdef HAVE_EPETRA_ARRAY_BOUNDS_CHECK
00500   if (RowIndex >= M_ || RowIndex < 0) 
00501     throw ReportError("Row index = " +toString(RowIndex) + 
00502                       " Out of Range 0 - " + toString(M_-1),-1);
00503   if (ColIndex >= N_ || ColIndex < 0) 
00504     throw ReportError("Column index = " +toString(ColIndex) + 
00505                       " Out of Range 0 - " + toString(N_-1),-2);
00506 #endif
00507   return(A_[ColIndex*LDA_ + RowIndex]);
00508 }
00509 //=========================================================================
00510 inline double* Epetra_SerialDenseMatrix::operator [] (int ColIndex) {
00511 #ifdef HAVE_EPETRA_ARRAY_BOUNDS_CHECK
00512   if (ColIndex >= N_ || ColIndex < 0) 
00513     throw ReportError("Column index = " +toString(ColIndex) + 
00514                       " Out of Range 0 - " + toString(N_-1),-2);
00515 #endif
00516   return(A_ + ColIndex*LDA_);
00517 }
00518 //=========================================================================
00519 inline const double* Epetra_SerialDenseMatrix::operator [] (int ColIndex) const {
00520 #ifdef HAVE_EPETRA_ARRAY_BOUNDS_CHECK
00521   if (ColIndex >= N_ || ColIndex < 0) 
00522     throw ReportError("Column index = " +toString(ColIndex) + 
00523                       " Out of Range 0 - " + toString(N_-1),-2);
00524 #endif
00525   return(A_+ ColIndex*LDA_);
00526 }
00527 //=========================================================================
00528 
00529 #endif /* EPETRA_SERIALDENSEMATRIX_H */

Generated on Wed May 12 21:41:05 2010 for Epetra Package Browser (Single Doxygen Collection) by  doxygen 1.4.7