Thyra_DefaultColumnwiseMultiVector.hpp

Go to the documentation of this file.
00001 // @HEADER
00002 // ***********************************************************************
00003 // 
00004 //    Thyra: Interfaces and Support for Abstract Numerical Algorithms
00005 //                 Copyright (2004) Sandia Corporation
00006 // 
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 // 
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //  
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //  
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00025 // 
00026 // ***********************************************************************
00027 // @HEADER
00028 
00029 #ifndef THYRA_MULTI_VECTOR_COLS_HPP
00030 #define THYRA_MULTI_VECTOR_COLS_HPP
00031 
00032 #include "Thyra_DefaultColumnwiseMultiVectorDecl.hpp"
00033 #include "Thyra_MultiVectorDefaultBase.hpp"
00034 #include "Thyra_VectorSpaceBase.hpp"
00035 #include "Thyra_VectorBase.hpp"
00036 #include "Thyra_SingleRhsLinearOpBase.hpp"
00037 #include "Thyra_MultiVectorBase.hpp"
00038 #include "Thyra_VectorSpaceFactoryBase.hpp"
00039 #include "Teuchos_TestForException.hpp"
00040 
00041 namespace Thyra {
00042 
00043 // Constructors/Initializers
00044 
00045 template<class Scalar>
00046 DefaultColumnwiseMultiVector<Scalar>::DefaultColumnwiseMultiVector()
00047 {}
00048 
00049 template<class Scalar>
00050 DefaultColumnwiseMultiVector<Scalar>::DefaultColumnwiseMultiVector(
00051   const Teuchos::RefCountPtr<VectorBase<Scalar> > &col_vec
00052   )
00053 {
00054   this->initialize(col_vec);
00055 }
00056 
00057 template<class Scalar>
00058 DefaultColumnwiseMultiVector<Scalar>::DefaultColumnwiseMultiVector(
00059   const  Teuchos::RefCountPtr<const VectorSpaceBase<Scalar> >   &range
00060   ,const  Teuchos::RefCountPtr<const VectorSpaceBase<Scalar> >  &domain
00061   ,const Teuchos::RefCountPtr<VectorBase<Scalar> >              col_vecs[]
00062   )
00063 {
00064   this->initialize(range,domain,col_vecs);
00065 }
00066 
00067 template<class Scalar>
00068 void DefaultColumnwiseMultiVector<Scalar>::initialize(
00069   const Teuchos::RefCountPtr<VectorBase<Scalar> > &col_vec
00070   )
00071 {
00072 #ifdef TEUCHOS_DEBUG
00073   const char err_msg[] = "DefaultColumnwiseMultiVector<Scalar>::initialize(...): Error!";
00074   TEST_FOR_EXCEPTION( col_vec.get() == NULL,           std::invalid_argument, err_msg ); 
00075   TEST_FOR_EXCEPTION( col_vec->space().get() == NULL,  std::invalid_argument, err_msg ); 
00076 #endif
00077   range_  = col_vec->space();
00078   domain_ = range_->smallVecSpcFcty()->createVecSpc(1);
00079   num_cols_ = 1;
00080   col_vecs_.resize(1);
00081   col_vecs_[0] = col_vec;
00082 }
00083   
00084 template<class Scalar>
00085 void DefaultColumnwiseMultiVector<Scalar>::initialize(
00086   const  Teuchos::RefCountPtr<const VectorSpaceBase<Scalar> >   &range
00087   ,const  Teuchos::RefCountPtr<const VectorSpaceBase<Scalar> >  &domain
00088   ,const Teuchos::RefCountPtr<VectorBase<Scalar> >              col_vecs[]
00089   )
00090 {
00091 #ifdef TEUCHOS_DEBUG
00092   const char err_msg[] = "DefaultColumnwiseMultiVector<Scalar>::initialize(...): Error!";
00093   TEST_FOR_EXCEPTION( range.get()   == NULL, std::invalid_argument, err_msg ); 
00094   TEST_FOR_EXCEPTION( domain.get()  == NULL, std::invalid_argument, err_msg ); 
00095   TEST_FOR_EXCEPTION( range->dim()  == 0,    std::invalid_argument, err_msg ); 
00096   TEST_FOR_EXCEPTION( domain->dim() == 0,    std::invalid_argument, err_msg );
00097   // ToDo: Check the compatibility of the vectors in col_vecs!
00098 #endif
00099   range_ = range;
00100   domain_ = domain;
00101   num_cols_ = domain->dim();
00102   col_vecs_.resize(num_cols_);
00103   if(col_vecs) {
00104     for( Index j = 0; j < num_cols_; ++j )
00105       col_vecs_[j] = col_vecs[j];
00106   }
00107   else {
00108     for( Index j = 0; j < num_cols_; ++j )
00109       col_vecs_[j] = createMember(range_);
00110   }
00111 }
00112 
00113 template<class Scalar>
00114 void DefaultColumnwiseMultiVector<Scalar>::set_uninitialized()
00115 {
00116   col_vecs_.resize(0);
00117   range_ = Teuchos::null;
00118   domain_ = Teuchos::null;
00119 }
00120 
00121 // Overridden from OpBase
00122 
00123 template<class Scalar>
00124 Teuchos::RefCountPtr<const VectorSpaceBase<Scalar> >
00125 DefaultColumnwiseMultiVector<Scalar>::range() const
00126 {
00127   return range_;
00128 }
00129 
00130 template<class Scalar>
00131 Teuchos::RefCountPtr<const VectorSpaceBase<Scalar> >
00132 DefaultColumnwiseMultiVector<Scalar>::domain() const
00133 {
00134   return domain_;
00135 }
00136 
00137 // Overridden from SingleRhsLinearOpBase
00138 
00139 template<class Scalar>
00140 bool DefaultColumnwiseMultiVector<Scalar>::opSupported(ETransp M_trans) const
00141 {
00142   typedef Teuchos::ScalarTraits<Scalar> ST;
00143   return ( ST::isComplex ? ( M_trans==NOTRANS || M_trans==CONJTRANS ) : true );
00144 }
00145 
00146 template<class Scalar>
00147 void DefaultColumnwiseMultiVector<Scalar>::apply(
00148   const ETransp                M_trans
00149   ,const VectorBase<Scalar>    &x
00150   ,VectorBase<Scalar>          *y
00151   ,const Scalar                alpha
00152   ,const Scalar                beta
00153   ) const
00154 {
00155 #ifdef TEUCHOS_DEBUG
00156   THYRA_ASSERT_LINEAR_OP_VEC_APPLY_SPACES("MultiVectorBase<Scalar>::apply()",*this,M_trans,x,y);
00157 #endif
00158   const Index nc = this->domain()->dim();
00159   // y *= beta
00160   Vt_S(y,beta);
00161   // y += alpha*op(M)*x
00162   if(M_trans == NOTRANS) {
00163     //
00164     // y += alpha*M*x = alpha*M.col(0)*x(0) + ... + alpha*M.col(nc-1)*x(nc-1)
00165     //
00166     // Extract an explicit view of x
00167     RTOpPack::ConstSubVectorView<Scalar> x_sub_vec;               
00168     x.acquireDetachedView(Range1D(),&x_sub_vec);
00169     // Loop through and add the multiple of each column
00170     for(Index j = 0; j < nc; ++j )
00171       Vp_StV( y, Scalar(alpha*x_sub_vec(j)), *this->col(j) );
00172     // Release the view of x
00173     x.releaseDetachedView(&x_sub_vec);
00174   }
00175   else {
00176     //
00177     //                   [ alpha*dot(M.col(0),x)    ]
00178     // y += alpha*M'*x = [ alpha*dot(M.col(1),x)    ]
00179     //                   [ ...                      ]
00180     //                   [ alpha*dot(M.col(nc-1),x) ]
00181     //
00182     // Extract an explicit view of y
00183     RTOpPack::SubVectorView<Scalar> y_sub_vec;               
00184     y->acquireDetachedView(Range1D(),&y_sub_vec);
00185     // Loop through and add to each element in y
00186     for(Index j = 0; j < nc; ++j )
00187       y_sub_vec(j) += alpha*dot(*this->col(j),x);
00188     // Commit explicit view of y
00189     y->commitDetachedView(&y_sub_vec);
00190   }
00191 }
00192 
00193 // Overridden from MultiVectorBase
00194 
00195 template<class Scalar>
00196 Teuchos::RefCountPtr<VectorBase<Scalar> >
00197 DefaultColumnwiseMultiVector<Scalar>::col(Index j)
00198 {
00199   TEST_FOR_EXCEPTION(
00200     !(  0 <= j  && j < num_cols_ ), std::logic_error
00201     ,"Error, j = " << j << " does not fall in the range [0,"<<(num_cols_-1)<< "]!"
00202     );
00203   return col_vecs_[j];
00204 }
00205 
00206 template<class Scalar>
00207 Teuchos::RefCountPtr<MultiVectorBase<Scalar> >
00208 DefaultColumnwiseMultiVector<Scalar>::subView( const Range1D& col_rng_in )
00209 {
00210   const Index numCols = domain_->dim();
00211   const Range1D col_rng = Teuchos::full_range(col_rng_in,0,numCols-1);
00212 #ifdef TEUCHOS_DEBUG
00213   TEST_FOR_EXCEPTION(
00214     !( col_rng.ubound() < numCols ), std::logic_error
00215     ,"DefaultColumnwiseMultiVector<Scalar>::subView(col_rng): Error, the input range col_rng = ["<<col_rng.lbound()<<","<<col_rng.ubound()<<"] "
00216     "is not in the range [0,"<<(numCols-1)<<"]!"
00217     );
00218 #endif
00219   return Teuchos::rcp(
00220     new DefaultColumnwiseMultiVector<Scalar>(
00221       range_
00222       ,domain_->smallVecSpcFcty()->createVecSpc(col_rng.size())
00223       ,&col_vecs_[col_rng.lbound()]
00224       )
00225     );
00226 }
00227   
00228 } // end namespace Thyra
00229 
00230 #endif // THYRA_MULTI_VECTOR_COLS_HPP

Generated on Thu Sep 18 12:33:02 2008 for Thyra Package Browser (Single Doxygen Collection) by doxygen 1.3.9.1