Belos Version of the Day
BelosAkxFactory.hpp
Go to the documentation of this file.
00001 //@HEADER
00002 // ************************************************************************
00003 //
00004 //                 Belos: Block Linear Solvers Package
00005 //                  Copyright 2004 Sandia Corporation
00006 //
00007 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00008 // the U.S. Government retains certain rights in this software.
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 Michael A. Heroux (maherou@sandia.gov)
00038 //
00039 // ************************************************************************
00040 //@HEADER
00041 
00042 #ifndef __Belos_AkxFactory_hpp
00043 #define __Belos_AkxFactory_hpp
00044 
00045 #include <BelosMonomialOpAkx.hpp>
00046 #include <Teuchos_ParameterList.hpp>
00047 
00048 #include <algorithm>
00049 #include <cctype>
00050 
00051 // TODO Include other OpAkx implementations, and other Akx
00052 // implementations (that don't implement OpAkx) as they are completed.
00053 
00054 namespace Belos {
00055 
00066   template<class Scalar, class MV>
00067   class AkxFactory {
00068   public:
00072     typedef Akx<Scalar, MV> akx_base_type;    
00073 
00077     AkxFactory();
00078 
00100     template<class OP>
00101     Teuchos::RCP<akx_base_type>
00102     makeOpAkx (const Teuchos::RCP<const OP>& A, 
00103          const Teuchos::RCP<const OP>& M_left, 
00104          const Teuchos::RCP<const OP>& M_right,
00105          const Teuchos::RCP<const Teuchos::ParameterList>& params);
00106 
00108     Teuchos::RCP<const Teuchos::ParameterList>
00109     getDefaultParameters ();
00110 
00116     const std::string& parameterListName () const {
00117       return listName_;
00118     }
00119 
00120   private:
00130     std::string
00131     validateBasisName (const std::string& basisName) const;
00132 
00134     const std::string listName_;
00135 
00137     std::vector<std::string> validBasisNames_;
00138     
00142     Teuchos::RCP<const Teuchos::ParameterList> defaultParams_;
00143   };
00144 
00145 
00146   template<class Scalar, class MV>
00147   AkxFactory<Scalar, MV>::AkxFactory() : listName_ ("Akx") 
00148   {
00149     typedef std::vector<std::string>::size_type size_type;
00150     //
00151     // TODO (mfh 14 Feb 2011) Support more basis types later.
00152     //
00153     const size_type numValidBasisNames = 1;
00154     const char* validBasisNames[] = {"Monomial"};
00155     validBasisNames_.resize (numValidBasisNames);
00156 
00157     for (size_type k = 0; k < numValidBasisNames; ++k)
00158       validBasisNames_.push_back (validBasisNames[k]);
00159   }
00160 
00161   template<class Scalar, class MV>
00162   Teuchos::RCP<typename AkxFactory<Scalar, MV>::akx_base_type>
00163   AkxFactory<Scalar, MV>::
00164   makeOpAkx (const Teuchos::RCP<const OP>& A, 
00165        const Teuchos::RCP<const OP>& M_left, 
00166        const Teuchos::RCP<const OP>& M_right,
00167        const Teuchos::RCP<const Teuchos::ParameterList>& params)
00168   {
00169     using Teuchos::ParameterList;
00170     using Teuchos::RCP;
00171     using Teuchos::Exceptions::InvalidParameter;
00172     using Teuchos::Exceptions::InvalidParameterName;
00173     using Teuchos::Exceptions::InvalidParameterType;
00174 
00175     RCP<const ParameterList> plist = 
00176       params.is_null() ? getDefaultParameters() : params;
00177     
00178     std::string basisType;
00179     try {
00180       basisType = plist->get<std::string> ("Basis Type");
00181     } catch (InvalidParameter& e) {
00182       // FIXME (mfh 14 Feb 2011, 03 Mar 2011)
00183       //
00184       // Just rethrow for now.  Later we might wrap and rethrow, or
00185       // implement more tolerant behavior.
00186       throw e; 
00187     }
00188     // Canonicalize (case, etc.) the string and validate it.  Throw an
00189     // exception if invalid, else return the canonicalized name.
00190     const std::string validName = validateBasisName (basisType);
00191 
00192     // Recommended (at least initially) basis length.
00193     int basisLength = 5; // a reasonable default
00194     try {
00195       basisLength = plist->get<int> ("Basis Length");
00196     } catch (InvalidParameter& e) {
00197       // Do nothing; let the default stay
00198     }
00199     TEST_FOR_EXCEPTION(basisLength < 1, std::invalid_argument,
00200            "The \"Basis Length\" parameter must be >= 1, "
00201            "but its value here is " << basisLength << ".");
00202     if (validName == "Monomial")
00203       {
00204   typedef MonomialOpAkx<Scalar, MV, OP> akx_impl_type;
00205   RCP<akx_impl_type> akx (new akx_impl_type (A, M_left, M_right, basisLength));
00206   return akx;
00207       }
00208     else
00209       {
00210   TEST_FOR_EXCEPTION(validName != "Monomial", std::logic_error,
00211          "We have not yet implemented the \"" 
00212          << validName << "\" basis.");
00213   throw std::logic_error("Should never get here!");
00214       }
00215   }
00216 
00217   template<class Scalar, class MV>
00218   Teuchos::RCP<const Teuchos::ParameterList>
00219   AkxFactory<Scalar, MV>::getDefaultParameters ()
00220   {
00221     using Teuchos::RCP;
00222     using Teuchos::ParameterList;
00223     using Teuchos::parameterList;
00224 
00225     if (defaultParams_.is_null())
00226       {
00227   RCP<ParameterList> plist = parameterList (parameterListName ());
00228   // Default basis type is monomial for now; this will likely
00229   // change to the Newton basis.
00230   const std::string defaultBasis ("Monomial");
00231   {
00232     std::ostringstream os;
00233     os << "Default basis type.  Valid choices include: {";
00234     for (std::vector<std::string>::size_type k = 0; 
00235          k < validBasisNames_.size(); ++k)
00236       {
00237         os << validBasisNames_[k];
00238         if (k < validBasisNames_.size() - 1)
00239     os << ", ";
00240       }
00241     os << "}.";
00242     plist->set ("Basis Type", defaultBasis, os.str());
00243   }
00244   // An initial basis length of 5 is a reasonable first guess.
00245   // Solvers should revise this dynamically, based on their
00246   // estimates of the condition number and rank of the generated
00247   // candidate basis.
00248   const int defaultBasisLength = 5;
00249   plist->set ("Basis Length", defaultBasisLength, 
00250         "Default recommended basis length.");
00251   defaultParams_ = plist;
00252       }
00253     return defaultParams_;
00254   }
00255 
00256   template<class Scalar, class MV>
00257   std::string
00258   AkxFactory<Scalar, MV>::validateBasisName (const std::string& basisName) const
00259   {
00260     // Canonicalize the basis name.  First, strip whitespace.
00261     // Second, capitalize the first letter and lowercase the rest.
00262     TEST_FOR_EXCEPTION(basisName.empty(), std::invalid_argument,
00263            "The matrix powers kernel basis name is an empty "
00264            "string.");
00265     const size_t npos = std::string::npos;
00266     size_t firstNonSpacePos = basisName.find_first_not_of (" \t\n");
00267     size_t lastNonSpacePos = basisName.find_last_not_of (" \t\n");
00268     TEST_FOR_EXCEPTION(firstNonSpacePos == npos, std::invalid_argument,
00269            "The matrix powers kernel basis name \"" << basisName 
00270            << "\" contains only whitespace.");
00271     // canonName must have length at least one.
00272     std::string canonName = 
00273       basisName.substr (firstNonSpacePos, lastNonSpacePos-firstNonSpacePos+1);
00274     // Capitalize the first letter and lowercase the rest, in place.
00275     canonName[0] = toupper (canonName[0]);
00276     std::transform (canonName.begin()+1, canonName.end(), 
00277         canonName.begin()+1, tolower);
00278     const bool foundIt = validBasisNames_.end() != 
00279       std::find (validBasisNames_.begin(), validBasisNames_.end(), canonName);
00280     TEST_FOR_EXCEPTION(! foundIt, std::invalid_argument,
00281            "Invalid basis name \"" << basisName << "\".");
00282     return canonName;
00283   }
00284 
00285 
00286 } // namespace Belos
00287 
00288 #endif // __Belos_AkxFactory_hpp
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines