Thyra Package Browser (Single Doxygen Collection) Version of the Day
Thyra_DefaultScaledAdjointLinearOp_def.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 
00030 #ifndef THYRA_DEFAULT_SCALED_ADJOINT_LINEAR_OP_DEF_HPP
00031 #define THYRA_DEFAULT_SCALED_ADJOINT_LINEAR_OP_DEF_HPP
00032 
00033 
00034 #include "Thyra_DefaultScaledAdjointLinearOp_decl.hpp"
00035 #include "Thyra_VectorSpaceBase.hpp"
00036 #include "Thyra_AssertOp.hpp"
00037 
00038 
00039 namespace Thyra {
00040 
00041 
00042 //Constructors/initializers/accessors
00043 
00044 
00045 template<class Scalar>
00046 void DefaultScaledAdjointLinearOp<Scalar>::initialize(
00047   const Scalar &scalar
00048   ,const EOpTransp &transp
00049   ,const Teuchos::RCP<LinearOpBase<Scalar> > &Op
00050   )
00051 {
00052   initializeImpl(scalar,transp,Op,false);
00053 }
00054 
00055 
00056 template<class Scalar>
00057 void DefaultScaledAdjointLinearOp<Scalar>::initialize(
00058   const Scalar &scalar
00059   ,const EOpTransp &transp
00060   ,const Teuchos::RCP<const LinearOpBase<Scalar> > &Op
00061   )
00062 {
00063   initializeImpl(scalar,transp,Op,false);
00064 }
00065 
00066 
00067 template<class Scalar>
00068 Teuchos::RCP<LinearOpBase<Scalar> >
00069 DefaultScaledAdjointLinearOp<Scalar>::getNonconstOp()
00070 {
00071   return getOpImpl().getNonconstObj();
00072 }
00073 
00074 
00075 template<class Scalar>
00076 Teuchos::RCP<const LinearOpBase<Scalar> >
00077 DefaultScaledAdjointLinearOp<Scalar>::getOp() const
00078 {
00079   return getOpImpl();
00080 }
00081 
00082 
00083 template<class Scalar>
00084 void DefaultScaledAdjointLinearOp<Scalar>::uninitialize()
00085 {
00086   typedef Teuchos::ScalarTraits<Scalar> ST;
00087   origOp_.uninitialize();
00088   overallScalar_ = ST::zero();
00089   overallTransp_ = NOTRANS;
00090   allScalarETransp_ = Teuchos::null;
00091 
00092 }
00093 
00094 
00095 // Overridden from Teuchos::Describable
00096 
00097  
00098 template<class Scalar>
00099 std::string DefaultScaledAdjointLinearOp<Scalar>::description() const
00100 {
00101   assertInitialized();
00102   typedef Teuchos::ScalarTraits<Scalar> ST;
00103   std::ostringstream oss;
00104   oss << Teuchos::Describable::description() << "{"
00105       << overallScalar() << ","<<toString(overallTransp())<<","
00106       << origOp_->description() << "}";
00107   return oss.str();
00108 }
00109 
00110 
00111 template<class Scalar>
00112 void DefaultScaledAdjointLinearOp<Scalar>::describe(
00113   Teuchos::FancyOStream &out_arg
00114   ,const Teuchos::EVerbosityLevel verbLevel
00115   ) const
00116 {
00117   typedef Teuchos::ScalarTraits<Scalar> ST;
00118   using Teuchos::RCP;
00119   using Teuchos::FancyOStream;
00120   using Teuchos::OSTab;
00121   assertInitialized();
00122   RCP<FancyOStream> out = rcp(&out_arg,false);
00123   OSTab tab(out);
00124   switch(verbLevel) {
00125     case Teuchos::VERB_DEFAULT:
00126     case Teuchos::VERB_LOW:
00127       *out << this->description() << std::endl;
00128       break;
00129     case Teuchos::VERB_MEDIUM:
00130     case Teuchos::VERB_HIGH:
00131     case Teuchos::VERB_EXTREME:
00132     {
00133       *out
00134         << Teuchos::Describable::description() << "{"
00135         << "rangeDim=" << this->range()->dim()
00136         << ",domainDim=" << this->domain()->dim() << "}\n";
00137       OSTab tab2(out);
00138       *out
00139         << "overallScalar="<< overallScalar() << std::endl
00140         << "overallTransp="<<toString(overallTransp()) << std::endl
00141         << "Constituent transformations:\n";
00142       for( int i = 0; i <= my_index_; ++i ) {
00143         const ScalarETransp<Scalar> &scalar_transp = (*allScalarETransp_)[my_index_-i];
00144         OSTab tab3(out,i+1);
00145         if(scalar_transp.scalar != ST::one() && scalar_transp.transp != NOTRANS)
00146           *out << "scalar="<<scalar_transp.scalar<<",transp="<<toString(scalar_transp.transp)<<std::endl;
00147         else if(scalar_transp.scalar != ST::one())
00148           *out << "scalar="<<scalar_transp.scalar<<std::endl;
00149         else if( scalar_transp.transp != NOTRANS )
00150           *out << "transp="<<toString(scalar_transp.transp)<<std::endl;
00151         else
00152           *out << "no-transformation\n";
00153       }
00154       tab.incrTab(my_index_+2);
00155       *out << "origOp = " << Teuchos::describe(*origOp_,verbLevel);
00156       break;
00157     }
00158     default:
00159       TEST_FOR_EXCEPT(true); // Should never be called!
00160   }
00161 }
00162 
00163 
00164 // Overridden from LinearOpBase
00165 
00166 
00167 template<class Scalar>
00168 Teuchos::RCP< const VectorSpaceBase<Scalar> >
00169 DefaultScaledAdjointLinearOp<Scalar>::range() const
00170 {
00171   assertInitialized();
00172   return ( real_trans(this->overallTransp())==NOTRANS
00173     ? this->getOrigOp()->range() : this->getOrigOp()->domain() );
00174 }
00175 
00176 
00177 template<class Scalar>
00178 Teuchos::RCP< const VectorSpaceBase<Scalar> >
00179 DefaultScaledAdjointLinearOp<Scalar>::domain() const
00180 {
00181   assertInitialized();
00182   return ( real_trans(this->overallTransp())==NOTRANS
00183     ? this->getOrigOp()->domain() : this->getOrigOp()->range() );
00184 }
00185 
00186 
00187 template<class Scalar>
00188 Teuchos::RCP<const LinearOpBase<Scalar> >
00189 DefaultScaledAdjointLinearOp<Scalar>::clone() const
00190 {
00191   return Teuchos::null; // Not supported yet but could be
00192 }
00193 
00194 
00195 // Overridden from ScaledAdointLinearOpBase
00196 
00197 
00198 template<class Scalar>
00199 Scalar DefaultScaledAdjointLinearOp<Scalar>::overallScalar() const
00200 {
00201   return overallScalar_;
00202 }
00203 
00204 
00205 template<class Scalar>
00206 EOpTransp DefaultScaledAdjointLinearOp<Scalar>::overallTransp() const
00207 {
00208   return overallTransp_;
00209 }
00210 
00211 
00212 template<class Scalar>
00213 Teuchos::RCP<LinearOpBase<Scalar> >
00214 DefaultScaledAdjointLinearOp<Scalar>::getNonconstOrigOp()
00215 {
00216   return origOp_.getNonconstObj();
00217 }
00218 
00219 
00220 template<class Scalar>
00221 Teuchos::RCP<const LinearOpBase<Scalar> >
00222 DefaultScaledAdjointLinearOp<Scalar>::getOrigOp() const
00223 {
00224   return origOp_;
00225 }
00226 
00227 
00228 // protected
00229 
00230 
00231 // Overridden from LinearOpBase
00232 
00233 
00234 template<class Scalar>
00235 bool DefaultScaledAdjointLinearOp<Scalar>::opSupportedImpl(EOpTransp M_trans) const
00236 {
00237   assertInitialized();
00238   return Thyra::opSupported(
00239     *this->getOrigOp(), trans_trans(this->overallTransp(), M_trans) );
00240 }
00241 
00242 
00243 template<class Scalar>
00244 void DefaultScaledAdjointLinearOp<Scalar>::applyImpl(
00245   const EOpTransp M_trans,
00246   const MultiVectorBase<Scalar> &X,
00247   const Ptr<MultiVectorBase<Scalar> > &Y,
00248   const Scalar alpha,
00249   const Scalar beta
00250   ) const
00251 {
00252   using Teuchos::as;
00253   assertInitialized();
00254   Thyra::apply(
00255     *this->getOrigOp(), trans_trans(M_trans, this->overallTransp()),
00256     X, Y, as<Scalar>(this->overallScalar()*alpha), beta
00257     );
00258 }
00259 
00260 
00261 // private
00262 
00263 
00264 template<class Scalar>
00265 void DefaultScaledAdjointLinearOp<Scalar>::initializeImpl(
00266   const Scalar &scalar
00267   ,const EOpTransp &transp
00268   ,const Teuchos::RCP<const LinearOpBase<Scalar> > &Op
00269   ,const bool isConst
00270   )
00271 {
00272   typedef Teuchos::ScalarTraits<Scalar> ST;
00273 #ifdef TEUCHOS_DEBUG
00274   TEST_FOR_EXCEPTION(
00275     Op.get()==NULL, std::invalid_argument
00276     ,"DefaultScaledAdjointLinearOp<"<<ST::name()<<">::initialize(scalar,transp,Op): "
00277     "Error!, Op.get()==NULL is not allowed!"
00278     );
00279 #endif // TEUCHOS_DEBUG
00280   Teuchos::RCP<const DefaultScaledAdjointLinearOp<Scalar> >
00281     saOp = Teuchos::rcp_dynamic_cast<const DefaultScaledAdjointLinearOp<Scalar> >(Op);
00282   if(saOp.get()) {
00283     origOp_ = saOp->origOp_;
00284     overallScalar_ = saOp->overallScalar_*scalar;
00285     overallTransp_ = trans_trans(saOp->overallTransp_,transp) ;
00286     my_index_ = saOp->my_index_ + 1;
00287     allScalarETransp_ = saOp->allScalarETransp_;
00288   }
00289   else {
00290     if(isConst)
00291       origOp_.initialize(Op);
00292     else
00293       origOp_.initialize(Teuchos::rcp_const_cast<LinearOpBase<Scalar> >(Op));
00294     overallScalar_ = scalar;
00295     overallTransp_ = transp;
00296     my_index_ = 0;
00297     allScalarETransp_ = Teuchos::rcp(new allScalarETransp_t());
00298   }
00299   allScalarETransp_->push_back(ScalarETransp<Scalar>(scalar,transp));
00300   // Set the object label
00301   std::string Op_label = Op->getObjectLabel();
00302   if(Op_label.length()==0)
00303     Op_label = "ANYM";
00304   std::ostringstream label;
00305   if(scalar!=ST::one())
00306     label << scalar << "*";
00307   switch(transp) {
00308     case NOTRANS:
00309       break; // No output
00310     case CONJ:
00311       label << "conj";
00312       break;
00313     case TRANS:
00314       label << "trans";
00315       break;
00316     case CONJTRANS:
00317       label << "adj";
00318       break;
00319     default:
00320       TEST_FOR_EXCEPT("Invalid EOpTransp value!");
00321   }
00322   label << "(" << Op_label << ")";
00323   this->setObjectLabel(label.str());
00324 }
00325 
00326 
00327 template<class Scalar>
00328 typename DefaultScaledAdjointLinearOp<Scalar>::CNLOC
00329 DefaultScaledAdjointLinearOp<Scalar>::getOpImpl() const
00330 {
00331   assertInitialized();
00332   if( my_index_ > 0 ) {
00333     const ScalarETransp<Scalar> &scalar_transp = allScalarETransp_->at(my_index_);
00334     Teuchos::RCP<DefaultScaledAdjointLinearOp<Scalar> >
00335       Op = Teuchos::rcp(new DefaultScaledAdjointLinearOp<Scalar>());
00336     Op->origOp_ = origOp_;
00337     Op->overallScalar_ = overallScalar_/scalar_transp.scalar;
00338     Op->overallTransp_ = trans_trans(overallTransp_,scalar_transp.transp);
00339     Op->my_index_ = my_index_ - 1;
00340     Op->allScalarETransp_ = allScalarETransp_;
00341     return CNLOC(
00342       Teuchos::rcp_implicit_cast<LinearOpBase<Scalar> >(Op)
00343       );
00344   }
00345   return origOp_;
00346 }
00347 
00348 
00349 } // namespace Thyra
00350 
00351 
00352 #endif  // THYRA_DEFAULT_SCALED_ADJOINT_LINEAR_OP_DEF_HPP
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines