fei_VectorTraits_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_VectorTraits_Epetra_h_
00010 #define _fei_VectorTraits_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_base.hpp or fei_Vector.hpp !!!
00017 //
00018 #include <fei_VectorTraits.hpp>
00019 #include <fei_Include_Trilinos.hpp>
00020 
00021 namespace fei {
00028   template<>
00029   struct VectorTraits<Epetra_MultiVector> {
00030     static const char* typeName()
00031       { return("Epetra_MultiVector"); }
00032 
00033     static int setValues(Epetra_MultiVector* vec, int firstLocalOffset,
00034                          double scalar, bool isSolnVector=false)
00035       {
00036         return( vec->PutScalar(scalar) );
00037       }
00038 
00039     //note that incoming indices are point-entry indices, not block-indices.
00040     static int putValuesIn(Epetra_MultiVector* vec,
00041                            int firstLocalOffset,
00042                            int numValues,
00043                            const int* indices,
00044                            const double* values,
00045                            bool sum_into,
00046                            bool isSolnVector=false,
00047                            int vectorIndex=0)
00048       {
00049         double* localVecValues = (*vec)[vectorIndex];
00050         if (sum_into) {
00051           for(int i=0; i<numValues; ++i) {
00052             localVecValues[indices[i]-firstLocalOffset] += values[i];
00053           }
00054         }
00055         else {
00056           for(int i=0; i<numValues; ++i) {
00057             localVecValues[indices[i]-firstLocalOffset] = values[i];
00058           }
00059         }
00060         return(0);
00061       }
00062 
00063     //note that incoming indices are point-entry indices, not block-indices.
00064     static int copyOut(Epetra_MultiVector* vec,
00065                        int firstLocalOffset,
00066                        int numValues, const int* indices, double* values,
00067                        bool isSolnVector=false,
00068                        int vectorIndex=0)
00069       {
00070         double* localVecValues = (*vec)[vectorIndex];
00071         for(int i=0; i<numValues; ++i) {
00072           values[i] = localVecValues[indices[i]-firstLocalOffset];
00073         }
00074 
00075         return(0);
00076       }
00077 
00078     static double* getLocalCoefsPtr(Epetra_MultiVector* vec,
00079                                     bool isSolnVector=false,
00080                                     int vectorIndex=0)
00081       {
00082         return((*vec)[vectorIndex]);
00083       }
00084 
00085     static int update(Epetra_MultiVector* vec,
00086                       double a,
00087                       const Epetra_MultiVector* x,
00088                       double b)
00089     {
00090       return( vec->Update(a, *x, b) );
00091     }
00092 
00093   };//struct VectorTraits<Epetra_MultiVector>
00094 }//namespace fei
00095 
00096 #endif //HAVE_FEI_EPETRA
00097 
00098 #endif // _fei_VectorTraits_Epetra_hpp_

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