Epetra Package Browser (Single Doxygen Collection) Development
Epetra_LongLongSerialDenseMatrix.cpp
Go to the documentation of this file.
00001 
00002 //@HEADER
00003 // ************************************************************************
00004 //
00005 //               Epetra: Linear Algebra Services Package
00006 //                 Copyright 2011 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_ConfigDefs.h"
00044 #include "Epetra_LongLongSerialDenseMatrix.h"
00045 #include "Epetra_Util.h"
00046 
00047 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00048 
00049 //=============================================================================
00050 Epetra_LongLongSerialDenseMatrix::Epetra_LongLongSerialDenseMatrix()
00051   : Epetra_Object("Epetra::LongLongSerialDenseMatrix"),
00052     CV_(Copy),
00053     A_Copied_(false),
00054     M_(0),
00055     N_(0),
00056     LDA_(0),
00057     A_(0)
00058 {
00059 }
00060 
00061 //=============================================================================
00062 Epetra_LongLongSerialDenseMatrix::Epetra_LongLongSerialDenseMatrix(int NumRows, int NumCols)
00063   : Epetra_Object("Epetra::LongLongSerialDenseMatrix"),
00064     CV_(Copy),
00065     A_Copied_(false),
00066     M_(0),
00067     N_(0),
00068     LDA_(0),
00069     A_(0)
00070 {
00071   if(NumRows < 0)
00072     throw ReportError("NumRows = " + toString(NumRows) + ". Should be >= 0", -1);
00073   if(NumCols < 0)
00074     throw ReportError("NumCols = " + toString(NumCols) + ". Should be >= 0", -1);
00075 
00076   int errorcode = Shape(NumRows, NumCols);
00077   if(errorcode != 0)
00078     throw ReportError("Shape returned non-zero (" + toString(errorcode) + ").", -2);
00079 }
00080 //=============================================================================
00081 Epetra_LongLongSerialDenseMatrix::Epetra_LongLongSerialDenseMatrix(Epetra_DataAccess CV_in, long long* A_in, int lda,
00082                                                          int NumRows, int NumCols)
00083   : Epetra_Object("Epetra::LongLongSerialDenseMatrix"),
00084     CV_(CV_in),
00085     A_Copied_(false),
00086     M_(NumRows),
00087     N_(NumCols),
00088     LDA_(lda),
00089     A_(A_in)
00090 {
00091   if(A_in == 0)
00092     throw ReportError("Null pointer passed as A_in parameter.", -3);
00093   if(NumRows < 0)
00094     throw ReportError("NumRows = " + toString(NumRows) + ". Should be >= 0", -1);
00095   if(NumCols < 0)
00096     throw ReportError("NumCols = " + toString(NumCols) + ". Should be >= 0", -1);
00097   if(lda < 0)
00098     throw ReportError("LDA = " + toString(lda) + ". Should be >= 0", -1);
00099 
00100   if(CV_in == Copy) {
00101     LDA_ = M_;
00102     const int newsize = LDA_ * N_;
00103     if(newsize > 0) {
00104       A_ = new long long[newsize];
00105       CopyMat(A_in, lda, M_, N_, A_, LDA_);
00106       A_Copied_ = true;
00107     }
00108     else {
00109       A_ = 0;
00110     }
00111   }
00112 }
00113 //=============================================================================
00114 Epetra_LongLongSerialDenseMatrix::Epetra_LongLongSerialDenseMatrix(const Epetra_LongLongSerialDenseMatrix& Source)
00115   : Epetra_Object(Source),
00116     CV_(Source.CV_),
00117     A_Copied_(false),
00118     M_(Source.M_),
00119     N_(Source.N_),
00120     LDA_(Source.LDA_),
00121     A_(Source.A_)
00122 {
00123   if(CV_ == Copy) {
00124     LDA_ = M_;
00125     const int newsize = LDA_ * N_;
00126     if(newsize > 0) {
00127       A_ = new long long[newsize];
00128       CopyMat(Source.A_, Source.LDA_, M_, N_, A_, LDA_);
00129       A_Copied_ = true;
00130     }
00131     else {
00132       A_ = 0;
00133       A_Copied_ = false;
00134     }
00135   }
00136 }
00137 //=============================================================================
00138 int Epetra_LongLongSerialDenseMatrix::Reshape(int NumRows, int NumCols) {
00139   if(NumRows < 0 || NumCols < 0)
00140     return(-1);
00141 
00142   long long* A_tmp = 0;
00143   const int newsize = NumRows * NumCols;
00144 
00145   if(newsize > 0) {
00146     // Allocate space for new matrix
00147     A_tmp = new long long[newsize];
00148     for(int k = 0; k < newsize; k++)
00149       A_tmp[k] = 0; // Zero out values
00150     int M_tmp = EPETRA_MIN(M_, NumRows);
00151     int N_tmp = EPETRA_MIN(N_, NumCols);
00152     if(A_ != 0)
00153       CopyMat(A_, LDA_, M_tmp, N_tmp, A_tmp, NumRows); // Copy principal submatrix of A to new A
00154   }
00155 
00156   CleanupData(); // Get rid of anything that might be already allocated
00157   M_ = NumRows;
00158   N_ = NumCols;
00159   LDA_ = M_;
00160   A_ = A_tmp; // Set pointer to new A
00161   A_Copied_ = (newsize>0);
00162   return(0);
00163 }
00164 //=============================================================================
00165 int Epetra_LongLongSerialDenseMatrix::Shape(int NumRows, int NumCols) {
00166   if(NumRows < 0 || NumCols < 0)
00167     return(-1);
00168 
00169   CleanupData(); // Get rid of anything that might be already allocated
00170   M_ = NumRows;
00171   N_ = NumCols;
00172   LDA_ = M_;
00173   const int newsize = LDA_ * N_;
00174   if(newsize > 0) {
00175     A_ = new long long[newsize];
00176 #ifdef EPETRA_HAVE_OMP
00177 #pragma omp parallel for
00178 #endif
00179     for(int k = 0; k < newsize; k++)
00180       A_[k] = 0; // Zero out values
00181     A_Copied_ = true;
00182   }
00183 
00184   return(0);
00185 }
00186 //=============================================================================
00187 Epetra_LongLongSerialDenseMatrix::~Epetra_LongLongSerialDenseMatrix()
00188 {
00189   CleanupData();
00190 }
00191 //=============================================================================
00192 void Epetra_LongLongSerialDenseMatrix::CleanupData()
00193 {
00194   if(A_Copied_)
00195     delete[] A_;
00196   A_ = 0;
00197   A_Copied_ = false;
00198   M_ = 0;
00199   N_ = 0;
00200   LDA_ = 0;
00201 }
00202 //=============================================================================
00203 Epetra_LongLongSerialDenseMatrix& Epetra_LongLongSerialDenseMatrix::operator = (const Epetra_LongLongSerialDenseMatrix& Source) {
00204   if(this == &Source)
00205     return(*this); // Special case of source same as target
00206   if((CV_ == View) && (Source.CV_ == View) && (A_ == Source.A_))
00207     return(*this); // Special case of both are views to same data.
00208 
00209   if(std::strcmp(Label(), Source.Label()) != 0)
00210     throw ReportError("operator= type mismatch (lhs = " + std::string(Label()) +
00211       ", rhs = " + std::string(Source.Label()) + ").", -5);
00212 
00213   if(Source.CV_ == View) {
00214     if(CV_ == Copy) { // C->V only
00215       CleanupData();
00216       CV_ = View;
00217     }
00218     M_ = Source.M_; // C->V and V->V
00219     N_ = Source.N_;
00220     LDA_ = Source.LDA_;
00221     A_ = Source.A_;
00222   }
00223   else {
00224     if(CV_ == View) { // V->C
00225       CV_ = Copy;
00226       M_ = Source.M_;
00227       N_ = Source.N_;
00228       LDA_ = Source.M_;
00229       const int newsize = LDA_ * N_;
00230       if(newsize > 0) {
00231         A_ = new long long[newsize];
00232         A_Copied_ = true;
00233       }
00234       else {
00235         A_ = 0;
00236         A_Copied_ = false;
00237       }
00238     }
00239     else { // C->C
00240       if((Source.M_ <= LDA_) && (Source.N_ == N_)) { // we don't need to reallocate
00241         M_ = Source.M_;
00242         N_ = Source.N_;
00243       }
00244       else { // we need to allocate more space (or less space)
00245         CleanupData();
00246         M_ = Source.M_;
00247         N_ = Source.N_;
00248         LDA_ = Source.M_;
00249         const int newsize = LDA_ * N_;
00250         if(newsize > 0) {
00251           A_ = new long long[newsize];
00252           A_Copied_ = true;
00253         }
00254       }
00255     }
00256     CopyMat(Source.A_, Source.LDA_, M_, N_, A_, LDA_); // V->C and C->C
00257   }
00258 
00259   return(*this);
00260 }
00261 
00262 //=============================================================================
00263 bool Epetra_LongLongSerialDenseMatrix::operator==(const Epetra_LongLongSerialDenseMatrix& rhs) const
00264 {
00265   if (M_ != rhs.M_ || N_ != rhs.N_) return(false);
00266 
00267   const long long* A_tmp = A_;
00268   const long long* rhsA = rhs.A_;
00269 
00270   for(int j=0; j<N_; ++j) {
00271     int offset = j*LDA_;
00272     int rhsOffset = j*rhs.LDA_;
00273     for(int i=0; i<M_; ++i) {
00274       if (A_tmp[offset+i] != rhsA[rhsOffset+i]) {
00275   return(false);
00276       }
00277     }
00278   }
00279 
00280   return(true);
00281 }
00282 
00283 //=============================================================================
00284 int Epetra_LongLongSerialDenseMatrix::MakeViewOf(const Epetra_LongLongSerialDenseMatrix& Source) {
00285   if(std::strcmp(Label(), Source.Label()) != 0)
00286     return(-1);
00287 
00288   if(CV_ == Copy) {
00289     CleanupData();
00290     CV_ = View;
00291   }
00292   M_ = Source.M_;
00293   N_ = Source.N_;
00294   LDA_ = Source.LDA_;
00295   A_ = Source.A_;
00296 
00297   return(0);
00298 }
00299 
00300 //=============================================================================
00301 void Epetra_LongLongSerialDenseMatrix::CopyMat(long long* Source, int Source_LDA, int NumRows, int NumCols,
00302                                           long long* Target, int Target_LDA)
00303 {
00304   int i, j;
00305   long long* targetPtr = Target;
00306   long long* sourcePtr = 0;
00307   for(j = 0; j < NumCols; j++) { // for each column
00308     targetPtr = Target + j * Target_LDA; // set pointers to correct stride
00309     sourcePtr = Source + j * Source_LDA;
00310     for(i = 0; i < NumRows; i++) // for each row
00311       *targetPtr++ = *sourcePtr++; // copy element (i,j) and increment pointer to (i,j+1)
00312   }
00313   return;
00314 }
00315 
00316 //=============================================================================
00317 long long Epetra_LongLongSerialDenseMatrix::OneNorm() {
00318   long long anorm = 0;
00319   long long* ptr = 0;
00320   for(int j = 0; j < N_; j++) {
00321     long long sum = 0;
00322     ptr = A_ + j*LDA_;
00323     for(int i = 0; i < M_; i++)
00324     {
00325       const long long val = *ptr++;
00326       sum += (val > 0 ? val : -val); // No std::abs(long long) on VS2005.
00327     }
00328     anorm = EPETRA_MAX(anorm, sum);
00329   }
00330   return(anorm);
00331 }
00332 
00333 //=============================================================================
00334 long long Epetra_LongLongSerialDenseMatrix::InfNorm() {
00335   long long anorm = 0;
00336   long long* ptr = 0;
00337   // Loop across columns in inner loop.  Most expensive memory access, but
00338   // requires no extra storage.
00339   for(int i = 0; i < M_; i++) {
00340     long long sum = 0;
00341     ptr = A_ + i;
00342     for(int j = 0; j < N_; j++) {
00343       const long long val = *ptr;
00344       sum += (val > 0 ? val : -val); // No std::abs(long long) on VS2005.
00345       ptr += LDA_;
00346     }
00347     anorm = EPETRA_MAX(anorm, sum);
00348   }
00349   return(anorm);
00350 }
00351 
00352 //=========================================================================
00353 void Epetra_LongLongSerialDenseMatrix::Print(std::ostream& os) const {
00354   if(CV_ == Copy)
00355     os << "Data access mode: Copy" << std::endl;
00356   else
00357     os << "Data access mode: View" << std::endl;
00358   if(A_Copied_)
00359     os << "A_Copied: yes" << std::endl;
00360   else
00361     os << "A_Copied: no" << std::endl;
00362   os << "Rows(M): " << M_ << std::endl;
00363   os << "Columns(N): " << N_ << std::endl;
00364   os << "LDA: " << LDA_ << std::endl;
00365   if(M_ == 0 || N_ == 0)
00366     os << "(matrix is empty, no values to display)" << std::endl;
00367   else
00368     for(int i = 0; i < M_; i++) {
00369       for(int j = 0; j < N_; j++){
00370         os << (*this)(i,j) << " ";
00371       }
00372       os << std::endl;
00373     }
00374 }
00375 
00376 //=========================================================================
00377 int Epetra_LongLongSerialDenseMatrix::Random() {
00378 
00379   Epetra_Util util;
00380 
00381   for(int j = 0; j < N_; j++) {
00382     long long* arrayPtr = A_ + (j * LDA_);
00383     for(int i = 0; i < M_; i++) {
00384       *arrayPtr++ = util.RandomInt();
00385     }
00386   }
00387 
00388   return(0);
00389 }
00390 
00391 #endif // EPETRA_NO_64BIT_GLOBAL_INDICES
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines