Thyra Version of the Day
Thyra_DefaultScaledAdjointLinearOp_def.hpp
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 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Roscoe A. Bartlett (bartlettra@ornl.gov) 
00038 // 
00039 // ***********************************************************************
00040 // @HEADER
00041 
00042 
00043 #ifndef THYRA_DEFAULT_SCALED_ADJOINT_LINEAR_OP_DEF_HPP
00044 #define THYRA_DEFAULT_SCALED_ADJOINT_LINEAR_OP_DEF_HPP
00045 
00046 
00047 #include "Thyra_DefaultScaledAdjointLinearOp_decl.hpp"
00048 #include "Thyra_VectorSpaceBase.hpp"
00049 #include "Thyra_AssertOp.hpp"
00050 
00051 
00052 namespace Thyra {
00053 
00054 
00055 //Constructors/initializers/accessors
00056 
00057 
00058 template<class Scalar>
00059 void DefaultScaledAdjointLinearOp<Scalar>::initialize(
00060   const Scalar &scalar
00061   ,const EOpTransp &transp
00062   ,const Teuchos::RCP<LinearOpBase<Scalar> > &Op
00063   )
00064 {
00065   initializeImpl(scalar,transp,Op,false);
00066 }
00067 
00068 
00069 template<class Scalar>
00070 void DefaultScaledAdjointLinearOp<Scalar>::initialize(
00071   const Scalar &scalar
00072   ,const EOpTransp &transp
00073   ,const Teuchos::RCP<const LinearOpBase<Scalar> > &Op
00074   )
00075 {
00076   initializeImpl(scalar,transp,Op,false);
00077 }
00078 
00079 
00080 template<class Scalar>
00081 Teuchos::RCP<LinearOpBase<Scalar> >
00082 DefaultScaledAdjointLinearOp<Scalar>::getNonconstOp()
00083 {
00084   return getOpImpl().getNonconstObj();
00085 }
00086 
00087 
00088 template<class Scalar>
00089 Teuchos::RCP<const LinearOpBase<Scalar> >
00090 DefaultScaledAdjointLinearOp<Scalar>::getOp() const
00091 {
00092   return getOpImpl();
00093 }
00094 
00095 
00096 template<class Scalar>
00097 void DefaultScaledAdjointLinearOp<Scalar>::uninitialize()
00098 {
00099   typedef Teuchos::ScalarTraits<Scalar> ST;
00100   origOp_.uninitialize();
00101   overallScalar_ = ST::zero();
00102   overallTransp_ = NOTRANS;
00103   allScalarETransp_ = Teuchos::null;
00104 
00105 }
00106 
00107 
00108 // Overridden from Teuchos::Describable
00109 
00110  
00111 template<class Scalar>
00112 std::string DefaultScaledAdjointLinearOp<Scalar>::description() const
00113 {
00114   assertInitialized();
00115   typedef Teuchos::ScalarTraits<Scalar> ST;
00116   std::ostringstream oss;
00117   oss << Teuchos::Describable::description() << "{"
00118       << overallScalar() << ","<<toString(overallTransp())<<","
00119       << origOp_->description() << "}";
00120   return oss.str();
00121 }
00122 
00123 
00124 template<class Scalar>
00125 void DefaultScaledAdjointLinearOp<Scalar>::describe(
00126   Teuchos::FancyOStream &out_arg
00127   ,const Teuchos::EVerbosityLevel verbLevel
00128   ) const
00129 {
00130   typedef Teuchos::ScalarTraits<Scalar> ST;
00131   using Teuchos::RCP;
00132   using Teuchos::FancyOStream;
00133   using Teuchos::OSTab;
00134   assertInitialized();
00135   RCP<FancyOStream> out = rcp(&out_arg,false);
00136   OSTab tab(out);
00137   switch(verbLevel) {
00138     case Teuchos::VERB_DEFAULT:
00139     case Teuchos::VERB_LOW:
00140       *out << this->description() << std::endl;
00141       break;
00142     case Teuchos::VERB_MEDIUM:
00143     case Teuchos::VERB_HIGH:
00144     case Teuchos::VERB_EXTREME:
00145     {
00146       *out
00147         << Teuchos::Describable::description() << "{"
00148         << "rangeDim=" << this->range()->dim()
00149         << ",domainDim=" << this->domain()->dim() << "}\n";
00150       OSTab tab2(out);
00151       *out
00152         << "overallScalar="<< overallScalar() << std::endl
00153         << "overallTransp="<<toString(overallTransp()) << std::endl
00154         << "Constituent transformations:\n";
00155       for( int i = 0; i <= my_index_; ++i ) {
00156         const ScalarETransp<Scalar> &scalar_transp = (*allScalarETransp_)[my_index_-i];
00157         OSTab tab3(out,i+1);
00158         if(scalar_transp.scalar != ST::one() && scalar_transp.transp != NOTRANS)
00159           *out << "scalar="<<scalar_transp.scalar<<",transp="<<toString(scalar_transp.transp)<<std::endl;
00160         else if(scalar_transp.scalar != ST::one())
00161           *out << "scalar="<<scalar_transp.scalar<<std::endl;
00162         else if( scalar_transp.transp != NOTRANS )
00163           *out << "transp="<<toString(scalar_transp.transp)<<std::endl;
00164         else
00165           *out << "no-transformation\n";
00166       }
00167       tab.incrTab(my_index_+2);
00168       *out << "origOp = " << Teuchos::describe(*origOp_,verbLevel);
00169       break;
00170     }
00171     default:
00172       TEUCHOS_TEST_FOR_EXCEPT(true); // Should never be called!
00173   }
00174 }
00175 
00176 
00177 // Overridden from LinearOpBase
00178 
00179 
00180 template<class Scalar>
00181 Teuchos::RCP< const VectorSpaceBase<Scalar> >
00182 DefaultScaledAdjointLinearOp<Scalar>::range() const
00183 {
00184   assertInitialized();
00185   return ( real_trans(this->overallTransp())==NOTRANS
00186     ? this->getOrigOp()->range() : this->getOrigOp()->domain() );
00187 }
00188 
00189 
00190 template<class Scalar>
00191 Teuchos::RCP< const VectorSpaceBase<Scalar> >
00192 DefaultScaledAdjointLinearOp<Scalar>::domain() const
00193 {
00194   assertInitialized();
00195   return ( real_trans(this->overallTransp())==NOTRANS
00196     ? this->getOrigOp()->domain() : this->getOrigOp()->range() );
00197 }
00198 
00199 
00200 template<class Scalar>
00201 Teuchos::RCP<const LinearOpBase<Scalar> >
00202 DefaultScaledAdjointLinearOp<Scalar>::clone() const
00203 {
00204   return Teuchos::null; // Not supported yet but could be
00205 }
00206 
00207 
00208 // Overridden from ScaledAdointLinearOpBase
00209 
00210 
00211 template<class Scalar>
00212 Scalar DefaultScaledAdjointLinearOp<Scalar>::overallScalar() const
00213 {
00214   return overallScalar_;
00215 }
00216 
00217 
00218 template<class Scalar>
00219 EOpTransp DefaultScaledAdjointLinearOp<Scalar>::overallTransp() const
00220 {
00221   return overallTransp_;
00222 }
00223 
00224 
00225 template<class Scalar>
00226 Teuchos::RCP<LinearOpBase<Scalar> >
00227 DefaultScaledAdjointLinearOp<Scalar>::getNonconstOrigOp()
00228 {
00229   return origOp_.getNonconstObj();
00230 }
00231 
00232 
00233 template<class Scalar>
00234 Teuchos::RCP<const LinearOpBase<Scalar> >
00235 DefaultScaledAdjointLinearOp<Scalar>::getOrigOp() const
00236 {
00237   return origOp_;
00238 }
00239 
00240 
00241 // protected
00242 
00243 
00244 // Overridden from LinearOpBase
00245 
00246 
00247 template<class Scalar>
00248 bool DefaultScaledAdjointLinearOp<Scalar>::opSupportedImpl(EOpTransp M_trans) const
00249 {
00250   assertInitialized();
00251   return Thyra::opSupported(
00252     *this->getOrigOp(), trans_trans(this->overallTransp(), M_trans) );
00253 }
00254 
00255 
00256 template<class Scalar>
00257 void DefaultScaledAdjointLinearOp<Scalar>::applyImpl(
00258   const EOpTransp M_trans,
00259   const MultiVectorBase<Scalar> &X,
00260   const Ptr<MultiVectorBase<Scalar> > &Y,
00261   const Scalar alpha,
00262   const Scalar beta
00263   ) const
00264 {
00265   using Teuchos::as;
00266   assertInitialized();
00267   Thyra::apply(
00268     *this->getOrigOp(), trans_trans(M_trans, this->overallTransp()),
00269     X, Y, as<Scalar>(this->overallScalar()*alpha), beta
00270     );
00271 }
00272 
00273 
00274 // private
00275 
00276 
00277 template<class Scalar>
00278 void DefaultScaledAdjointLinearOp<Scalar>::initializeImpl(
00279   const Scalar &scalar
00280   ,const EOpTransp &transp
00281   ,const Teuchos::RCP<const LinearOpBase<Scalar> > &Op
00282   ,const bool isConst
00283   )
00284 {
00285   typedef Teuchos::ScalarTraits<Scalar> ST;
00286 #ifdef TEUCHOS_DEBUG
00287   TEUCHOS_TEST_FOR_EXCEPTION(
00288     Op.get()==NULL, std::invalid_argument
00289     ,"DefaultScaledAdjointLinearOp<"<<ST::name()<<">::initialize(scalar,transp,Op): "
00290     "Error!, Op.get()==NULL is not allowed!"
00291     );
00292 #endif // TEUCHOS_DEBUG
00293   Teuchos::RCP<const DefaultScaledAdjointLinearOp<Scalar> >
00294     saOp = Teuchos::rcp_dynamic_cast<const DefaultScaledAdjointLinearOp<Scalar> >(Op);
00295   if(saOp.get()) {
00296     origOp_ = saOp->origOp_;
00297     overallScalar_ = saOp->overallScalar_*scalar;
00298     overallTransp_ = trans_trans(saOp->overallTransp_,transp) ;
00299     my_index_ = saOp->my_index_ + 1;
00300     allScalarETransp_ = saOp->allScalarETransp_;
00301   }
00302   else {
00303     if(isConst)
00304       origOp_.initialize(Op);
00305     else
00306       origOp_.initialize(Teuchos::rcp_const_cast<LinearOpBase<Scalar> >(Op));
00307     overallScalar_ = scalar;
00308     overallTransp_ = transp;
00309     my_index_ = 0;
00310     allScalarETransp_ = Teuchos::rcp(new allScalarETransp_t());
00311   }
00312   allScalarETransp_->push_back(ScalarETransp<Scalar>(scalar,transp));
00313   // Set the object label
00314   std::string Op_label = Op->getObjectLabel();
00315   if(Op_label.length()==0)
00316     Op_label = "ANYM";
00317   std::ostringstream label;
00318   if(scalar!=ST::one())
00319     label << scalar << "*";
00320   switch(transp) {
00321     case NOTRANS:
00322       break; // No output
00323     case CONJ:
00324       label << "conj";
00325       break;
00326     case TRANS:
00327       label << "trans";
00328       break;
00329     case CONJTRANS:
00330       label << "adj";
00331       break;
00332     default:
00333       TEUCHOS_TEST_FOR_EXCEPT("Invalid EOpTransp value!");
00334   }
00335   label << "(" << Op_label << ")";
00336   this->setObjectLabel(label.str());
00337 }
00338 
00339 
00340 template<class Scalar>
00341 typename DefaultScaledAdjointLinearOp<Scalar>::CNLOC
00342 DefaultScaledAdjointLinearOp<Scalar>::getOpImpl() const
00343 {
00344   assertInitialized();
00345   if( my_index_ > 0 ) {
00346     const ScalarETransp<Scalar> &scalar_transp = allScalarETransp_->at(my_index_);
00347     Teuchos::RCP<DefaultScaledAdjointLinearOp<Scalar> >
00348       Op = Teuchos::rcp(new DefaultScaledAdjointLinearOp<Scalar>());
00349     Op->origOp_ = origOp_;
00350     Op->overallScalar_ = overallScalar_/scalar_transp.scalar;
00351     Op->overallTransp_ = trans_trans(overallTransp_,scalar_transp.transp);
00352     Op->my_index_ = my_index_ - 1;
00353     Op->allScalarETransp_ = allScalarETransp_;
00354     return CNLOC(
00355       Teuchos::rcp_implicit_cast<LinearOpBase<Scalar> >(Op)
00356       );
00357   }
00358   return origOp_;
00359 }
00360 
00361 
00362 } // namespace Thyra
00363 
00364 
00365 #endif  // THYRA_DEFAULT_SCALED_ADJOINT_LINEAR_OP_DEF_HPP
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines