Belos Package Browser (Single Doxygen Collection) Development
MyBetterOperator.hpp
Go to the documentation of this file.
00001 /*
00002 //@HEADER
00003 // ************************************************************************
00004 //
00005 //                 Belos: Block Linear Solvers Package
00006 //                  Copyright 2004 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 
00044 #ifndef MY_BETTER_OPERATOR_HPP
00045 #define MY_BETTER_OPERATOR_HPP
00046 
00047 #include "BelosConfigDefs.hpp"
00048 #include "BelosOperator.hpp"
00049 #include "MyMultiVec.hpp"
00050 #include "Teuchos_BLAS.hpp"
00051 
00053 
00065 template <class ScalarType>
00066 class MyBetterOperator : public Belos::Operator<ScalarType> 
00067 {
00068 
00069 public:
00070 
00072   MyBetterOperator(const int nrows, const int *colptr,
00073                    const int nnz, const int *rowin, const ScalarType *vals)
00074   : _nr(nrows), _nnz(nnz), _cptr(nrows+1), _rind(nnz), _vals(nnz)
00075   {
00076     std::copy<const int*,IntIter>(colptr,colptr+nrows+1,_cptr.begin());
00077     std::copy<const int*,IntIter>(rowin,rowin+nnz,_rind.begin());
00078     std::copy<const ScalarType*,STIter>(vals,vals+nnz,_vals.begin());
00079   }
00080 
00082   ~MyBetterOperator()
00083   { }
00084 
00086   void Apply(const Belos::MultiVec<ScalarType>& X, 
00087              Belos::MultiVec<ScalarType>& Y,
00088              Belos::ETrans trans = Belos::NOTRANS) const
00089   {
00090     const MyMultiVec<ScalarType>* MyX;
00091     MyX = dynamic_cast<const MyMultiVec<ScalarType>*>(&X); 
00092     assert (MyX != 0);
00093     
00094     MyMultiVec<ScalarType>* MyY;
00095     MyY = dynamic_cast<MyMultiVec<ScalarType>*>(&Y); 
00096     assert (MyY != 0);
00097 
00098     // Initialize output std::vector to zero.
00099     MyY->MvInit( Teuchos::ScalarTraits<ScalarType>::zero() );
00100     
00101     assert (X.GetNumberVecs() == Y.GetNumberVecs());
00102     assert (X.GetVecLength() == Y.GetVecLength());
00103     
00104     int nv = X.GetNumberVecs();
00105 
00106     // Apply operator
00107     int IA1, IA2, ri;
00108     ScalarType aval;
00109     int i,j,v;
00110     for (j=0; j<_nr; j++) {
00111       IA1 = _cptr[j]-1;
00112       IA2 = _cptr[j+1]-1;
00113       for (i=IA1; i<IA2; i++) {
00114   ri = _rind[i]-1;
00115         aval = _vals[i];
00116         for (v=0; v<nv; v++) {
00117           (*MyY)[v][ri] += aval*(*MyX)[v][j];
00118         }
00119       }
00120     }
00121   }
00122 
00123   void Print( std::ostream& os ) {
00124     for (int j=0; j<_nr; j++) {
00125       int IA1 = _cptr[j]-1;
00126       int IA2 = _cptr[j+1]-1;
00127       for (int i=IA1; i<IA2; i++) {
00128         os << "("<<_rind[i]-1<<","<<j<<")\t"<<_vals[i]<< std::endl;
00129       }
00130     } 
00131   }
00132 
00133   private:
00134   typedef typename std::vector<ScalarType>::iterator STIter;
00135   typedef std::vector<int>::iterator        IntIter;
00137   int _nr, _nnz;
00139   std::vector<int> _cptr;
00141   std::vector<int> _rind;
00143   std::vector<ScalarType> _vals;
00144 };
00145 
00146 #endif //MY_BETTER_OPERATOR_HPP
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines