Epetra Package Browser (Single Doxygen Collection) Development
Epetra_IntSerialDenseMatrix.cpp
Go to the documentation of this file.
00001 
00002 //@HEADER
00003 // ************************************************************************
00004 // 
00005 //               Epetra: Linear Algebra Services Package 
00006 //                 Copyright 2001 Sandia Corporation
00007 // 
00008 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00009 // the U.S. Government retains certain rights in this software.
00010 //
00011 // Redistribution and use in source and binary forms, with or without
00012 // modification, are permitted provided that the following conditions are
00013 // met:
00014 //
00015 // 1. Redistributions of source code must retain the above copyright
00016 // notice, this list of conditions and the following disclaimer.
00017 //
00018 // 2. Redistributions in binary form must reproduce the above copyright
00019 // notice, this list of conditions and the following disclaimer in the
00020 // documentation and/or other materials provided with the distribution.
00021 //
00022 // 3. Neither the name of the Corporation nor the names of the
00023 // contributors may be used to endorse or promote products derived from
00024 // this software without specific prior written permission.
00025 //
00026 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00027 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00028 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00029 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00030 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00031 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00032 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00033 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00034 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00035 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00036 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00037 //
00038 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00039 // 
00040 // ************************************************************************
00041 //@HEADER
00042 
00043 #include "Epetra_IntSerialDenseMatrix.h"
00044 #include "Epetra_Util.h"
00045 
00046 //=============================================================================
00047 Epetra_IntSerialDenseMatrix::Epetra_IntSerialDenseMatrix()
00048   : Epetra_Object("Epetra::IntSerialDenseMatrix"),
00049     CV_(Copy),
00050     A_Copied_(false),
00051     M_(0),
00052     N_(0),
00053     LDA_(0),
00054     A_(0)
00055 {
00056 }
00057 
00058 //=============================================================================
00059 Epetra_IntSerialDenseMatrix::Epetra_IntSerialDenseMatrix(int NumRows, int NumCols)
00060   : Epetra_Object("Epetra::IntSerialDenseMatrix"),
00061     CV_(Copy),
00062     A_Copied_(false),
00063     M_(0),
00064     N_(0),
00065     LDA_(0),
00066     A_(0)
00067 {
00068   if(NumRows < 0)
00069     throw ReportError("NumRows = " + toString(NumRows) + ". Should be >= 0", -1);
00070   if(NumCols < 0)
00071     throw ReportError("NumCols = " + toString(NumCols) + ". Should be >= 0", -1);
00072 
00073   int errorcode = Shape(NumRows, NumCols);
00074   if(errorcode != 0)
00075     throw ReportError("Shape returned non-zero (" + toString(errorcode) + ").", -2);
00076 }
00077 //=============================================================================
00078 Epetra_IntSerialDenseMatrix::Epetra_IntSerialDenseMatrix(Epetra_DataAccess CV, int* A, int LDA, 
00079                                                          int NumRows, int NumCols)
00080   : Epetra_Object("Epetra::IntSerialDenseMatrix"),
00081     CV_(CV),
00082     A_Copied_(false),
00083     M_(NumRows),
00084     N_(NumCols),
00085     LDA_(LDA),
00086     A_(A)
00087 {
00088   if(A == 0)
00089     throw ReportError("Null pointer passed as A parameter.", -3);
00090   if(NumRows < 0)
00091     throw ReportError("NumRows = " + toString(NumRows) + ". Should be >= 0", -1);
00092   if(NumCols < 0)
00093     throw ReportError("NumCols = " + toString(NumCols) + ". Should be >= 0", -1);
00094   if(LDA < 0)
00095     throw ReportError("LDA = " + toString(LDA) + ". Should be >= 0", -1);
00096 
00097   if(CV == Copy) {
00098     LDA_ = M_;
00099     const int newsize = LDA_ * N_;
00100     if(newsize > 0) {
00101       A_ = new int[newsize];
00102       CopyMat(A, LDA, M_, N_, A_, LDA_);
00103       A_Copied_ = true;
00104     }
00105     else {
00106       A_ = 0;
00107     }
00108   }
00109 }
00110 //=============================================================================
00111 Epetra_IntSerialDenseMatrix::Epetra_IntSerialDenseMatrix(const Epetra_IntSerialDenseMatrix& Source)
00112   : Epetra_Object(Source),
00113     CV_(Source.CV_),
00114     A_Copied_(false),
00115     M_(Source.M_),
00116     N_(Source.N_),
00117     LDA_(Source.LDA_),
00118     A_(Source.A_)
00119 {
00120   if(CV_ == Copy) {
00121     LDA_ = M_;
00122     const int newsize = LDA_ * N_;
00123     if(newsize > 0) {
00124       A_ = new int[newsize];
00125       CopyMat(Source.A_, Source.LDA_, M_, N_, A_, LDA_);
00126       A_Copied_ = true;
00127     }
00128     else {
00129       A_ = 0;
00130       A_Copied_ = false;
00131     }
00132   }
00133 }
00134 //=============================================================================
00135 int Epetra_IntSerialDenseMatrix::Reshape(int NumRows, int NumCols) {
00136   if(NumRows < 0 || NumCols < 0)
00137     return(-1);
00138 
00139   int* A_tmp = 0;
00140   const int newsize = NumRows * NumCols;
00141 
00142   if(newsize > 0) {
00143     // Allocate space for new matrix
00144     A_tmp = new int[newsize];
00145     for(int k = 0; k < newsize; k++) 
00146       A_tmp[k] = 0; // Zero out values
00147     int M_tmp = EPETRA_MIN(M_, NumRows);
00148     int N_tmp = EPETRA_MIN(N_, NumCols);
00149     if(A_ != 0) 
00150       CopyMat(A_, LDA_, M_tmp, N_tmp, A_tmp, NumRows); // Copy principal submatrix of A to new A
00151   }
00152   
00153   CleanupData(); // Get rid of anything that might be already allocated  
00154   M_ = NumRows;
00155   N_ = NumCols;
00156   LDA_ = M_;
00157   A_ = A_tmp; // Set pointer to new A
00158   A_Copied_ = (newsize>0);
00159   return(0);
00160 }
00161 //=============================================================================
00162 int Epetra_IntSerialDenseMatrix::Shape(int NumRows, int NumCols) {
00163   if(NumRows < 0 || NumCols < 0)
00164     return(-1);
00165 
00166   CleanupData(); // Get rid of anything that might be already allocated
00167   M_ = NumRows;
00168   N_ = NumCols;
00169   LDA_ = M_;
00170   const int newsize = LDA_ * N_;
00171   if(newsize > 0) {
00172     A_ = new int[newsize];
00173 #ifdef Epetra_HAVE_OMP
00174 #pragma omp parallel for
00175 #endif
00176     for(int k = 0; k < newsize; k++)
00177       A_[k] = 0; // Zero out values
00178     A_Copied_ = true;
00179   }
00180 
00181   return(0);
00182 }
00183 //=============================================================================
00184 Epetra_IntSerialDenseMatrix::~Epetra_IntSerialDenseMatrix()
00185 {
00186   CleanupData();
00187 }
00188 //=============================================================================
00189 void Epetra_IntSerialDenseMatrix::CleanupData()
00190 {
00191   if(A_Copied_)
00192     delete[] A_; 
00193   A_ = 0; 
00194   A_Copied_ = false;
00195   M_ = 0;
00196   N_ = 0;
00197   LDA_ = 0;
00198 }
00199 //=============================================================================
00200 Epetra_IntSerialDenseMatrix& Epetra_IntSerialDenseMatrix::operator = (const Epetra_IntSerialDenseMatrix& Source) {
00201   if(this == &Source)
00202     return(*this); // Special case of source same as target
00203   if((CV_ == View) && (Source.CV_ == View) && (A_ == Source.A_))
00204     return(*this); // Special case of both are views to same data.
00205 
00206   if(std::strcmp(Label(), Source.Label()) != 0)
00207     throw ReportError("operator= type mismatch (lhs = " + string(Label()) + 
00208                       ", rhs = " + string(Source.Label()) + ").", -5);
00209   
00210   if(Source.CV_ == View) {
00211     if(CV_ == Copy) { // C->V only
00212       CleanupData();
00213       CV_ = View;
00214     }
00215     M_ = Source.M_; // C->V and V->V
00216     N_ = Source.N_;
00217     LDA_ = Source.LDA_;
00218     A_ = Source.A_;
00219   }
00220   else {
00221     if(CV_ == View) { // V->C
00222       CV_ = Copy;
00223       M_ = Source.M_;
00224       N_ = Source.N_;
00225       LDA_ = Source.M_;
00226       const int newsize = LDA_ * N_;
00227       if(newsize > 0) {
00228         A_ = new int[newsize];
00229         A_Copied_ = true;
00230       }
00231       else {
00232         A_ = 0;
00233         A_Copied_ = false;
00234       }
00235     }
00236     else { // C->C
00237       if((Source.M_ <= LDA_) && (Source.N_ == N_)) { // we don't need to reallocate
00238         M_ = Source.M_;
00239         N_ = Source.N_;
00240       }
00241       else { // we need to allocate more space (or less space)
00242         CleanupData();
00243         M_ = Source.M_;
00244         N_ = Source.N_;
00245         LDA_ = Source.M_;
00246         const int newsize = LDA_ * N_;
00247         if(newsize > 0) {
00248           A_ = new int[newsize];
00249           A_Copied_ = true;
00250         }
00251       }
00252     }
00253     CopyMat(Source.A_, Source.LDA_, M_, N_, A_, LDA_); // V->C and C->C
00254   }
00255   
00256   return(*this);
00257 }
00258 
00259 //=============================================================================
00260 bool Epetra_IntSerialDenseMatrix::operator==(const Epetra_IntSerialDenseMatrix& rhs) const
00261 {
00262   if (M_ != rhs.M_ || N_ != rhs.N_) return(false);
00263 
00264   const int* A = A_;
00265   const int* rhsA = rhs.A_;
00266 
00267   for(int j=0; j<N_; ++j) {
00268     int offset = j*LDA_;
00269     int rhsOffset = j*rhs.LDA_;
00270     for(int i=0; i<M_; ++i) {
00271       if (A[offset+i] != rhsA[rhsOffset+i]) {
00272   return(false);
00273       }
00274     }
00275   }
00276 
00277   return(true);
00278 }
00279 
00280 //=============================================================================
00281 int Epetra_IntSerialDenseMatrix::MakeViewOf(const Epetra_IntSerialDenseMatrix& Source) {
00282   if(std::strcmp(Label(), Source.Label()) != 0)
00283     return(-1);
00284 
00285   if(CV_ == Copy) {
00286     CleanupData();
00287     CV_ = View;
00288   }
00289   M_ = Source.M_;
00290   N_ = Source.N_;
00291   LDA_ = Source.LDA_;
00292   A_ = Source.A_;
00293 
00294   return(0);
00295 }
00296 
00297 //=============================================================================
00298 void Epetra_IntSerialDenseMatrix::CopyMat(int* Source, int Source_LDA, int NumRows, int NumCols, 
00299                                           int* Target, int Target_LDA) 
00300 {
00301   int i, j;
00302   int* targetPtr = Target;
00303   int* sourcePtr = 0;
00304   for(j = 0; j < NumCols; j++) { // for each column
00305     targetPtr = Target + j * Target_LDA; // set pointers to correct stride
00306     sourcePtr = Source + j * Source_LDA;
00307     for(i = 0; i < NumRows; i++) // for each row
00308       *targetPtr++ = *sourcePtr++; // copy element (i,j) and increment pointer to (i,j+1)
00309   }
00310   return;
00311 }
00312 
00313 //=============================================================================
00314 int Epetra_IntSerialDenseMatrix::OneNorm() {
00315   int anorm = 0;
00316   int* ptr = 0;
00317   for(int j = 0; j < N_; j++) {
00318     int sum = 0;
00319     ptr = A_ + j*LDA_;
00320     for(int i = 0; i < M_; i++) 
00321       sum += std::abs(*ptr++);
00322     anorm = EPETRA_MAX(anorm, sum);
00323   }
00324   return(anorm);
00325 }
00326 
00327 //=============================================================================
00328 int Epetra_IntSerialDenseMatrix::InfNorm() {  
00329   int anorm = 0;
00330   int* ptr = 0;
00331   // Loop across columns in inner loop.  Most expensive memory access, but 
00332   // requires no extra storage.
00333   for(int i = 0; i < M_; i++) {
00334     int sum = 0;
00335     ptr = A_ + i;
00336     for(int j = 0; j < N_; j++) {
00337       sum += std::abs(*ptr);
00338       ptr += LDA_;
00339     }
00340     anorm = EPETRA_MAX(anorm, sum);
00341   }
00342   return(anorm);
00343 }
00344 
00345 //=========================================================================
00346 void Epetra_IntSerialDenseMatrix::Print(ostream& os) const {
00347   if(CV_ == Copy)
00348     os << "Data access mode: Copy" << endl;
00349   else
00350     os << "Data access mode: View" << endl;
00351   if(A_Copied_)
00352     os << "A_Copied: yes" << endl;
00353   else
00354     os << "A_Copied: no" << endl;
00355   os << "Rows(M): " << M_ << endl;
00356   os << "Columns(N): " << N_ << endl;
00357   os << "LDA: " << LDA_ << endl;
00358   if(M_ == 0 || N_ == 0)
00359     os << "(matrix is empty, no values to display)" << endl;
00360   else
00361     for(int i = 0; i < M_; i++) {
00362       for(int j = 0; j < N_; j++){
00363         os << (*this)(i,j) << " ";
00364       }
00365       os << endl;
00366     }
00367 }
00368 
00369 //=========================================================================
00370 int Epetra_IntSerialDenseMatrix::Random() {
00371 
00372   Epetra_Util util;
00373 
00374   for(int j = 0; j < N_; j++) {
00375     int* arrayPtr = A_ + (j * LDA_);
00376     for(int i = 0; i < M_; i++) {
00377       *arrayPtr++ = util.RandomInt();
00378     }
00379   }
00380   
00381   return(0);
00382 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines