Thyra_BelosLinearOpWithSolveFactory.hpp

00001 
00002 #ifndef THYRA_BELOS_LINEAR_OP_WITH_SOLVE_FACTORY_HPP
00003 #define THYRA_BELOS_LINEAR_OP_WITH_SOLVE_FACTORY_HPP
00004 
00005 #include "Thyra_BelosLinearOpWithSolveFactoryDecl.hpp"
00006 #include "Thyra_BelosLinearOpWithSolve.hpp"
00007 #include "Thyra_ScaledAdjointLinearOpBase.hpp"
00008 #include "BelosBlockGmresSolMgr.hpp"
00009 #include "BelosPseudoBlockGmresSolMgr.hpp"
00010 #include "BelosBlockCGSolMgr.hpp"
00011 #include "BelosThyraAdapter.hpp"
00012 #include "Teuchos_VerboseObjectParameterListHelpers.hpp"
00013 #include "Teuchos_StandardParameterEntryValidators.hpp"
00014 #include "Teuchos_ParameterList.hpp"
00015 #include "Teuchos_dyn_cast.hpp"
00016 
00017 namespace Thyra {
00018 
00019 // Parameter names for Paramter List
00020 
00021 template<class Scalar>
00022 const std::string BelosLinearOpWithSolveFactory<Scalar>::SolverType_name = "Solver Type";
00023 template<class Scalar>
00024 const std::string BelosLinearOpWithSolveFactory<Scalar>::SolverType_default = "Block GMRES";
00025 template<class Scalar>
00026 const std::string BelosLinearOpWithSolveFactory<Scalar>::SolverTypes_name = "Solver Types";
00027 template<class Scalar>
00028 const std::string BelosLinearOpWithSolveFactory<Scalar>::BlockGMRES_name = "Block GMRES";
00029 template<class Scalar>
00030 const std::string  BelosLinearOpWithSolveFactory<Scalar>::PseudoBlockGMRES_name = "Pseudo Block GMRES";
00031 template<class Scalar>
00032 const std::string  BelosLinearOpWithSolveFactory<Scalar>::BlockCG_name = "Block CG";
00033 
00034 // Constructors/initializers/accessors
00035 
00036 template<class Scalar>
00037 BelosLinearOpWithSolveFactory<Scalar>::BelosLinearOpWithSolveFactory()
00038   :solverType_(SOLVER_TYPE_PSEUDO_BLOCK_GMRES)
00039 {
00040   updateThisValidParamList();
00041 }
00042 
00043 template<class Scalar>
00044 BelosLinearOpWithSolveFactory<Scalar>::BelosLinearOpWithSolveFactory(
00045   const Teuchos::RCP<PreconditionerFactoryBase<Scalar> >  &precFactory
00046   )
00047   :solverType_(SOLVER_TYPE_PSEUDO_BLOCK_GMRES)
00048 {
00049   this->setPreconditionerFactory(precFactory);
00050 }
00051 
00052 // Overridden from LinearOpWithSolveFactoryBase
00053 
00054 template<class Scalar>
00055 bool BelosLinearOpWithSolveFactory<Scalar>::acceptsPreconditionerFactory() const
00056 {
00057   return true;
00058 }
00059 
00060 template<class Scalar>
00061 void BelosLinearOpWithSolveFactory<Scalar>::setPreconditionerFactory(
00062   const Teuchos::RCP<PreconditionerFactoryBase<Scalar> > &precFactory,
00063   const std::string &precFactoryName
00064   )
00065 {
00066   TEST_FOR_EXCEPT(!precFactory.get());
00067   Teuchos::RCP<const Teuchos::ParameterList>
00068     precFactoryValidPL = precFactory->getValidParameters();
00069   const std::string _precFactoryName =
00070     ( precFactoryName != ""
00071       ? precFactoryName
00072       : ( precFactoryValidPL.get() ? precFactoryValidPL->name() : "GENERIC PRECONDITIONER FACTORY" )
00073       );
00074   precFactory_ = precFactory;
00075   precFactoryName_ = _precFactoryName;
00076   updateThisValidParamList();
00077 }
00078 
00079 template<class Scalar>
00080 Teuchos::RCP<PreconditionerFactoryBase<Scalar> >
00081 BelosLinearOpWithSolveFactory<Scalar>::getPreconditionerFactory() const
00082 {
00083   return precFactory_;
00084 }
00085 
00086 template<class Scalar>
00087 void BelosLinearOpWithSolveFactory<Scalar>::unsetPreconditionerFactory(
00088   Teuchos::RCP<PreconditionerFactoryBase<Scalar> >  *precFactory
00089   ,std::string                                              *precFactoryName
00090   )
00091 {
00092   if(precFactory) *precFactory = precFactory_;
00093   if(precFactoryName) *precFactoryName = precFactoryName_;
00094   precFactory_ = Teuchos::null;
00095   precFactoryName_ = "";
00096   updateThisValidParamList();
00097 }
00098 
00099 template<class Scalar>
00100 bool BelosLinearOpWithSolveFactory<Scalar>::isCompatible(
00101   const LinearOpSourceBase<Scalar> &fwdOpSrc
00102   ) const
00103 {
00104   if(precFactory_.get())
00105     return precFactory_->isCompatible(fwdOpSrc);
00106   return true; // Without a preconditioner, we are compatible with all linear operators!
00107 }
00108 
00109 template<class Scalar>
00110 Teuchos::RCP<LinearOpWithSolveBase<Scalar> >
00111 BelosLinearOpWithSolveFactory<Scalar>::createOp() const
00112 {
00113   return Teuchos::rcp(new BelosLinearOpWithSolve<Scalar>());
00114 }
00115 
00116 template<class Scalar>
00117 void BelosLinearOpWithSolveFactory<Scalar>::initializeOp(
00118   const Teuchos::RCP<const LinearOpSourceBase<Scalar> >    &fwdOpSrc
00119   ,LinearOpWithSolveBase<Scalar>                                   *Op
00120   ,const ESupportSolveUse                                          supportSolveUse
00121   ) const
00122 {
00123   using Teuchos::null;
00124   initializeOpImpl(fwdOpSrc,null,null,false,Op,supportSolveUse);
00125 }
00126 
00127 template<class Scalar>
00128 void BelosLinearOpWithSolveFactory<Scalar>::initializeAndReuseOp(
00129   const Teuchos::RCP<const LinearOpSourceBase<Scalar> >    &fwdOpSrc
00130   ,LinearOpWithSolveBase<Scalar>                                   *Op
00131   ) const
00132 {
00133   using Teuchos::null;
00134   initializeOpImpl(fwdOpSrc,null,null,true,Op,SUPPORT_SOLVE_UNSPECIFIED);
00135 }
00136 
00137 template<class Scalar>
00138 bool BelosLinearOpWithSolveFactory<Scalar>::supportsPreconditionerInputType(
00139   const EPreconditionerInputType precOpType
00140   ) const
00141 {
00142   if(precFactory_.get())
00143     return true;
00144   return (precOpType==PRECONDITIONER_INPUT_TYPE_AS_OPERATOR);
00145 }
00146 
00147 template<class Scalar>
00148 void BelosLinearOpWithSolveFactory<Scalar>::initializePreconditionedOp(
00149   const Teuchos::RCP<const LinearOpSourceBase<Scalar> >       &fwdOpSrc
00150   ,const Teuchos::RCP<const PreconditionerBase<Scalar> >      &prec
00151   ,LinearOpWithSolveBase<Scalar>                                      *Op
00152   ,const ESupportSolveUse                                             supportSolveUse
00153   ) const
00154 {
00155   using Teuchos::null;
00156   initializeOpImpl(fwdOpSrc,null,prec,false,Op,supportSolveUse);
00157 }
00158 
00159 template<class Scalar>
00160 void BelosLinearOpWithSolveFactory<Scalar>::initializeApproxPreconditionedOp(
00161   const Teuchos::RCP<const LinearOpSourceBase<Scalar> >      &fwdOpSrc
00162   ,const Teuchos::RCP<const LinearOpSourceBase<Scalar> >     &approxFwdOpSrc
00163   ,LinearOpWithSolveBase<Scalar>                                      *Op
00164   ,const ESupportSolveUse                                             supportSolveUse
00165   ) const
00166 {
00167   using Teuchos::null;
00168   initializeOpImpl(fwdOpSrc,approxFwdOpSrc,null,false,Op,supportSolveUse);
00169 }
00170 
00171 template<class Scalar>
00172 void BelosLinearOpWithSolveFactory<Scalar>::uninitializeOp(
00173   LinearOpWithSolveBase<Scalar>                               *Op
00174   ,Teuchos::RCP<const LinearOpSourceBase<Scalar> >    *fwdOpSrc
00175   ,Teuchos::RCP<const PreconditionerBase<Scalar> >    *prec
00176   ,Teuchos::RCP<const LinearOpSourceBase<Scalar> >    *approxFwdOpSrc
00177   ,ESupportSolveUse                                           *supportSolveUse
00178   ) const
00179 {
00180 #ifdef TEUCHOS_DEBUG
00181   TEST_FOR_EXCEPT(Op==NULL);
00182 #endif
00183   BelosLinearOpWithSolve<Scalar>
00184     &belosOp = Teuchos::dyn_cast<BelosLinearOpWithSolve<Scalar> >(*Op);
00185   Teuchos::RCP<const LinearOpSourceBase<Scalar> > 
00186     _fwdOpSrc = belosOp.extract_fwdOpSrc();
00187   Teuchos::RCP<const PreconditionerBase<Scalar> >
00188     _prec = ( belosOp.isExternalPrec() ? belosOp.extract_prec() : Teuchos::null );
00189   // Note: above we only extract the preconditioner if it was passed in
00190   // externally.  Otherwise, we need to hold on to it so that we can reuse it
00191   // in the next initialization.
00192   Teuchos::RCP<const LinearOpSourceBase<Scalar> >
00193     _approxFwdOpSrc = belosOp.extract_approxFwdOpSrc();
00194   ESupportSolveUse
00195     _supportSolveUse = belosOp.supportSolveUse();
00196   if(fwdOpSrc) *fwdOpSrc = _fwdOpSrc;
00197   if(prec) *prec = _prec;
00198   if(approxFwdOpSrc) *approxFwdOpSrc = _approxFwdOpSrc;
00199   if(supportSolveUse) *supportSolveUse = _supportSolveUse;
00200 }
00201 
00202 // Overridden from ParameterListAcceptor
00203 
00204 template<class Scalar>
00205 void BelosLinearOpWithSolveFactory<Scalar>::setParameterList(
00206   Teuchos::RCP<Teuchos::ParameterList> const& paramList
00207   )
00208 {
00209   TEST_FOR_EXCEPT(paramList.get()==NULL);
00210   paramList->validateParametersAndSetDefaults(*this->getValidParameters(),1);
00211   paramList_ = paramList;
00212   solverType_ = Teuchos::getIntegralValue<ESolverType>(*paramList_,SolverType_name);
00213   Teuchos::readVerboseObjectSublist(&*paramList_,this);
00214 }
00215 
00216 template<class Scalar>
00217 Teuchos::RCP<Teuchos::ParameterList>
00218 BelosLinearOpWithSolveFactory<Scalar>::getNonconstParameterList()
00219 {
00220   return paramList_;
00221 }
00222 
00223 template<class Scalar>
00224 Teuchos::RCP<Teuchos::ParameterList>
00225 BelosLinearOpWithSolveFactory<Scalar>::unsetParameterList()
00226 {
00227   Teuchos::RCP<Teuchos::ParameterList> _paramList = paramList_;
00228   paramList_ = Teuchos::null;
00229   return _paramList;
00230 }
00231 
00232 template<class Scalar>
00233 Teuchos::RCP<const Teuchos::ParameterList>
00234 BelosLinearOpWithSolveFactory<Scalar>::getParameterList() const
00235 {
00236   return paramList_;
00237 }
00238 
00239 template<class Scalar>
00240 Teuchos::RCP<const Teuchos::ParameterList>
00241 BelosLinearOpWithSolveFactory<Scalar>::getValidParameters() const
00242 {
00243   return thisValidParamList_;
00244 }
00245 
00246 // Public functions overridden from Teuchos::Describable
00247 
00248 template<class Scalar>
00249 std::string BelosLinearOpWithSolveFactory<Scalar>::description() const
00250 {
00251   std::ostringstream oss;
00252   oss << "Thyra::BelosLinearOpWithSolveFactory";
00253   //oss << "{";
00254   // ToDo: Fill this in some!
00255   //oss << "}";
00256   return oss.str();
00257 }
00258 
00259 // private
00260 
00261 template<class Scalar>
00262 Teuchos::RCP<const Teuchos::ParameterList>
00263 BelosLinearOpWithSolveFactory<Scalar>::generateAndGetValidParameters()
00264 {
00265   using Teuchos::tuple;
00266   using Teuchos::setStringToIntegralParameter;
00267   typedef MultiVectorBase<Scalar> MV_t;
00268   typedef LinearOpBase<Scalar> LO_t;
00269   static Teuchos::RCP<Teuchos::ParameterList> validParamList;
00270   if(validParamList.get()==NULL) {
00271     validParamList = Teuchos::rcp(new Teuchos::ParameterList("BelosLinearOpWithSolveFactory"));
00272     setStringToIntegralParameter<ESolverType>(
00273       SolverType_name, SolverType_default,
00274       "Type of linear solver algorithm to use.",
00275       tuple<std::string>(
00276         "Block GMRES",
00277         "Pseudo Block GMRES",
00278         "Block CG",
00279         "GCRODR"
00280         ),
00281       tuple<std::string>(
00282         "Performs block and single single-RHS GMRES as well as\n"
00283         "flexible GMRES by setting options in the \"Block GMRES\" sublist.",
00284 
00285         "GMRES solver that performs single-RHS GMRES on multiple RHSs taking\n"
00286         "advantage of operator multi-vector multiplication and the amortization\n"
00287         "of global communication.  Individual linear systems are deflated out as\n"
00288         "they are solved.",
00289 
00290         "CG solver that performs block and single-RHS CG.",
00291 
00292         "GMRES solver that performs subspace recycling between RHS and linear systems."
00293         ),
00294       tuple<ESolverType>(
00295         SOLVER_TYPE_BLOCK_GMRES,
00296         SOLVER_TYPE_PSEUDO_BLOCK_GMRES,
00297         SOLVER_TYPE_BLOCK_CG,
00298         SOLVER_TYPE_GCRODR
00299         ),
00300       &*validParamList
00301       );
00302     Teuchos::ParameterList
00303       &solverTypesSL = validParamList->sublist(SolverTypes_name);
00304     {
00305       Belos::BlockGmresSolMgr<Scalar,MV_t,LO_t> mgr;
00306       solverTypesSL.sublist(BlockGMRES_name).setParameters(
00307         *mgr.getValidParameters()
00308         );
00309     }
00310     {
00311       Belos::PseudoBlockGmresSolMgr<Scalar,MV_t,LO_t> mgr;
00312       solverTypesSL.sublist(PseudoBlockGMRES_name).setParameters(
00313         *mgr.getValidParameters()
00314         );
00315     }
00316     {
00317       Belos::BlockCGSolMgr<Scalar,MV_t,LO_t> mgr;
00318       solverTypesSL.sublist(BlockCG_name).setParameters(
00319         *mgr.getValidParameters()
00320         );
00321     }
00322   }
00323   return validParamList;
00324 }
00325 
00326 template<class Scalar>
00327 void BelosLinearOpWithSolveFactory<Scalar>::updateThisValidParamList()
00328 {
00329   thisValidParamList_ = Teuchos::rcp(
00330     new Teuchos::ParameterList(*generateAndGetValidParameters())
00331     );
00332   Teuchos::setupVerboseObjectSublist(&*thisValidParamList_);
00333 }
00334 
00335 template<class Scalar>
00336 void BelosLinearOpWithSolveFactory<Scalar>::initializeOpImpl(
00337   const Teuchos::RCP<const LinearOpSourceBase<Scalar> >       &fwdOpSrc
00338   ,const Teuchos::RCP<const LinearOpSourceBase<Scalar> >      &approxFwdOpSrc
00339   ,const Teuchos::RCP<const PreconditionerBase<Scalar> >      &prec_in
00340   ,const bool                                                         reusePrec
00341   ,LinearOpWithSolveBase<Scalar>                                      *Op
00342   ,const ESupportSolveUse                                             supportSolveUse
00343   ) const
00344 {
00345 
00346   using Teuchos::RCP;
00347   using Teuchos::rcp;
00348   using Teuchos::set_extra_data;
00349   typedef Teuchos::ScalarTraits<Scalar> ST;
00350   typedef typename ST::magnitudeType ScalarMag;
00351   typedef MultiVectorBase<Scalar>    MV_t;
00352   typedef LinearOpBase<Scalar>       LO_t;
00353 
00354   const Teuchos::RCP<Teuchos::FancyOStream> out       = this->getOStream();
00355   const Teuchos::EVerbosityLevel                    verbLevel = this->getVerbLevel();
00356   Teuchos::OSTab tab(out);
00357   if(out.get() && static_cast<int>(verbLevel) >= static_cast<int>(Teuchos::VERB_LOW))
00358     *out << "\nEntering Thyra::BelosLinearOpWithSolveFactory<"<<ST::name()<<">::initializeOpImpl(...) ...\n";
00359 
00360   typedef Teuchos::VerboseObjectTempState<PreconditionerFactoryBase<Scalar> > VOTSPF;
00361   VOTSPF precFactoryOutputTempState(precFactory_,out,verbLevel);
00362   
00363   TEST_FOR_EXCEPT(Op==NULL);
00364   TEST_FOR_EXCEPT(fwdOpSrc.get()==NULL);
00365   TEST_FOR_EXCEPT(fwdOpSrc->getOp().get()==NULL);
00366   Teuchos::RCP<const LinearOpBase<Scalar> >
00367     fwdOp = fwdOpSrc->getOp(),
00368     approxFwdOp = ( approxFwdOpSrc.get() ? approxFwdOpSrc->getOp() : Teuchos::null );
00369 
00370   //
00371   // Get the BelosLinearOpWithSolve interface
00372   //
00373   BelosLinearOpWithSolve<Scalar>
00374     *belosOp = &Teuchos::dyn_cast<BelosLinearOpWithSolve<Scalar> >(*Op);
00375   //
00376   // Get/Create the preconditioner
00377   //
00378   RCP<PreconditionerBase<Scalar> >         myPrec = Teuchos::null;
00379   RCP<const PreconditionerBase<Scalar> >   prec = Teuchos::null;
00380   if(prec_in.get()) {
00381     // Use an externally defined preconditioner
00382     prec = prec_in;
00383   }
00384   else {
00385     // Try and generate a preconditioner on our own
00386     if(precFactory_.get()) {
00387       myPrec =
00388         ( !belosOp->isExternalPrec()
00389           ? Teuchos::rcp_const_cast<PreconditionerBase<Scalar> >(belosOp->extract_prec())
00390           : Teuchos::null
00391           );
00392       bool hasExistingPrec = false;
00393       if(myPrec.get()) {
00394         hasExistingPrec = true;
00395         // ToDo: Get the forward operator and validate that it is the same
00396         // operator that is used here!
00397       }
00398       else {
00399         hasExistingPrec = false;
00400         myPrec = precFactory_->createPrec();
00401       }
00402       if( hasExistingPrec && reusePrec ) {
00403         // Just reuse the existing preconditioner again!
00404       }
00405       else {
00406         // Update the preconditioner
00407         if(approxFwdOp.get())
00408           precFactory_->initializePrec(approxFwdOpSrc,&*myPrec);
00409         else
00410           precFactory_->initializePrec(fwdOpSrc,&*myPrec);
00411       }
00412       prec = myPrec;
00413     }
00414   }
00415   //
00416   // Uninitialize the current solver object
00417   //
00418   bool oldIsExternalPrec = false;
00419   RCP<Belos::LinearProblem<Scalar,MV_t,LO_t> >     oldLP = Teuchos::null;
00420   RCP<Belos::SolverManager<Scalar,MV_t,LO_t> >     oldIterSolver = Teuchos::null;
00421   RCP<const LinearOpSourceBase<Scalar> >           oldFwdOpSrc = Teuchos::null;
00422   RCP<const LinearOpSourceBase<Scalar> >           oldApproxFwdOpSrc = Teuchos::null;   
00423   ESupportSolveUse                                 oldSupportSolveUse = SUPPORT_SOLVE_UNSPECIFIED;
00424 
00425   belosOp->uninitialize( &oldLP,
00426                          NULL,
00427                          &oldIterSolver,
00428                          &oldFwdOpSrc,
00429                          NULL,
00430                          &oldIsExternalPrec,
00431                          &oldApproxFwdOpSrc,
00432                          &oldSupportSolveUse );
00433   //
00434   // Create the Belos linear problem
00435   // NOTE:  If one exists already, reuse it.
00436   typedef Belos::LinearProblem<Scalar,MV_t,LO_t> LP_t;
00437   RCP<LP_t> lp;
00438   if (oldLP != Teuchos::null) {
00439     lp = oldLP;
00440   }
00441   else {
00442     lp = rcp(new LP_t());
00443   }
00444   //
00445   // Set the operator
00446   //
00447   lp->setOperator(fwdOp);
00448   //
00449   // Set the preconditioner
00450   //
00451   if(prec.get()) {
00452     RCP<const LinearOpBase<Scalar> > unspecified = prec->getUnspecifiedPrecOp();
00453     RCP<const LinearOpBase<Scalar> > left        = prec->getLeftPrecOp();
00454     RCP<const LinearOpBase<Scalar> > right       = prec->getRightPrecOp();
00455     TEST_FOR_EXCEPTION(
00456       !( left.get() || right.get() || unspecified.get() ), std::logic_error
00457       ,"Error, at least one preconditoner linear operator objects must be set!"
00458       );
00459     if(unspecified.get()) {
00460       lp->setRightPrec(unspecified);
00461       // ToDo: Allow user to determine whether this should be placed on the
00462       // left or on the right through a parameter in the parameter list!
00463     }
00464     else {
00465       // Set a left, right or split preconditioner
00466       TEST_FOR_EXCEPTION(
00467         left.get(),std::logic_error
00468         ,"Error, we can not currently handle a left preconditioner!"
00469         );
00470       lp->setRightPrec(right);
00471     }
00472   }
00473   if(myPrec.get()) {
00474     set_extra_data<RCP<PreconditionerBase<Scalar> > >(myPrec,"Belos::InternalPrec",
00475                     &lp, Teuchos::POST_DESTROY, false);
00476   }
00477   else if(prec.get()) {
00478     set_extra_data<RCP<const PreconditionerBase<Scalar> > >(prec,"Belos::ExternalPrec",
00479                     &lp, Teuchos::POST_DESTROY, false);
00480   }
00481   //
00482   // Generate the parameter list.
00483   //
00484   typedef Belos::SolverManager<Scalar,MV_t,LO_t> IterativeSolver_t;
00485   RCP<IterativeSolver_t> iterativeSolver = Teuchos::null;
00486   RCP<Teuchos::ParameterList> solverPL = Teuchos::rcp( new Teuchos::ParameterList() );
00487   
00488   switch(solverType_) {
00489     case SOLVER_TYPE_BLOCK_GMRES: 
00490     {
00491       // Set the PL
00492       if(paramList_.get()) {
00493         Teuchos::ParameterList &solverTypesPL = paramList_->sublist(SolverTypes_name);
00494         Teuchos::ParameterList &gmresPL = solverTypesPL.sublist(BlockGMRES_name);
00495         solverPL = Teuchos::rcp( &gmresPL, false );
00496       }
00497       // Create the solver
00498       if (oldIterSolver != Teuchos::null) {
00499         iterativeSolver = oldIterSolver;
00500         iterativeSolver->setProblem( lp );
00501         iterativeSolver->setParameters( solverPL );
00502       } 
00503       else {
00504         iterativeSolver = rcp(new Belos::BlockGmresSolMgr<Scalar,MV_t,LO_t>(lp,solverPL));
00505       }
00506       break;
00507     }
00508     case SOLVER_TYPE_PSEUDO_BLOCK_GMRES:
00509     {
00510       // Set the PL
00511       if(paramList_.get()) {
00512         Teuchos::ParameterList &solverTypesPL = paramList_->sublist(SolverTypes_name);
00513         Teuchos::ParameterList &pbgmresPL = solverTypesPL.sublist(PseudoBlockGMRES_name);
00514         solverPL = Teuchos::rcp( &pbgmresPL, false );
00515       }
00516       // 
00517       // Create the solver
00518       // 
00519       if (oldIterSolver != Teuchos::null) {
00520         iterativeSolver = oldIterSolver;
00521         iterativeSolver->setProblem( lp );
00522         iterativeSolver->setParameters( solverPL );
00523       }
00524       else {
00525         iterativeSolver = rcp(new Belos::PseudoBlockGmresSolMgr<Scalar,MV_t,LO_t>(lp,solverPL));
00526       }
00527       break;
00528     }
00529     case SOLVER_TYPE_BLOCK_CG:
00530     {
00531       // Set the PL
00532       if(paramList_.get()) {
00533         Teuchos::ParameterList &solverTypesPL = paramList_->sublist(SolverTypes_name);
00534         Teuchos::ParameterList &cgPL = solverTypesPL.sublist(BlockCG_name);
00535         solverPL = Teuchos::rcp( &cgPL, false );
00536       }
00537       // Create the solver
00538       if (oldIterSolver != Teuchos::null) {
00539         iterativeSolver = oldIterSolver;
00540         iterativeSolver->setProblem( lp );
00541         iterativeSolver->setParameters( solverPL );
00542       }
00543       else {
00544         iterativeSolver = rcp(new Belos::BlockCGSolMgr<Scalar,MV_t,LO_t>(lp,solverPL));
00545       }
00546       break;
00547     }
00548     case SOLVER_TYPE_GCRODR:
00549     {
00550       TEST_FOR_EXCEPT("ToDo: Support GCRODR"); 
00551       break;
00552     }
00553     default:
00554     {
00555       TEST_FOR_EXCEPT(true);
00556     }
00557   }
00558 
00559   //
00560   // Initialize the LOWS object
00561   //
00562   belosOp->initialize(
00563     lp,solverPL,iterativeSolver
00564     ,fwdOpSrc,prec,myPrec.get()==NULL,approxFwdOpSrc,supportSolveUse
00565     );
00566   belosOp->setOStream(out);
00567   belosOp->setVerbLevel(verbLevel);
00568 #ifdef TEUCHOS_DEBUG
00569   if(paramList_.get()) {
00570     // Make sure we read the list correctly
00571     paramList_->validateParameters(*this->getValidParameters(),1); // Validate 0th and 1st level deep
00572   }
00573 #endif
00574   if(out.get() && static_cast<int>(verbLevel) >= static_cast<int>(Teuchos::VERB_LOW))
00575     *out << "\nLeaving Thyra::BelosLinearOpWithSolveFactory<"<<ST::name()<<">::initializeOpImpl(...) ...\n";
00576   
00577 }
00578 
00579 } // namespace Thyra
00580 
00581 #endif // THYRA_BELOS_LINEAR_OP_WITH_SOLVE_FACTORY_HPP

Generated on Wed Jul 22 12:56:30 2009 for Stratimikos by  doxygen 1.5.8