FEI Version of the Day
fei_AztecDMSR_Matrix.hpp
00001 /*
00002 // @HEADER
00003 // ************************************************************************
00004 //             FEI: Finite Element Interface to Linear Solvers
00005 //                  Copyright (2005) Sandia Corporation.
00006 //
00007 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation, the
00008 // U.S. Government retains certain rights in this software.
00009 //
00010 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Alan Williams (william@sandia.gov) 
00038 //
00039 // ************************************************************************
00040 // @HEADER
00041 */
00042 
00043 #ifndef _AztecDMSR_Matrix_h_
00044 #define _AztecDMSR_Matrix_h_
00045 
00046 #ifdef HAVE_FEI_AZTECOO
00047 
00048 
00049 //
00050 // This class is a wrapper for the Aztec DMSR matrix data structure.
00051 //
00052 // Important usage notes:
00053 //
00054 // * The 'oneBased' argument to the constructor indicates whether the
00055 //   matrix should use 1-based indices (row numbers and column indices) in
00056 //   the input and output arguments to its interfaces (e.g., getRow),
00057 //   with the exception of the update list -- keep reading.
00058 //   'oneBased' should be 1 for 1-based indices, 0 for 0-based.
00059 //   Here's the confusing part -- the update list should contain 0-based
00060 //   indices, regardless of the value of 'oneBased'.  That's because the
00061 //   update list gets used internally by Aztec functions that only work
00062 //   in 0-based numbers.
00063 //
00064 // * The 'rowLengths' array, input argument to the configure function,
00065 //   must contain the lengths of each row, *NOT* including the
00066 //   coefficient on the diagonal.
00067 //
00068 #include <az_aztec.h>
00069 #include <fei_SharedPtr.hpp>
00070 #include <fei_Aztec_Map.hpp>
00071 #include "fei_iostream.hpp"
00072 #include "fei_fstream.hpp"
00073 #include "fei_sstream.hpp"
00074 
00075 namespace fei_trilinos {
00076 
00077 class Aztec_LSVector;
00078 
00079 class AztecDMSR_Matrix {
00080     
00081   public:
00082     // Constructor.
00083     AztecDMSR_Matrix(fei::SharedPtr<Aztec_Map> map);
00084 
00085     //Copy constructor
00086     AztecDMSR_Matrix(const AztecDMSR_Matrix& src);
00087 
00088     virtual ~AztecDMSR_Matrix ();
00089 
00090     // Mathematical functions.
00091     void matvec(const Aztec_LSVector& x, Aztec_LSVector& y) const;
00092 
00093     void put(double s);
00094     void getDiagonal(Aztec_LSVector& diagVector) const;
00095 
00096     fei::SharedPtr<Aztec_Map> getAztec_Map() const {return(amap_);};
00097 
00098     int rowLength(int row) const;
00099     
00100     // ... to read matrix.
00101     void getRow(int row, int& length, double *coefs, int *colInd) const;
00102     void getRow(int row, int& length, double *coefs) const;
00103     void getRow(int row, int& length, int *colInd) const;
00104 
00106     int setDiagEntry(int row, double value);
00107 
00109     double getDiagEntry(int row) const;
00110 
00111     // ... to write matrix.
00112     int putRow(int row, int len, const double *coefs, 
00113                        const int *colInd);
00114 
00115     int sumIntoRow(int numRows, const int* rows,
00116                  int numCols, const int* cols,
00117                  const double* const* coefs);
00118 
00119     int sumIntoRow(int row, int len, const double *coefs, 
00120                            const int *colInd);
00121 
00122     int addScaledMatrix(double scalar, const AztecDMSR_Matrix& source);
00123 
00124     void scale(double scalar);
00125 
00128     int getOffDiagRowPointers(int row, int*& colIndices, double*& coefs,
00129             int& offDiagRowLength);
00130 
00131     void allocate(int *rowLengths);
00132 
00133     //inform about structure, including column-indices, so that val and bindx
00134     //can be allocated *and* so that bindx can be populated.
00135     void allocate(int *rowLengths, const int* const* colIndices);
00136 
00137     //inform that data fill is complete, so AZ_transform can be called.
00138     void fillComplete();
00139 
00140     bool isFilled() const {return(isFilled_);};
00141     void setFilled(bool flag) {isFilled_ = flag;};
00142     bool isAllocated() const {return(isAllocated_);};
00143     void setAllocated(bool flag) {isAllocated_ = flag;};
00144 
00145     void copyStructure(AztecDMSR_Matrix& source);
00146 
00147     bool readFromFile(const char *filename);
00148     bool writeToFile(const char *fileName) const;
00149     bool rowMax() const {return true;};
00150     double rowMax(int row) const;
00151  
00152     int getNumNonZeros() {return(nnzeros_);};
00153 
00154     double* getBeginPointer() { return val; }
00155 
00156     int getOffset(int row, int col)
00157     {
00158       int localRow;
00159       if (!amap_->inUpdate(row,localRow)){
00160         throw std::runtime_error("row not found");
00161       }
00162 
00163       if (row == col) return localRow;
00164 
00165       int* row_ptr = &bindx[bindx[localRow]];
00166       int* end_row = &bindx[bindx[localRow+1]];
00167 
00168       int col_offset = 0;
00169       for(; row_ptr != end_row; ++row_ptr) {
00170         if (amap_->getTransformedEqn(*row_ptr) == col) break;
00171         ++col_offset;
00172       }
00173       if (row_ptr == end_row){
00174         FEI_OSTRINGSTREAM osstr;
00175         osstr << "Col "<<col << " not found for row "<<row;
00176         throw std::runtime_error(osstr.str());
00177       }
00178       return bindx[localRow] + col_offset;
00179     }
00180 
00181     //Aztec-specific functions:
00182 
00183     AZ_MATRIX* getAZ_MATRIX_PTR() {return(Amat_);};
00184 
00185   private:
00186     void messageAbort(const char* mesg);
00187     int insert(int item, int offset, int* list, int& len, int allocLen);
00188     int insert(double item, int offset, double* list, int& len, int allocLen);
00189     void expand_array(int*& array, int& arraylen, int newlen);
00190     void expand_array(double*& array, int& arraylen, int newlen);
00191 
00192     bool isFilled_;
00193     bool isAllocated_;
00194 
00195     int localOffset_;
00196     int localSize_;
00197 
00198     fei::SharedPtr<Aztec_Map> amap_;
00199 
00200     AZ_MATRIX* Amat_;
00201 
00202     bool arraysAllocated_;
00203     double *val;
00204     int *bindx;
00205     int *rowLengths_;
00206     int nnzeros_; //val and bindx are of length nnzeros_+1
00207 
00208     int N_update_;
00209 
00210     int* tmp_array_;
00211     int tmp_array_len_;
00212     double* dtmp_array_;
00213     int dtmp_array_len_;
00214 
00215     bool azTransformed_;
00216 };
00217 
00218 }//namespace fei_trilinos
00219 
00220 #endif //HAVE_FEI_AZTECOO
00221 
00222 #endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Friends