fei_MatrixTraits_Epetra.hpp

00001 /*--------------------------------------------------------------------*/
00002 /*    Copyright 2005 Sandia Corporation.                              */
00003 /*    Under the terms of Contract DE-AC04-94AL85000, there is a       */
00004 /*    non-exclusive license for use of this work by or on behalf      */
00005 /*    of the U.S. Government.  Export of this program may require     */
00006 /*    a license from the United States Government.                    */
00007 /*--------------------------------------------------------------------*/
00008 
00009 #ifndef _fei_MatrixTraits_Epetra_h_
00010 #define _fei_MatrixTraits_Epetra_h_
00011 
00012 #include <fei_trilinos_macros.hpp>
00013 
00014 #ifdef HAVE_FEI_EPETRA
00015 
00016 //
00017 //IMPORTANT NOTE: Make sure that wherever this file is included from, it
00018 //appears BEFORE any include of fei_Vector_Impl.hpp or fei_Matrix_Impl.hpp !!!
00019 //
00020 #include <fei_MatrixTraits.hpp>
00021 #include <snl_fei_BlockMatrixTraits.hpp>
00022 #include <fei_VectorTraits_Epetra.hpp>
00023 #include <fei_Include_Trilinos.hpp>
00024 #include <fei_Vector_Impl.hpp>
00025 
00026 namespace fei {
00033   template<>
00034   struct MatrixTraits<Epetra_CrsMatrix> {
00035     static const char* typeName()
00036       { return("Epetra_CrsMatrix"); }
00037 
00038     static int setValues(Epetra_CrsMatrix* mat, double scalar)
00039       {
00040         return( mat->PutScalar(scalar) );
00041       }
00042 
00043     static int getNumLocalRows(Epetra_CrsMatrix* mat, int& numRows)
00044     {
00045       numRows = mat->NumMyRows();
00046       return(0);
00047     }
00048 
00049     static int getRowLength(Epetra_CrsMatrix* mat, int row, int& length)
00050       {
00051         length = mat->NumGlobalEntries(row);
00052         if (length < 0) return(-1);
00053         return( 0 );
00054       }
00055 
00056     static int copyOutRow(Epetra_CrsMatrix* mat,
00057                       int row, int len, double* coefs, int* indices)
00058       {
00059         int dummy;
00060         return(mat->ExtractGlobalRowCopy(row, len, dummy, coefs, indices));
00061       }
00062 
00063     static int putValuesIn(Epetra_CrsMatrix* mat,
00064                      int numRows, const int* rows,
00065                      int numCols, const int* cols,
00066                      const double* const* values,
00067                            bool sum_into)
00068       {
00069         if (sum_into) {
00070           for(int i=0; i<numRows; ++i) {
00071             int err = mat->SumIntoGlobalValues(rows[i], numCols,
00072                                                (double*)values[i],
00073                                                (int*)cols);
00074             if (err != 0) {
00075               return(err);
00076             }
00077           }
00078         }
00079         else {
00080           for(int i=0; i<numRows; ++i) {
00081             int err = mat->ReplaceGlobalValues(rows[i], numCols,
00082                                                (double*)values[i],
00083                                                (int*)cols);
00084             if (err != 0) {
00085               return(err);
00086             }
00087           }
00088         }
00089         return(0);
00090       }
00091 
00092     static int globalAssemble(Epetra_CrsMatrix* mat)
00093     {
00094       if (!mat->Filled()) {
00095         int err = mat->FillComplete();
00096         if (err != 0) {
00097           FEI_CERR << "MatrixTraits<Epetra_CrsMatrix>::globalAssemble"
00098                    << " ERROR in mat->FillComplete" << FEI_ENDL;
00099           return(-1);
00100         }
00101       }
00102 
00103       if (!mat->StorageOptimized()) {
00104         mat->OptimizeStorage();
00105       }
00106 
00107       return( 0 );
00108     }
00109 
00110     static int matvec(Epetra_CrsMatrix* mat,
00111                       fei::Vector* x,
00112                       fei::Vector* y)
00113     {
00114       fei::Vector_Impl<Epetra_MultiVector>* evx =
00115         dynamic_cast<fei::Vector_Impl<Epetra_MultiVector>* >(x);
00116       fei::Vector_Impl<Epetra_MultiVector>* evy =
00117         dynamic_cast<fei::Vector_Impl<Epetra_MultiVector>* >(y);
00118 
00119       if (evx == NULL || evy == NULL) {
00120         return(-1);
00121       }
00122 
00123       Epetra_MultiVector* ex = evx->getUnderlyingVector();
00124       Epetra_MultiVector* ey = evy->getUnderlyingVector();
00125 
00126       return( mat->Multiply(false, *ex, *ey) );
00127     }
00128 
00129   };//struct MatrixTraits<Epetra_CrsMatrix>
00130 }//namespace fei
00131 
00132 namespace snl_fei {
00139   template<>
00140   struct BlockMatrixTraits<Epetra_VbrMatrix> {
00141     static const char* typeName()
00142       { return("Epetra_VbrMatrix"); }
00143 
00144     static int putScalar(Epetra_VbrMatrix* mat, double scalar)
00145       {
00146         return( mat->PutScalar(scalar) );
00147       }
00148 
00149     static int getRowLength(Epetra_VbrMatrix* mat, int row, int& length)
00150       {
00151         length = mat->NumGlobalBlockEntries(row);
00152         return(0);
00153       }
00154 
00155     static int getPointRowLength(Epetra_VbrMatrix* mat, int row, int& length)
00156     {
00157       const Epetra_Map& map = mat->RowMatrixRowMap();
00158       int minLocalRow = map.MinMyGID();
00159       int localRow = row - minLocalRow;
00160       int error = mat->NumMyRowEntries(localRow, length);
00161       return(error);
00162     }
00163 
00164     static int copyOutRow(Epetra_VbrMatrix* mat,
00165                           int row, int numBlkCols,
00166                           int rowDim,
00167                           int* blkCols,
00168                           int* colDims,
00169                           double* coefs,
00170                           int coefsLen,
00171                           int& blkRowLength)
00172       {
00173         int checkRowDim;
00174         int error = mat->BeginExtractGlobalBlockRowCopy(row, numBlkCols,
00175                                                         checkRowDim,
00176                                                         blkRowLength,
00177                                                         blkCols, colDims);
00178         if (error != 0 || checkRowDim != rowDim || blkRowLength != numBlkCols) {
00179           return(error);
00180         }
00181 
00182         int offset = 0;
00183         for(int i=0; i<numBlkCols; ++i) {
00184           if (offset >= coefsLen) {
00185             cerr << "BlockMatrixTraits::copyOutRow ran off end of coefs array."
00186                  << endl;
00187             return(-2);
00188           }
00189           int numValues = rowDim*colDims[i];
00190           error = mat->ExtractEntryCopy(numValues, &(coefs[offset]),
00191                                         rowDim, false);
00192           if (error != 0) {
00193             return(error);
00194           }
00195           offset += numValues;
00196         }
00197 
00198         return(0);
00199       }
00200 
00201     static int copyOutPointRow(Epetra_VbrMatrix* mat,
00202                                int firstLocalOffset,
00203                                int row,
00204                                int len,
00205                                double* coefs,
00206                                int* indices,
00207                                int& rowLength)
00208       {
00209         int error = mat->ExtractMyRowCopy(row-firstLocalOffset,
00210                                           len, rowLength,
00211                                           coefs, indices);
00212 
00213         const Epetra_Map& colmap = mat->RowMatrixColMap();
00214         for(int i=0; i<len; ++i) {
00215           indices[i] = colmap.GID(indices[i]);
00216         }
00217 
00218         return(error);
00219       }
00220 
00221     static int sumIn(Epetra_VbrMatrix* mat,
00222                      int blockRow,
00223                      int rowDim,
00224                      int numBlockCols,
00225                      const int* blockCols,
00226                      const int* colDims,
00227                      int LDA,
00228                      const double* values)
00229     {
00230       int err, *nc_cols = const_cast<int*>(blockCols);
00231       double* nc_values = const_cast<double*>(values);
00232 
00233       err = mat->BeginSumIntoGlobalValues(blockRow, numBlockCols, nc_cols);
00234       if (err != 0) return(err);
00235 
00236       int voffset = 0;
00237       for(int j=0; j<numBlockCols; ++j) {
00238         err = mat->SubmitBlockEntry(&(nc_values[voffset]), LDA,
00239                                     rowDim, colDims[j]);
00240         if (err != 0) return(err);
00241 
00242         voffset += colDims[j]*LDA;
00243       }
00244 
00245       err = mat->EndSubmitEntries();
00246       if (err != 0) return(err);
00247 
00248       return(0);
00249     }
00250 
00251     static int copyIn(Epetra_VbrMatrix* mat,
00252                       int blockRow,
00253                       int rowDim,
00254                       int numBlockCols,
00255                       const int* blockCols,
00256                       const int* colDims,
00257                       int LDA,
00258                       const double* values)
00259     {
00260       int* nc_cols = const_cast<int*>(blockCols);
00261       double* nc_values = const_cast<double*>(values);
00262 
00263       int err = mat->BeginReplaceGlobalValues(blockRow, numBlockCols, nc_cols);
00264       if (err != 0) return(err);
00265 
00266       int voffset = 0;
00267       for(int j=0; j<numBlockCols; ++j) {
00268         err = mat->SubmitBlockEntry(&(nc_values[voffset]), LDA,
00269                                     rowDim, colDims[j]);
00270         if (err != 0) return(err);
00271 
00272         voffset += colDims[j]*LDA;
00273       }
00274 
00275       err = mat->EndSubmitEntries();
00276       if (err != 0) return(err);
00277 
00278       return(0);
00279     }
00280 
00281     static int sumIn(Epetra_VbrMatrix* mat,
00282                      int row,
00283                      int rowDim,
00284                      int numCols,
00285                      const int* cols,
00286                      const int* LDAs,
00287                      const int* colDims,
00288                      const double* const* values)
00289       {
00290         int* nc_cols = const_cast<int*>(cols);
00291         double** nc_values = const_cast<double**>(values);
00292         int err = mat->BeginSumIntoGlobalValues(row,numCols,nc_cols);
00293         if (err != 0) {
00294           return(err);
00295         }
00296         for(int i=0; i<numCols; ++i) {
00297           err = mat->SubmitBlockEntry(nc_values[i], LDAs[i], rowDim, colDims[i]);
00298           if (err != 0) {
00299             return(err);
00300           }
00301         }
00302         err = mat->EndSubmitEntries();
00303 
00304         return(err);
00305       }
00306 
00307     static int copyIn(Epetra_VbrMatrix* mat,
00308                       int row,
00309                       int rowDim,
00310                       int numCols,
00311                       const int* cols,
00312                       const int* LDAs,
00313                       const int* colDims,
00314                       const double* const* values)
00315       {
00316         int* nc_cols = const_cast<int*>(cols);
00317         double** nc_values = const_cast<double**>(values);
00318         int err = mat->BeginReplaceGlobalValues(row, numCols, nc_cols);
00319         if (err != 0) return(err);
00320         for(int i=0; i<numCols; ++i) {
00321           err = mat->SubmitBlockEntry(nc_values[i], LDAs[i], rowDim, colDims[i]);
00322           if (err != 0) return(err);
00323         }
00324         err = mat->EndSubmitEntries();
00325 
00326         return(err);
00327       }
00328 
00329     static int globalAssemble(Epetra_VbrMatrix* mat)
00330     {
00331       return( mat->FillComplete() );
00332     }
00333   };//struct BlockMatrixTraits<Epetra_VbrMatrix>
00334 }//namespace snl_fei
00335 #endif //HAVE_FEI_EPETRA
00336 #endif // _fei_MatrixTraits_Epetra_hpp_

Generated on Tue Jul 13 09:27:45 2010 for FEI by  doxygen 1.4.7