Thyra Version of the Day
Thyra_DefaultBlockedTriangularLinearOpWithSolveFactory_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 #ifndef THYRA_DEFAULT_BLOCKED_TRIANGULAR_LINEAR_OP_WITH_SOLVE_FACTORY_HPP
00043 #define THYRA_DEFAULT_BLOCKED_TRIANGULAR_LINEAR_OP_WITH_SOLVE_FACTORY_HPP
00044 
00045 
00046 #include "Thyra_DefaultBlockedTriangularLinearOpWithSolveFactory_decl.hpp"
00047 #include "Thyra_LinearOpWithSolveBase.hpp"
00048 #include "Thyra_LinearOpWithSolveFactoryHelpers.hpp"
00049 #include "Thyra_PhysicallyBlockedLinearOpBase.hpp"
00050 #include "Thyra_PhysicallyBlockedLinearOpWithSolveBase.hpp" // Interface
00051 #include "Thyra_DefaultBlockedTriangularLinearOpWithSolve.hpp" // Implementation
00052 #include "Thyra_DefaultLinearOpSource.hpp"
00053 
00054 
00055 namespace Thyra {
00056 
00057 
00058 // Overridden from Constructors/Initializers/Accessors
00059 
00060   
00061 template<class Scalar>
00062 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::DefaultBlockedTriangularLinearOpWithSolveFactory(
00063   const RCP<LinearOpWithSolveFactoryBase<Scalar> > &lowsf
00064   )
00065 {
00066 #ifdef TEUCHOS_DEBUG
00067   TEUCHOS_TEST_FOR_EXCEPT(is_null(lowsf));
00068 #endif
00069   lowsf_.initialize(lowsf);
00070 }
00071 
00072   
00073 template<class Scalar>
00074 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::DefaultBlockedTriangularLinearOpWithSolveFactory(
00075   const RCP<const LinearOpWithSolveFactoryBase<Scalar> > &lowsf
00076   )
00077 {
00078 #ifdef TEUCHOS_DEBUG
00079   TEUCHOS_TEST_FOR_EXCEPT(is_null(lowsf));
00080 #endif
00081   lowsf_.initialize(lowsf);
00082 }
00083 
00084 
00085 template<class Scalar>
00086 RCP<LinearOpWithSolveFactoryBase<Scalar> >
00087 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::getUnderlyingLOWSF()
00088 {
00089   return lowsf_.getNonconstObj();
00090 }
00091 
00092 
00093 template<class Scalar>
00094 RCP<const LinearOpWithSolveFactoryBase<Scalar> >
00095 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::getUnderlyingLOWSF() const
00096 {
00097   return lowsf_.getConstObj();
00098 }
00099 
00100 
00101 // Overridden from Teuchos::Describable
00102 
00103 
00104 template<class Scalar>
00105 std::string
00106 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::description() const
00107 {
00108   std::ostringstream oss;
00109   oss << this->Teuchos::Describable::description()
00110       << "{"
00111       << "lowsf=";
00112   if (!is_null(lowsf_.getConstObj()))
00113     oss << lowsf_.getConstObj()->description();
00114   else
00115     oss << "NULL";
00116   oss << "}";
00117   return oss.str();
00118 }
00119 
00120 
00121 // Overridden from ParameterListAcceptor
00122 
00123 
00124 template<class Scalar>
00125 void
00126 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::setParameterList(
00127   RCP<ParameterList> const& paramList
00128   )
00129 {
00130   lowsf_.getNonconstObj()->setParameterList(paramList);
00131 }
00132 
00133 
00134 template<class Scalar>
00135 RCP<ParameterList>
00136 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::getNonconstParameterList()
00137 {
00138   return lowsf_.getNonconstObj()->getNonconstParameterList();
00139 }
00140 
00141 
00142 template<class Scalar>
00143 RCP<ParameterList> 
00144 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::unsetParameterList()
00145 {
00146   return lowsf_.getNonconstObj()->unsetParameterList();
00147 }
00148 
00149 
00150 template<class Scalar>
00151 RCP<const ParameterList>
00152 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::getParameterList() const
00153 {
00154   return lowsf_.getConstObj()->getParameterList();
00155 }
00156 
00157 
00158 template<class Scalar>
00159 RCP<const ParameterList>
00160 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::getValidParameters() const
00161 {
00162   return lowsf_.getConstObj()->getValidParameters();
00163 }
00164 
00165 
00166 // Overridden from LinearOpWithSolveFactoyBase
00167 
00168 
00169 template<class Scalar>
00170 bool
00171 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::acceptsPreconditionerFactory() const
00172 {
00173   return false;
00174 }
00175 
00176 
00177 template<class Scalar>
00178 void
00179 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::setPreconditionerFactory(
00180   const RCP<PreconditionerFactoryBase<Scalar> > &precFactory,
00181   const std::string &precFactoryName
00182   )
00183 {
00184   TEUCHOS_TEST_FOR_EXCEPTION(true,std::logic_error,
00185     "Error, we don't support a preconditioner factory!");
00186 }
00187 
00188 
00189 template<class Scalar>
00190 RCP<PreconditionerFactoryBase<Scalar> >
00191 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::getPreconditionerFactory() const
00192 {
00193   return Teuchos::null;
00194 }
00195 
00196 
00197 template<class Scalar>
00198 void DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::unsetPreconditionerFactory(
00199   RCP<PreconditionerFactoryBase<Scalar> > *precFactory,
00200   std::string *precFactoryName
00201   )
00202 {
00203   TEUCHOS_TEST_FOR_EXCEPTION(true,std::logic_error,
00204     "Error, we don't support a preconditioner factory!");
00205 }
00206 
00207 
00208 template<class Scalar>
00209 bool
00210 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::isCompatible(
00211   const LinearOpSourceBase<Scalar> &fwdOpSrc
00212   ) const
00213 {
00214   TEUCHOS_TEST_FOR_EXCEPT(true);
00215   return false;
00216 }
00217 
00218 
00219 template<class Scalar>
00220 RCP<LinearOpWithSolveBase<Scalar> >
00221 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::createOp() const
00222 {
00223   return defaultBlockedTriangularLinearOpWithSolve<Scalar>();
00224 }
00225 
00226 
00227 template<class Scalar>
00228 void
00229 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::initializeOp(
00230   const RCP<const LinearOpSourceBase<Scalar> > &fwdOpSrc,
00231   LinearOpWithSolveBase<Scalar> *Op,
00232   const ESupportSolveUse supportSolveUse
00233   ) const
00234 {
00235 
00236   using Teuchos::dyn_cast;
00237   using Teuchos::rcp_dynamic_cast;
00238 
00239 #ifdef TEUCHOS_DEBUG
00240   TEUCHOS_TEST_FOR_EXCEPT(0==Op);
00241 #endif
00242 
00243   // Set the verbosity settings for the wrapped LOWSF object!
00244   lowsf_.getConstObj()->setOStream(this->getOStream());
00245   lowsf_.getConstObj()->setVerbLevel(this->getVerbLevel());
00246 
00247   // Get the block interface to get at the blocks
00248   typedef PhysicallyBlockedLinearOpBase<Scalar> PBLOB;
00249   const RCP<const PBLOB> blo =
00250     rcp_dynamic_cast<const PBLOB>(fwdOpSrc->getOp().assert_not_null());
00251   
00252   // Dynamic cast to get the DefaultBlockedTriangularLinearOpWithSolveBase
00253   // interface that we will fill.
00254 
00255   typedef DefaultBlockedTriangularLinearOpWithSolve<Scalar> DBTLOWS;
00256   DBTLOWS &btlows = dyn_cast<DBTLOWS>(*Op);
00257 
00258   // Determine if this is the first time through or if we have already
00259   // initialized before.  This will be needed to allow efficient reuse of the
00260   // LOWSB objects for the diagonal blocks.
00261   const bool firstTime = is_null(btlows.range());
00262 
00263   // If this is the first time through, we need to fill and create the block
00264   // structure
00265   if (firstTime)
00266     btlows.beginBlockFill(blo->productRange(),blo->productDomain());
00267 
00268   const int N = blo->productRange()->numBlocks();
00269   for ( int k = 0; k < N; ++k ) {
00270     const RCP<const LinearOpBase<Scalar> > fwdOp_k =
00271       blo->getBlock(k,k).assert_not_null();
00272     if (firstTime) {
00273       // This is the first time through so reate and initialize a new LOWSB
00274       // object for each block
00275       btlows.setNonconstLOWSBlock( k, k,
00276         linearOpWithSolve<Scalar>(*lowsf_.getConstObj(),fwdOp_k)
00277         );
00278     }
00279     else {
00280       // This is not the first time through so we need to just reinitiallize
00281       // the object that is already created.  This allows us to efficiently
00282       // reuse precreated structure and storage.
00283       RCP<LinearOpWithSolveBase<Scalar> >
00284         invOp_k = btlows.getNonconstLOWSBlock(k,k).assert_not_null();
00285       Thyra::initializeOp<Scalar>(*lowsf_.getConstObj(), fwdOp_k, invOp_k.ptr());
00286     }
00287   }
00288 
00289   // If this is the first time through, then we need to finalize the block
00290   // structure.
00291   if (firstTime)
00292     btlows.endBlockFill();
00293 
00294   // After the block structure has been setup, set the off-diagonal blocks.
00295   // Note that this also sets the diagonal blocks but these are ignored since
00296   // the LOWSB blocks created above override these.
00297   btlows.setBlocks(blo);
00298 
00299   // Set the verbosity settings
00300   btlows.setOStream(this->getOStream());
00301   btlows.setVerbLevel(this->getVerbLevel());
00302 
00303 }
00304 
00305 
00306 template<class Scalar>
00307 void
00308 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::initializeAndReuseOp(
00309   const RCP<const LinearOpSourceBase<Scalar> > &fwdOpSrc,
00310   LinearOpWithSolveBase<Scalar> *Op
00311   ) const
00312 {
00313   TEUCHOS_TEST_FOR_EXCEPT(true);
00314 }
00315 
00316 
00317 template<class Scalar>
00318 void
00319 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::uninitializeOp(
00320   LinearOpWithSolveBase<Scalar> *Op,
00321   RCP<const LinearOpSourceBase<Scalar> > *fwdOpSrc,
00322   RCP<const PreconditionerBase<Scalar> > *prec,
00323   RCP<const LinearOpSourceBase<Scalar> > *approxFwdOpSrc,
00324   ESupportSolveUse *supportSolveUse
00325   ) const
00326 {
00327   using Teuchos::dyn_cast;
00328   using Teuchos::rcp_implicit_cast;
00329   using Teuchos::rcp_dynamic_cast;
00330   typedef DefaultBlockedTriangularLinearOpWithSolve<Scalar> DBTLOWS;
00331   TEUCHOS_TEST_FOR_EXCEPT(0==Op);
00332   DBTLOWS &btlowsOp = dyn_cast<DBTLOWS>(*Op);
00333   if (fwdOpSrc) {
00334     const RCP<const LinearOpBase<Scalar> > fwdOp = btlowsOp.getBlocks();
00335     if (!is_null(fwdOp))
00336       *fwdOpSrc = defaultLinearOpSource<Scalar>(fwdOp);
00337     else
00338       *fwdOpSrc = Teuchos::null;
00339   }
00340   if (prec) *prec = Teuchos::null;
00341   if (approxFwdOpSrc) *approxFwdOpSrc = Teuchos::null;
00342 }
00343 
00344 
00345 template<class Scalar>
00346 bool
00347 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::supportsPreconditionerInputType(
00348   const EPreconditionerInputType precOpType
00349   ) const
00350 {
00351   // We don't support any external preconditioners!
00352   return false;
00353   // 20071006: rabartl: Note: We could support external preconditioners but it
00354   // will take some work.  We would have to extract out the individual
00355   // preconditioners from each block.  This would be pretty easy to do but I
00356   // am not going to do this until we have to.
00357 }
00358 
00359 
00360 template<class Scalar>
00361 void
00362 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::initializePreconditionedOp(
00363   const RCP<const LinearOpSourceBase<Scalar> > &fwdOpSrc,
00364   const RCP<const PreconditionerBase<Scalar> > &prec,
00365   LinearOpWithSolveBase<Scalar> *Op,
00366   const ESupportSolveUse supportSolveUse
00367   ) const
00368 {
00369   TEUCHOS_TEST_FOR_EXCEPTION(true,std::logic_error,
00370     "Error, we don't support an external preconditioner!");
00371 }
00372 
00373 
00374 template<class Scalar>
00375 void
00376 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::initializeApproxPreconditionedOp(
00377   const RCP<const LinearOpSourceBase<Scalar> > &fwdOpSrc,
00378   const RCP<const LinearOpSourceBase<Scalar> > &approxFwdOpSrc,
00379   LinearOpWithSolveBase<Scalar> *Op,
00380   const ESupportSolveUse supportSolveUse
00381   ) const
00382 {
00383   TEUCHOS_TEST_FOR_EXCEPTION(true,std::logic_error,
00384     "Error, we don't support an external preconditioner!");
00385 }
00386 
00387 
00388 // protected
00389 
00390 
00391 template<class Scalar>
00392 void
00393 DefaultBlockedTriangularLinearOpWithSolveFactory<Scalar>::informUpdatedVerbosityState() const
00394 {
00395   lowsf_.getConstObj()->setVerbLevel(this->getVerbLevel());
00396   lowsf_.getConstObj()->setOStream(this->getOStream());
00397 }
00398 
00399 
00400 } // namespace Thyra
00401 
00402 
00403 #endif // THYRA_DEFAULT_BLOCKED_TRIANGULAR_LINEAR_OP_WITH_SOLVE_FACTORY_HPP
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines