Thyra Version of the Day
Thyra_SolveSupportTypes.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_SOLVE_SUPPORT_TYPES_HPP
00043 #define THYRA_SOLVE_SUPPORT_TYPES_HPP
00044 
00045 #include "Thyra_OperatorVectorTypes.hpp"
00046 #include "Teuchos_ParameterList.hpp"
00047 #include "Teuchos_FancyOStream.hpp"
00048 #include "Teuchos_Describable.hpp"
00049 
00050 
00051 namespace Thyra {
00052 
00053 
00060 enum ESolveMeasureNormType {
00062   SOLVE_MEASURE_ONE,
00064   SOLVE_MEASURE_NORM_RESIDUAL,
00066   SOLVE_MEASURE_NORM_SOLUTION,
00068   SOLVE_MEASURE_NORM_INIT_RESIDUAL,
00070   SOLVE_MEASURE_NORM_RHS
00071 };
00072 
00073 
00078 inline
00079 const std::string toString(const ESolveMeasureNormType solveMeasureNormType)
00080 {
00081   switch(solveMeasureNormType) {
00082     case SOLVE_MEASURE_ONE:
00083       return "SOLVE_MEASURE_ONE";
00084     case SOLVE_MEASURE_NORM_RESIDUAL:
00085       return "SOLVE_MEASURE_NORM_RESIDUAL";
00086     case SOLVE_MEASURE_NORM_SOLUTION:
00087       return "SOLVE_MEASURE_NORM_SOLUTION";
00088     case SOLVE_MEASURE_NORM_INIT_RESIDUAL:
00089       return "SOLVE_MEASURE_NORM_INIT_RESIDUAL";
00090     case SOLVE_MEASURE_NORM_RHS:
00091       return "SOLVE_MEASURE_NORM_RHS";
00092     default:
00093       TEUCHOS_TEST_FOR_EXCEPT(true);
00094   }
00095   return NULL; // Never be called!
00096 }
00097 
00098 
00112 struct SolveMeasureType {
00114   ESolveMeasureNormType  numerator;
00116   ESolveMeasureNormType  denominator;
00118   SolveMeasureType()
00119     :numerator(SOLVE_MEASURE_ONE),denominator(SOLVE_MEASURE_ONE)
00120     {}
00122   SolveMeasureType(ESolveMeasureNormType _numerator, ESolveMeasureNormType _denominator)
00123     :numerator(_numerator),denominator(_denominator)
00124     {}
00126   void set(ESolveMeasureNormType _numerator, ESolveMeasureNormType _denominator)
00127     { numerator = _numerator; denominator = _denominator; }
00131   bool useDefault() const
00132     { return ( numerator==SOLVE_MEASURE_ONE && denominator==SOLVE_MEASURE_ONE ); }
00134   bool operator()(ESolveMeasureNormType numerator_in,
00135     ESolveMeasureNormType denominator_in
00136     ) const
00137     { return ( numerator==numerator_in && denominator==denominator_in ); }
00139   bool contains(ESolveMeasureNormType measure) const
00140     { return ( numerator==measure || denominator==measure ); }
00141 };
00142 
00143 
00148 inline
00149 std::ostream& operator<<(std::ostream &out, const SolveMeasureType &solveMeasureType)
00150 {
00151   out << "("<<toString(solveMeasureType.numerator)
00152       << "/"<<toString(solveMeasureType.denominator)<<")";
00153   return out;
00154 }
00155 
00156 
00160 template<class Scalar>
00161 class ReductionFunctional : public Teuchos::Describable {
00162 public:
00163 
00166 
00175   typename ScalarTraits<Scalar>::magnitudeType
00176   reduce( const VectorBase<Scalar> &v ) const
00177     {
00178 #ifdef THYRA_DEBUG
00179       TEUCHOS_TEST_FOR_EXCEPTION(!isCompatible(v), Exceptions::IncompatibleVectorSpaces,
00180         "Error, the vector v="<<v.description()<<" is not compatiable with"
00181         " *this="<<this->description()<<"!");
00182 #endif
00183       return reduceImpl(v);
00184     }      
00185 
00189   bool isCompatible( const VectorBase<Scalar> &v ) const
00190     { return isCompatibleImpl(v); }
00191 
00193 
00194 protected:
00195 
00198 
00200   virtual typename ScalarTraits<Scalar>::magnitudeType
00201   reduceImpl( const VectorBase<Scalar> &v ) const = 0;
00202 
00204   virtual bool isCompatibleImpl( const VectorBase<Scalar> &v ) const = 0;
00205 
00207 
00208 };
00209 
00210 
00311 template <class Scalar>
00312 struct SolveCriteria {
00314   typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType ScalarMag;
00316   static ScalarMag unspecifiedTolerance() { return ScalarMag(-1.0); }
00319   SolveMeasureType solveMeasureType;
00322   ScalarMag requestedTol;
00328   RCP<ParameterList> extraParameters;
00331   RCP<const ReductionFunctional<Scalar> > numeratorReductionFunc;
00334   RCP<const ReductionFunctional<Scalar> > denominatorReductionFunc;
00336   SolveCriteria()
00337     : requestedTol(unspecifiedTolerance())
00338     {}
00340   SolveCriteria(
00341     SolveMeasureType solveMeasureType_in,
00342     ScalarMag requestedTol_in,
00343     const RCP<ParameterList> &extraParameters_in = Teuchos::null,
00344     const RCP<ReductionFunctional<Scalar> > &numeratorReductionFunc_in = Teuchos::null,
00345     const RCP<ReductionFunctional<Scalar> > &denominatorReductionFunc_in = Teuchos::null
00346     )
00347     : solveMeasureType(solveMeasureType_in),
00348       requestedTol(requestedTol_in), 
00349       extraParameters(extraParameters_in),
00350       numeratorReductionFunc(numeratorReductionFunc_in),
00351       denominatorReductionFunc(denominatorReductionFunc_in)
00352     {}
00353 };
00354 
00355 
00360 template<class Scalar>
00361 std::ostream& operator<<(std::ostream &out, const SolveCriteria<Scalar> &solveCriteria)
00362 {
00363   out << typeName(solveCriteria) << "{";
00364   out << "solveMeasureType="<<solveCriteria.solveMeasureType;
00365   out << ", requestedTol="<<solveCriteria.requestedTol;
00366   if (nonnull(solveCriteria.extraParameters)) {
00367     out << ", extraParameters="<<solveCriteria.extraParameters;
00368   }
00369   if (nonnull(solveCriteria.numeratorReductionFunc)) {
00370     out << ", numeratorReductionFunc="<<solveCriteria.numeratorReductionFunc->description();
00371   }
00372   if (nonnull(solveCriteria.denominatorReductionFunc)) {
00373     out << ", denominatorReductionFunc="<<solveCriteria.denominatorReductionFunc->description();
00374   }
00375   out << "}";
00376   return out;
00377 }
00378 
00379 
00384 template <class Scalar>
00385 struct THYRA_DEPRECATED BlockSolveCriteria {
00387   SolveCriteria<Scalar> solveCriteria;
00389   int                     numRhs;
00391   BlockSolveCriteria()
00392     : solveCriteria(), numRhs(1)
00393     {}
00395   BlockSolveCriteria( const SolveCriteria<Scalar> &_solveCriteria, int _numRhs )
00396     : solveCriteria(_solveCriteria), numRhs(_numRhs)
00397     {}
00398 };
00399 
00400 
00405 class CatastrophicSolveFailure : public std::runtime_error
00406 {public: CatastrophicSolveFailure(const std::string& what_arg) : std::runtime_error(what_arg) {}};
00407 
00408 
00413 enum ESolveStatus {
00414   SOLVE_STATUS_CONVERGED        
00415   ,SOLVE_STATUS_UNCONVERGED     
00416   ,SOLVE_STATUS_UNKNOWN         
00417 };
00418 
00419 
00424 inline
00425 const std::string toString(const ESolveStatus solveStatus)
00426 {
00427   switch(solveStatus) {
00428     case SOLVE_STATUS_CONVERGED:    return "SOLVE_STATUS_CONVERGED";
00429     case SOLVE_STATUS_UNCONVERGED:  return "SOLVE_STATUS_UNCONVERGED";
00430     case SOLVE_STATUS_UNKNOWN:      return "SOLVE_STATUS_UNKNOWN";
00431     default: TEUCHOS_TEST_FOR_EXCEPT(true);
00432   }
00433   return ""; // Never be called!
00434 }
00435 
00436 
00443 template <class Scalar>
00444 struct SolveStatus {
00446   typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType ScalarMag;
00448   static ScalarMag unknownTolerance() { return ScalarMag(-1); }
00450   ESolveStatus solveStatus;
00454   ScalarMag achievedTol;
00456   std::string message;
00459   RCP<ParameterList> extraParameters;
00461   SolveStatus()
00462     :solveStatus(SOLVE_STATUS_UNKNOWN), achievedTol(unknownTolerance())
00463     {}
00466   static std::string achievedTolToString( const ScalarMag &achievedTol )
00467     {
00468       if(achievedTol==unknownTolerance()) return "unknownTolerance()";
00469       std::ostringstream oss; oss << achievedTol; return oss.str();
00470     }
00471 };
00472 
00473 
00478 template <class Scalar>
00479 std::ostream& operator<<( std::ostream& out_arg, const SolveStatus<Scalar> &solveStatus )
00480 {
00481   RCP<Teuchos::FancyOStream>
00482     out = Teuchos::getFancyOStream(Teuchos::rcp(&out_arg,false));
00483   Teuchos::OSTab tab(out);
00484   *out
00485     << "solveStatus = " << toString(solveStatus.solveStatus) << std::endl
00486     << "achievedTol = " << SolveStatus<Scalar>::achievedTolToString(solveStatus.achievedTol) << std::endl;
00487   *out << "message:";
00488   if (solveStatus.message.length()) {
00489     Teuchos::OSTab tab2(out);
00490     *out << "\n" << solveStatus.message << "\n";
00491   }
00492   *out << "extraParameters:";
00493   if(solveStatus.extraParameters.get()) {
00494     *out << "\n";
00495     Teuchos::OSTab tab3(out);
00496     solveStatus.extraParameters->print(*out, 1000, true);
00497   }
00498   else {
00499     *out << " NONE\n";
00500   }
00501   return out_arg;
00502 }
00503 
00504 
00510 enum ESupportSolveUse {
00511   SUPPORT_SOLVE_UNSPECIFIED  
00512   ,SUPPORT_SOLVE_FORWARD_ONLY  
00513   ,SUPPORT_SOLVE_TRANSPOSE_ONLY  
00514   ,SUPPORT_SOLVE_FORWARD_AND_TRANSPOSE  
00515 };
00516 
00517 
00522 enum EPreconditionerInputType {
00523   PRECONDITIONER_INPUT_TYPE_AS_OPERATOR  
00524   ,PRECONDITIONER_INPUT_TYPE_AS_MATRIX   
00525 };
00526 
00527 
00532 template <class Scalar>
00533 void accumulateSolveStatusInit(
00534   const Ptr<SolveStatus<Scalar> > &overallSolveStatus
00535   )
00536 {
00537   overallSolveStatus->solveStatus = SOLVE_STATUS_CONVERGED;
00538 }
00539 
00540 
00557 template <class Scalar>
00558 void accumulateSolveStatus(
00559   const SolveCriteria<Scalar>, // ToDo: Never used, need to take this out!
00560   const SolveStatus<Scalar> &solveStatus,
00561   const Ptr<SolveStatus<Scalar> > &overallSolveStatus
00562   )
00563 {
00564   switch(solveStatus.solveStatus) {
00565     case SOLVE_STATUS_UNCONVERGED:
00566     {
00567       // First, if we see any unconverged solve status, then the entire block is
00568       // unconverged!
00569       overallSolveStatus->solveStatus = SOLVE_STATUS_UNCONVERGED;
00570       overallSolveStatus->message = solveStatus.message;
00571       overallSolveStatus->extraParameters = solveStatus.extraParameters;
00572       break;
00573     }
00574     case SOLVE_STATUS_UNKNOWN:
00575     {
00576       // Next, if any solve status is unknown, then if the overall solve
00577       // status says converged, then we have to mark it as unknown.  Note that
00578       // unknown could mean that the system is actually converged!
00579       switch(overallSolveStatus->solveStatus) {
00580         case SOLVE_STATUS_CONVERGED:
00581           overallSolveStatus->solveStatus = SOLVE_STATUS_UNKNOWN;
00582           break;
00583         case SOLVE_STATUS_UNCONVERGED:
00584         case SOLVE_STATUS_UNKNOWN:
00585           // If we get here then the overall solve status is either unknown
00586           // already or says unconverged and this will not change here!
00587           overallSolveStatus->message = solveStatus.message;
00588           overallSolveStatus->extraParameters = solveStatus.extraParameters;
00589           break;
00590         default:
00591           TEUCHOS_TEST_FOR_EXCEPT(true); // Corrupted enum?
00592       }
00593       break;
00594     }
00595     case SOLVE_STATUS_CONVERGED:
00596     {
00597       // If we get here then the overall solve status is either unknown,
00598       // unconverged, or converged and this will not change here!
00599       if(overallSolveStatus->message == "")
00600         overallSolveStatus->message = solveStatus.message;
00601       break;
00602     }
00603     default:
00604       TEUCHOS_TEST_FOR_EXCEPT(true); // Corrupted enum?
00605   }
00606   // Update the achieved tolerence to the maximum returned
00607   if( solveStatus.achievedTol > overallSolveStatus->achievedTol ) {
00608     overallSolveStatus->achievedTol = solveStatus.achievedTol;
00609   }
00610   // Set a message if none is set
00611   if(overallSolveStatus->message == "")
00612     overallSolveStatus->message = solveStatus.message;
00613   // Set the extra parameters if none is set
00614   if(overallSolveStatus->extraParameters.get()==NULL)
00615     overallSolveStatus->extraParameters = solveStatus.extraParameters;
00616 }
00617 
00618 
00623 template <class Scalar>
00624 THYRA_DEPRECATED
00625 void accumulateSolveStatus(
00626   const SolveCriteria<Scalar>, // ToDo: Never used, need to take this out!
00627   const SolveStatus<Scalar> &solveStatus,
00628   SolveStatus<Scalar> *overallSolveStatus
00629   )
00630 {
00631   accumulateSolveStatus(
00632     SolveCriteria<Scalar>(),
00633     solveStatus, Teuchos::ptr(overallSolveStatus)
00634     );
00635 }
00636 
00637 
00638 } // namespace Thyra
00639 
00640 
00641 #endif // THYRA_SOLVE_SUPPORT_TYPES_HPP
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines