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

Generated on Wed May 12 21:30:41 2010 for FEI by  doxygen 1.4.7