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 terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00009 // license for use of this work by or on behalf of the U.S. Government.
00010 // 
00011 // This library is free software; you can redistribute it and/or modify
00012 // it under the terms of the GNU Lesser General Public License as
00013 // published by the Free Software Foundation; either version 2.1 of the
00014 // License, or (at your option) any later version.
00015 //  
00016 // This library is distributed in the hope that it will be useful, but
00017 // WITHOUT ANY WARRANTY; without even the implied warranty of
00018 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00019 // Lesser General Public License for more details.
00020 //  
00021 // You should have received a copy of the GNU Lesser General Public
00022 // License along with this library; if not, write to the Free Software
00023 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00024 // USA
00025 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00026 // 
00027 // ************************************************************************
00028 //@HEADER
00029 
00030 #include "Epetra_IntSerialDenseMatrix.h"
00031 #include "Epetra_Util.h"
00032 
00033 //=============================================================================
00034 Epetra_IntSerialDenseMatrix::Epetra_IntSerialDenseMatrix()
00035   : Epetra_Object("Epetra::IntSerialDenseMatrix"),
00036     CV_(Copy),
00037     A_Copied_(false),
00038     M_(0),
00039     N_(0),
00040     LDA_(0),
00041     A_(0)
00042 {
00043 }
00044 
00045 //=============================================================================
00046 Epetra_IntSerialDenseMatrix::Epetra_IntSerialDenseMatrix(int NumRows, int NumCols)
00047   : Epetra_Object("Epetra::IntSerialDenseMatrix"),
00048     CV_(Copy),
00049     A_Copied_(false),
00050     M_(0),
00051     N_(0),
00052     LDA_(0),
00053     A_(0)
00054 {
00055   if(NumRows < 0)
00056     throw ReportError("NumRows = " + toString(NumRows) + ". Should be >= 0", -1);
00057   if(NumCols < 0)
00058     throw ReportError("NumCols = " + toString(NumCols) + ". Should be >= 0", -1);
00059 
00060   int errorcode = Shape(NumRows, NumCols);
00061   if(errorcode != 0)
00062     throw ReportError("Shape returned non-zero (" + toString(errorcode) + ").", -2);
00063 }
00064 //=============================================================================
00065 Epetra_IntSerialDenseMatrix::Epetra_IntSerialDenseMatrix(Epetra_DataAccess CV, int* A, int LDA, 
00066                                                          int NumRows, int NumCols)
00067   : Epetra_Object("Epetra::IntSerialDenseMatrix"),
00068     CV_(CV),
00069     A_Copied_(false),
00070     M_(NumRows),
00071     N_(NumCols),
00072     LDA_(LDA),
00073     A_(A)
00074 {
00075   if(A == 0)
00076     throw ReportError("Null pointer passed as A parameter.", -3);
00077   if(NumRows < 0)
00078     throw ReportError("NumRows = " + toString(NumRows) + ". Should be >= 0", -1);
00079   if(NumCols < 0)
00080     throw ReportError("NumCols = " + toString(NumCols) + ". Should be >= 0", -1);
00081   if(LDA < 0)
00082     throw ReportError("LDA = " + toString(LDA) + ". Should be >= 0", -1);
00083 
00084   if(CV == Copy) {
00085     LDA_ = M_;
00086     const int newsize = LDA_ * N_;
00087     if(newsize > 0) {
00088       A_ = new int[newsize];
00089       CopyMat(A, LDA, M_, N_, A_, LDA_);
00090       A_Copied_ = true;
00091     }
00092     else {
00093       A_ = 0;
00094     }
00095   }
00096 }
00097 //=============================================================================
00098 Epetra_IntSerialDenseMatrix::Epetra_IntSerialDenseMatrix(const Epetra_IntSerialDenseMatrix& Source)
00099   : Epetra_Object(Source),
00100     CV_(Source.CV_),
00101     A_Copied_(false),
00102     M_(Source.M_),
00103     N_(Source.N_),
00104     LDA_(Source.LDA_),
00105     A_(Source.A_)
00106 {
00107   if(CV_ == Copy) {
00108     LDA_ = M_;
00109     const int newsize = LDA_ * N_;
00110     if(newsize > 0) {
00111       A_ = new int[newsize];
00112       CopyMat(Source.A_, Source.LDA_, M_, N_, A_, LDA_);
00113       A_Copied_ = true;
00114     }
00115     else {
00116       A_ = 0;
00117       A_Copied_ = false;
00118     }
00119   }
00120 }
00121 //=============================================================================
00122 int Epetra_IntSerialDenseMatrix::Reshape(int NumRows, int NumCols) {
00123   if(NumRows < 0 || NumCols < 0)
00124     return(-1);
00125 
00126   int* A_tmp = 0;
00127   const int newsize = NumRows * NumCols;
00128 
00129   if(newsize > 0) {
00130     // Allocate space for new matrix
00131     A_tmp = new int[newsize];
00132     for(int k = 0; k < newsize; k++) 
00133       A_tmp[k] = 0; // Zero out values
00134     int M_tmp = EPETRA_MIN(M_, NumRows);
00135     int N_tmp = EPETRA_MIN(N_, NumCols);
00136     if(A_ != 0) 
00137       CopyMat(A_, LDA_, M_tmp, N_tmp, A_tmp, NumRows); // Copy principal submatrix of A to new A
00138   }
00139   
00140   CleanupData(); // Get rid of anything that might be already allocated  
00141   M_ = NumRows;
00142   N_ = NumCols;
00143   LDA_ = M_;
00144   A_ = A_tmp; // Set pointer to new A
00145   A_Copied_ = (newsize>0);
00146   return(0);
00147 }
00148 //=============================================================================
00149 int Epetra_IntSerialDenseMatrix::Shape(int NumRows, int NumCols) {
00150   if(NumRows < 0 || NumCols < 0)
00151     return(-1);
00152 
00153   CleanupData(); // Get rid of anything that might be already allocated
00154   M_ = NumRows;
00155   N_ = NumCols;
00156   LDA_ = M_;
00157   const int newsize = LDA_ * N_;
00158   if(newsize > 0) {
00159     A_ = new int[newsize];
00160     for(int k = 0; k < newsize; k++)
00161       A_[k] = 0; // Zero out values
00162     A_Copied_ = true;
00163   }
00164 
00165   return(0);
00166 }
00167 //=============================================================================
00168 Epetra_IntSerialDenseMatrix::~Epetra_IntSerialDenseMatrix()
00169 {
00170   CleanupData();
00171 }
00172 //=============================================================================
00173 void Epetra_IntSerialDenseMatrix::CleanupData()
00174 {
00175   if(A_Copied_)
00176     delete[] A_; 
00177   A_ = 0; 
00178   A_Copied_ = false;
00179   M_ = 0;
00180   N_ = 0;
00181   LDA_ = 0;
00182 }
00183 //=============================================================================
00184 Epetra_IntSerialDenseMatrix& Epetra_IntSerialDenseMatrix::operator = (const Epetra_IntSerialDenseMatrix& Source) {
00185   if(this == &Source)
00186     return(*this); // Special case of source same as target
00187   if((CV_ == View) && (Source.CV_ == View) && (A_ == Source.A_))
00188     return(*this); // Special case of both are views to same data.
00189 
00190   if(strcmp(Label(), Source.Label()) != 0)
00191     throw ReportError("operator= type mismatch (lhs = " + string(Label()) + 
00192                       ", rhs = " + string(Source.Label()) + ").", -5);
00193   
00194   if(Source.CV_ == View) {
00195     if(CV_ == Copy) { // C->V only
00196       CleanupData();
00197       CV_ = View;
00198     }
00199     M_ = Source.M_; // C->V and V->V
00200     N_ = Source.N_;
00201     LDA_ = Source.LDA_;
00202     A_ = Source.A_;
00203   }
00204   else {
00205     if(CV_ == View) { // V->C
00206       CV_ = Copy;
00207       M_ = Source.M_;
00208       N_ = Source.N_;
00209       LDA_ = Source.M_;
00210       const int newsize = LDA_ * N_;
00211       if(newsize > 0) {
00212         A_ = new int[newsize];
00213         A_Copied_ = true;
00214       }
00215       else {
00216         A_ = 0;
00217         A_Copied_ = false;
00218       }
00219     }
00220     else { // C->C
00221       if((Source.M_ <= LDA_) && (Source.N_ == N_)) { // we don't need to reallocate
00222         M_ = Source.M_;
00223         N_ = Source.N_;
00224       }
00225       else { // we need to allocate more space (or less space)
00226         CleanupData();
00227         M_ = Source.M_;
00228         N_ = Source.N_;
00229         LDA_ = Source.M_;
00230         const int newsize = LDA_ * N_;
00231         if(newsize > 0) {
00232           A_ = new int[newsize];
00233           A_Copied_ = true;
00234         }
00235       }
00236     }
00237     CopyMat(Source.A_, Source.LDA_, M_, N_, A_, LDA_); // V->C and C->C
00238   }
00239   
00240   return(*this);
00241 }
00242 
00243 //=============================================================================
00244 bool Epetra_IntSerialDenseMatrix::operator==(const Epetra_IntSerialDenseMatrix& rhs) const
00245 {
00246   if (M_ != rhs.M_ || N_ != rhs.N_) return(false);
00247 
00248   const int* A = A_;
00249   const int* rhsA = rhs.A_;
00250 
00251   for(int j=0; j<N_; ++j) {
00252     int offset = j*LDA_;
00253     int rhsOffset = j*rhs.LDA_;
00254     for(int i=0; i<M_; ++i) {
00255       if (A[offset+i] != rhsA[rhsOffset+i]) {
00256   return(false);
00257       }
00258     }
00259   }
00260 
00261   return(true);
00262 }
00263 
00264 //=============================================================================
00265 int Epetra_IntSerialDenseMatrix::MakeViewOf(const Epetra_IntSerialDenseMatrix& Source) {
00266   if(strcmp(Label(), Source.Label()) != 0)
00267     return(-1);
00268 
00269   if(CV_ == Copy) {
00270     CleanupData();
00271     CV_ = View;
00272   }
00273   M_ = Source.M_;
00274   N_ = Source.N_;
00275   LDA_ = Source.LDA_;
00276   A_ = Source.A_;
00277 
00278   return(0);
00279 }
00280 
00281 //=============================================================================
00282 void Epetra_IntSerialDenseMatrix::CopyMat(int* Source, int Source_LDA, int NumRows, int NumCols, 
00283                                           int* Target, int Target_LDA) 
00284 {
00285   int i, j;
00286   int* targetPtr = Target;
00287   int* sourcePtr = 0;
00288   for(j = 0; j < NumCols; j++) { // for each column
00289     targetPtr = Target + j * Target_LDA; // set pointers to correct stride
00290     sourcePtr = Source + j * Source_LDA;
00291     for(i = 0; i < NumRows; i++) // for each row
00292       *targetPtr++ = *sourcePtr++; // copy element (i,j) and increment pointer to (i,j+1)
00293   }
00294   return;
00295 }
00296 
00297 //=============================================================================
00298 int Epetra_IntSerialDenseMatrix::OneNorm() {
00299   int anorm = 0;
00300   int* ptr = 0;
00301   for(int j = 0; j < N_; j++) {
00302     int sum = 0;
00303     ptr = A_ + j*LDA_;
00304     for(int i = 0; i < M_; i++) 
00305       sum += std::abs(*ptr++);
00306     anorm = EPETRA_MAX(anorm, sum);
00307   }
00308   return(anorm);
00309 }
00310 
00311 //=============================================================================
00312 int Epetra_IntSerialDenseMatrix::InfNorm() {  
00313   int anorm = 0;
00314   int* ptr = 0;
00315   // Loop across columns in inner loop.  Most expensive memory access, but 
00316   // requires no extra storage.
00317   for(int i = 0; i < M_; i++) {
00318     int sum = 0;
00319     ptr = A_ + i;
00320     for(int j = 0; j < N_; j++) {
00321       sum += std::abs(*ptr);
00322       ptr += LDA_;
00323     }
00324     anorm = EPETRA_MAX(anorm, sum);
00325   }
00326   return(anorm);
00327 }
00328 
00329 //=========================================================================
00330 void Epetra_IntSerialDenseMatrix::Print(ostream& os) const {
00331   if(CV_ == Copy)
00332     os << "Data access mode: Copy" << endl;
00333   else
00334     os << "Data access mode: View" << endl;
00335   if(A_Copied_)
00336     os << "A_Copied: yes" << endl;
00337   else
00338     os << "A_Copied: no" << endl;
00339   os << "Rows(M): " << M_ << endl;
00340   os << "Columns(N): " << N_ << endl;
00341   os << "LDA: " << LDA_ << endl;
00342   if(M_ == 0 || N_ == 0)
00343     os << "(matrix is empty, no values to display)" << endl;
00344   else
00345     for(int i = 0; i < M_; i++) {
00346       for(int j = 0; j < N_; j++){
00347         os << (*this)(i,j) << " ";
00348       }
00349       os << endl;
00350     }
00351 }
00352 
00353 //=========================================================================
00354 int Epetra_IntSerialDenseMatrix::Random() {
00355 
00356   Epetra_Util util;
00357 
00358   for(int j = 0; j < N_; j++) {
00359     int* arrayPtr = A_ + (j * LDA_);
00360     for(int i = 0; i < M_; i++) {
00361       *arrayPtr++ = util.RandomInt();
00362     }
00363   }
00364   
00365   return(0);
00366 }

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