Teko Version of the Day
Teko_InverseLibrary.cpp
00001 /*
00002 // @HEADER
00003 // 
00004 // ***********************************************************************
00005 // 
00006 //      Teko: A package for block and physics based preconditioning
00007 //                  Copyright 2010 Sandia Corporation 
00008 //  
00009 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00010 // the U.S. Government retains certain rights in this software.
00011 //  
00012 // Redistribution and use in source and binary forms, with or without
00013 // modification, are permitted provided that the following conditions are
00014 // met:
00015 //  
00016 // 1. Redistributions of source code must retain the above copyright
00017 // notice, this list of conditions and the following disclaimer.
00018 //  
00019 // 2. Redistributions in binary form must reproduce the above copyright
00020 // notice, this list of conditions and the following disclaimer in the
00021 // documentation and/or other materials provided with the distribution.
00022 //  
00023 // 3. Neither the name of the Corporation nor the names of the
00024 // contributors may be used to endorse or promote products derived from
00025 // this software without specific prior written permission. 
00026 //  
00027 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00028 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00029 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00030 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00031 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00032 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00033 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00034 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00035 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING 
00036 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00037 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00038 //  
00039 // Questions? Contact Eric C. Cyr (eccyr@sandia.gov)
00040 // 
00041 // ***********************************************************************
00042 // 
00043 // @HEADER
00044 
00045 */
00046 
00047 #include "Teko_InverseLibrary.hpp"
00048 
00049 #include "Teko_SolveInverseFactory.hpp"
00050 #include "Teko_PreconditionerInverseFactory.hpp"
00051 #include "Teko_BlockPreconditionerFactory.hpp"
00052 
00053 #include "Teko_NeumannSeriesPreconditionerFactory.hpp"
00054 #include "Teuchos_AbstractFactoryStd.hpp"
00055 
00056 #include <algorithm>
00057 
00058 using Teuchos::RCP;
00059 using Teuchos::rcp;
00060 
00061 namespace Teko {
00062 
00066 void addToStratimikosBuilder(Stratimikos::DefaultLinearSolverBuilder & builder)
00067 {
00068    typedef Thyra::PreconditionerFactoryBase<double> PrecFactory;
00069 
00070    RCP<const Teuchos::AbstractFactory<Thyra::PreconditionerFactoryBase<double> > > factory;
00071      
00072    factory = Teuchos::abstractFactoryStd<PrecFactory,Teko::NeumannSeriesPreconditionerFactory<double> >();
00073    builder.setPreconditioningStrategyFactory(factory,"Neumann Series");
00074 }
00075 
00076 InverseLibrary::InverseLibrary()
00077 {
00078    Teko_DEBUG_SCOPE("InverseLibrary::InverseLibrary", 10);
00079 
00080    // setup some valid Stratimikos parameters
00082 
00083    // set valid solve factory names
00084    stratValidSolver_.push_back("Belos"); 
00085    stratValidSolver_.push_back("Amesos"); 
00086    stratValidSolver_.push_back("AztecOO"); 
00087 
00088    // set valid preconditioner factory name
00089    stratValidPrecond_.push_back("ML"); 
00090    stratValidPrecond_.push_back("Ifpack"); 
00091    stratValidPrecond_.push_back("Neumann Series"); 
00092 
00093    // set valid Teko preconditioner factory names
00094    PreconditionerFactory::getPreconditionerFactoryNames(blockValidPrecond_);
00095 
00096    Teko_DEBUG_MSG_BEGIN(10)
00097       DEBUG_STREAM << "Loaded \"block\" preconditioners = ";
00098       for(std::size_t i=0;i<blockValidPrecond_.size();i++)
00099          DEBUG_STREAM << blockValidPrecond_[i] << ", ";
00100       DEBUG_STREAM << std::endl;
00101    Teko_DEBUG_MSG_END()
00102 }
00103 
00105 void InverseLibrary::addInverse(const std::string & label,const Teuchos::ParameterList & pl)
00106 {
00107    // strip out the label
00108    const std::string type = pl.get<std::string>("Type");
00109 
00110    // copy the parameter list so we can modify it
00111    Teuchos::ParameterList settingsList;
00112    settingsList.set(type,pl);
00113    settingsList.sublist(type).remove("Type");
00114 
00115    // is this a Stratimikos preconditioner or solver
00116    if(std::find(stratValidPrecond_.begin(),stratValidPrecond_.end(),type)!=stratValidPrecond_.end()) {
00117       // this is a Stratimikos preconditioner factory
00118       addStratPrecond(label,type,settingsList);
00119    }
00120    else if(std::find(stratValidSolver_.begin(),stratValidSolver_.end(),type)!=stratValidSolver_.end()) {
00121       // this is a Stratimikos preconditioner factory
00122       addStratSolver(label,type,settingsList);
00123    }
00124    else if(std::find(blockValidPrecond_.begin(),blockValidPrecond_.end(),type)!=blockValidPrecond_.end()) {
00125       // this is a Teko preconditioner factory
00126       addBlockPrecond(label,type,settingsList);
00127    }
00128    else {
00129       Teuchos::FancyOStream & os = *Teko::getOutputStream();
00130       os << "ERROR: Could not find inverse type \"" << type 
00131          << "\" required by inverse name \"" << label << "\"" << std::endl;
00132       TEUCHOS_ASSERT(false);
00133    }
00134 }
00135 
00137 void InverseLibrary::addStratSolver(const std::string & label,const std::string & type,const Teuchos::ParameterList & pl)
00138 {
00139    // add some additional parameters onto the list
00140    RCP<Teuchos::ParameterList> stratList = rcp(new Teuchos::ParameterList());
00141    stratList->set("Linear Solver Type",type);
00142    stratList->set("Linear Solver Types",pl);
00143    stratList->set("Preconditioner Type","None");
00144 
00145    stratSolver_[label] = stratList;
00146 }
00147 
00149 void InverseLibrary::addStratPrecond(const std::string & label,const std::string & type,const Teuchos::ParameterList & pl)
00150 {
00151    // add some additional parameters onto the list
00152    RCP<Teuchos::ParameterList> stratList = rcp(new Teuchos::ParameterList());
00153    stratList->set("Preconditioner Type",type);
00154    stratList->set("Preconditioner Types",pl);
00155 
00156    stratPrecond_[label] = stratList;
00157 }
00158 
00160 void InverseLibrary::addBlockPrecond(const std::string & label,const std::string & type,const Teuchos::ParameterList & pl)
00161 {
00162    // add some additional parameters onto the list
00163    RCP<Teuchos::ParameterList> blockList = rcp(new Teuchos::ParameterList());
00164    blockList->set("Preconditioner Type",type);
00165    blockList->set("Preconditioner Settings",pl.sublist(type));
00166 
00167    // add the Teko preconditioner parameter list into the library
00168    blockPrecond_[label] = blockList;
00169 }
00170 
00178 Teuchos::RCP<const Teuchos::ParameterList> InverseLibrary::getParameterList(const std::string & label) const
00179 {
00180    std::map<std::string,RCP<const Teuchos::ParameterList> >::const_iterator itr;
00181    
00182    // check preconditioners
00183    itr = stratPrecond_.find(label);
00184    if(itr!=stratPrecond_.end()) return itr->second;
00185     
00186    // check solvers
00187    itr = stratSolver_.find(label);
00188    if(itr!=stratSolver_.end()) return itr->second;
00189    
00190    // check solvers
00191    itr = blockPrecond_.find(label);
00192    if(itr!=blockPrecond_.end()) return itr->second;
00193 
00194    return Teuchos::null;
00195 }
00196 
00198 Teuchos::RCP<InverseFactory> InverseLibrary::getInverseFactory(const std::string & label) const
00199 {
00200    Teko_DEBUG_SCOPE("InverseLibrary::getInverseFactory",10);
00201 
00202    std::map<std::string,RCP<const Teuchos::ParameterList> >::const_iterator itr;
00203 
00204    bool isStratSolver=false,isStratPrecond=false,isBlockPrecond=false;
00205 
00206    // is this a Stratimikos solver?
00207    itr = stratPrecond_.find(label);
00208    isStratPrecond = itr!=stratPrecond_.end();
00209 
00210    // is this a Stratimikos preconditioner?
00211    if(not isStratPrecond) {
00212       itr = stratSolver_.find(label);
00213       isStratSolver = itr!=stratSolver_.end();
00214    }
00215 
00216    // must be a "block" preconditioner
00217    if(not (isStratSolver || isStratPrecond)) {
00218       itr = blockPrecond_.find(label);
00219       isBlockPrecond = itr!=blockPrecond_.end();
00220    }
00221 
00222    Teko_DEBUG_MSG("Inverse \"" << label << "\" is of type " 
00223              << "strat prec = " << isStratPrecond << ", "
00224              << "strat solv = " << isStratSolver << ", " 
00225              << "block prec = " << isBlockPrecond,3);
00226 
00227    // Must be one of Strat solver, strat preconditioner, block preconditioner
00228    if(not (isStratSolver || isStratPrecond || isBlockPrecond)) {
00229       RCP<Teuchos::FancyOStream> out = Teuchos::VerboseObjectBase::getDefaultOStream();
00230 
00231       *out << "InverseLibrary::getInverseFactory could not find \"" << label << "\" ... aborting\n";
00232       *out << "Choose one of: " << std::endl; 
00233 
00234       *out << "   Stratimikos preconditioners = ";
00235       for(itr=stratPrecond_.begin();itr!=stratPrecond_.end();++itr)
00236          *out << "      \"" << itr->first << "\"\n";
00237       *out << std::endl;
00238 
00239       *out << "   Stratimikos solvers = ";
00240       for(itr=stratSolver_.begin();itr!=stratSolver_.end();++itr)
00241          *out << "      \"" << itr->first << "\"\n";
00242       *out << std::endl;
00243 
00244       *out << "   Block preconditioners = ";
00245       for(itr=blockPrecond_.begin();itr!=blockPrecond_.end();++itr)
00246          *out << "      \"" << itr->first << "\"\n";
00247       *out << std::endl;
00248 
00249       TEUCHOS_ASSERT(isStratSolver || isStratPrecond || isBlockPrecond);
00250    }
00251    
00252    RCP<const Teuchos::ParameterList> pl = itr->second;
00253 
00254    // build inverse factory
00255    if(isStratPrecond) {
00256       // remove required parameters
00257       RCP<Teuchos::ParameterList> plCopy = rcp(new Teuchos::ParameterList(*pl));
00258       std::string type = plCopy->get<std::string>("Preconditioner Type");
00259       RCP<Teuchos::ParameterList> xtraParams;
00260       if(plCopy->sublist("Preconditioner Types").sublist(type).isParameter("Required Parameters")) {
00261          xtraParams = rcp(new Teuchos::ParameterList(
00262                plCopy->sublist("Preconditioner Types").sublist(type).sublist("Required Parameters"))); 
00263          plCopy->sublist("Preconditioner Types").sublist(type).remove("Required Parameters"); 
00264       }
00265 
00266       // print some debuggin info
00267       Teko_DEBUG_MSG_BEGIN(10);
00268          DEBUG_STREAM << "Printing parameter list: " << std::endl; 
00269          Teko_DEBUG_PUSHTAB(); plCopy->print(DEBUG_STREAM); Teko_DEBUG_POPTAB();
00270 
00271          if(xtraParams!=Teuchos::null) {
00272             DEBUG_STREAM << "Printing extra parameters: " << std::endl; 
00273             Teko_DEBUG_PUSHTAB(); xtraParams->print(DEBUG_STREAM); Teko_DEBUG_POPTAB();
00274          }
00275       Teko_DEBUG_MSG_END();
00276 
00277       Stratimikos::DefaultLinearSolverBuilder strat;
00278       addToStratimikosBuilder(strat);
00279       strat.setParameterList(plCopy);
00280 
00281       // try to build a preconditioner factory
00282       RCP<Thyra::PreconditionerFactoryBase<double> > precFact = strat.createPreconditioningStrategy(type);
00283 
00284       // string must map to a preconditioner
00285       RCP<Teko::PreconditionerInverseFactory> precInvFact 
00286             = rcp(new PreconditionerInverseFactory(precFact,xtraParams,getRequestHandler()));
00287       precInvFact->setupParameterListFromRequestHandler();
00288       return precInvFact;
00289    }
00290    else if(isStratSolver) {
00291       RCP<Teuchos::ParameterList> solveList = rcp(new Teuchos::ParameterList(*pl));
00292       std::string type = solveList->get<std::string>("Linear Solver Type");
00293 
00294       // get preconditioner name, remove "Use Preconditioner" parameter
00295       Teuchos::ParameterList & solveSettings = solveList->sublist("Linear Solver Types").sublist(type);
00296       std::string precKeyWord = "Use Preconditioner";
00297       std::string precName = "None";
00298       if(solveSettings.isParameter(precKeyWord)) {
00299          precName = solveSettings.get<std::string>(precKeyWord);
00300          solveSettings.remove(precKeyWord);
00301       }
00302 
00303       // build Thyra preconditioner factory
00304       RCP<Thyra::PreconditionerFactoryBase<double> > precFactory;
00305       if(precName!="None") {
00306          // we will manually set the preconditioner, so set this to null
00307          solveList->set<std::string>("Preconditioner Type","None");
00308          
00309          // build inverse that preconditioner corresponds to
00310          RCP<PreconditionerInverseFactory> precInvFactory 
00311                = Teuchos::rcp_dynamic_cast<PreconditionerInverseFactory>(getInverseFactory(precName));
00312 
00313          // extract preconditioner factory from preconditioner _inverse_ factory
00314          precFactory = precInvFactory->getPrecFactory();
00315       }
00316 
00317       Stratimikos::DefaultLinearSolverBuilder strat;
00318       addToStratimikosBuilder(strat);
00319       strat.setParameterList(solveList);
00320 
00321       // try to build a solver factory
00322       RCP<Thyra::LinearOpWithSolveFactoryBase<double> > solveFact = strat.createLinearSolveStrategy(type);
00323       if(precFactory!=Teuchos::null)
00324          solveFact->setPreconditionerFactory(precFactory,precName);
00325 
00326       // if its around, build a InverseFactory
00327       return rcp(new SolveInverseFactory(solveFact));
00328    }
00329    else if(isBlockPrecond) {
00330       try {
00331          std::string type = pl->get<std::string>("Preconditioner Type");
00332          const Teuchos::ParameterList & settings = pl->sublist("Preconditioner Settings");
00333    
00334          // build preconditioner factory from the string
00335          RCP<PreconditionerFactory> precFact 
00336                = PreconditionerFactory::buildPreconditionerFactory(type,settings,Teuchos::rcpFromRef(*this));
00337     
00338          TEUCHOS_ASSERT(precFact!=Teuchos::null);
00339    
00340          // return the inverse factory object
00341          return rcp(new PreconditionerInverseFactory(precFact,getRequestHandler()));   
00342       }
00343       catch(std::exception & e) {
00344          RCP<Teuchos::FancyOStream> out = Teko::getOutputStream();
00345          
00346          *out << "Teko: \"getInverseFactory\" failed, Parameter List =\n";
00347          pl->print(*out);
00348 
00349          *out << "*** THROWN EXCEPTION ***\n";
00350          *out << e.what() << std::endl;
00351          *out << "************************\n";
00352          
00353          throw e;
00354       }
00355    }
00356 
00357    TEUCHOS_ASSERT(false);
00358 }
00359 
00361 void InverseLibrary::PrintAvailableInverses(std::ostream & os) const
00362 {
00363    std::map<std::string,Teuchos::RCP<const Teuchos::ParameterList> >::const_iterator itr; 
00364 
00365    os << "Stratimikos Solvers: " << std::endl;
00366    os << "********************************" << std::endl;
00367    for(itr=stratSolver_.begin();itr!=stratSolver_.end();++itr) {
00368       os << "name = \"" << itr->first << "\"" << std::endl;
00369       itr->second->print(os);
00370       os << std::endl;
00371    }
00372 
00373    os << "Stratimikos Preconditioners: " << std::endl;
00374    os << "********************************" << std::endl;
00375    for(itr=stratPrecond_.begin();itr!=stratPrecond_.end();++itr) {
00376       os << "name = \"" << itr->first << "\"" << std::endl;
00377       itr->second->print(os);
00378       os << std::endl;
00379    }
00380 
00381    os << "Teko Preconditioners: " << std::endl;
00382    os << "********************************" << std::endl;
00383    for(itr=blockPrecond_.begin();itr!=blockPrecond_.end();++itr) {
00384       os << "name = \"" << itr->first << "\"" << std::endl;
00385       itr->second->print(os);
00386       os << std::endl;
00387    }
00388 }
00389 
00399 RCP<InverseLibrary> InverseLibrary::buildFromParameterList(const Teuchos::ParameterList & pl,bool useStratDefaults)
00400 {
00401    // build from Stratimikos or allocate a new inverse library
00402    RCP<InverseLibrary> invLib;
00403    if(useStratDefaults)
00404       invLib = InverseLibrary::buildFromStratimikos();
00405    else
00406       invLib = rcp(new InverseLibrary());
00407 
00408    // to convert the void* like entry
00409    Teuchos::ParameterList * temp = 0;
00410 
00411    // loop over all entries in parameter list
00412    Teuchos::ParameterList::ConstIterator itr;
00413    for(itr=pl.begin();itr!=pl.end();++itr) {
00414       // get current entry
00415       std::string label             = itr->first;
00416       Teuchos::ParameterList & list = itr->second.getValue(temp);
00417       
00418       // add to library
00419       invLib->addInverse(label,list);
00420    }
00421    
00422    return invLib;
00423 }
00424 
00434 Teuchos::RCP<InverseLibrary> InverseLibrary::buildFromStratimikos(const Stratimikos::DefaultLinearSolverBuilder & strat)
00435 {
00436    RCP<InverseLibrary> invLib = rcp(new InverseLibrary());
00437 
00438    // get default inveres in Stratimikos
00439    RCP<Teuchos::ParameterList> pl = rcp(new Teuchos::ParameterList(*strat.getValidParameters()));
00440    Teuchos::ParameterList lst(pl->sublist("Linear Solver Types"));
00441    Teuchos::ParameterList pft(pl->sublist("Preconditioner Types"));
00442 
00443    Teuchos::ParameterList::ConstIterator itr;
00444    Teuchos::ParameterList * temp = 0;
00445 
00446    // loop over all entries in solver list
00447    for(itr=lst.begin();itr!=lst.end();++itr) {
00448       // get current entry
00449       std::string label             = itr->first;
00450       Teuchos::ParameterList & list = itr->second.getValue(temp);
00451       list.set("Type",label);
00452       
00453       // add to library
00454       invLib->addInverse(label,list);
00455    }
00456 
00457    // loop over all entries in preconditioner list
00458    for(itr=pft.begin();itr!=pft.end();++itr) {
00459       // get current entry
00460       std::string label             = itr->first;
00461       Teuchos::ParameterList & list = itr->second.getValue(temp);
00462       list.set("Type",label);
00463       
00464       // add to library
00465       invLib->addInverse(label,list);
00466    }
00467 
00468    return invLib;
00469 }
00470 
00471 } // end namespace Teko
 All Classes Files Functions Variables