Thyra_SolveSupportTypes.hpp

// @HEADER
// ***********************************************************************
// 
//    Thyra: Interfaces and Support for Abstract Numerical Algorithms
//                 Copyright (2004) Sandia Corporation
// 
// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
// license for use of this work by or on behalf of the U.S. Government.
// 
// This library is free software; you can redistribute it and/or modify
// it under the terms of the GNU Lesser General Public License as
// published by the Free Software Foundation; either version 2.1 of the
// License, or (at your option) any later version.
//  
// This library is distributed in the hope that it will be useful, but
// WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// Lesser General Public License for more details.
//  
// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
// USA
// Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
// 
// ***********************************************************************
// @HEADER

#ifndef THYRA_SOLVE_SUPPORT_TYPES_HPP
#define THYRA_SOLVE_SUPPORT_TYPES_HPP

#include "Thyra_OperatorVectorTypes.hpp"

namespace Thyra {

enum ESolveTolType {
  SOLVE_TOL_DEFAULT                     
  ,SOLVE_TOL_REL_RESIDUAL_NORM          
  ,SOLVE_TOL_REL_SOLUTION_ERR_NORM      
};

inline
const char* toString(const ESolveTolType solveTolType)
{
  switch(solveTolType) {
    case SOLVE_TOL_DEFAULT: return"SOLVE_TOL_DEFAULT";
    case SOLVE_TOL_REL_RESIDUAL_NORM: return "SOLVE_TOL_REL_RESIDUAL_NORM";
    case SOLVE_TOL_REL_SOLUTION_ERR_NORM: return "SOLVE_TOL_REL_SOLUTION_ERR_NORM";
    default: TEST_FOR_EXCEPT(true);
  }
  return ""; // Never be called!
}

template <class Scalar>
struct SolveCriteria {
  typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType ScalarMag;
  static ScalarMag unspecifiedTolerance() { return ScalarMag(-1); }
  static int unspecifiedMaxIterations() { return -1; }
  ESolveTolType    solveTolType;
  ScalarMag        requestedTol;
  int maxIterations;
  SolveCriteria()
    :solveTolType(SOLVE_TOL_DEFAULT)
    ,requestedTol(unspecifiedTolerance())
     ,maxIterations(unspecifiedMaxIterations())
    {}
  SolveCriteria(ESolveTolType _solveTolType, ScalarMag _requestedTol, int _maxIterations = unspecifiedMaxIterations())
    : solveTolType(_solveTolType), requestedTol(_requestedTol), maxIterations(_maxIterations)
    {}
};

template <class Scalar>
struct BlockSolveCriteria {
  SolveCriteria<Scalar>   solveCriteria;
  int                     numRhs;
  BlockSolveCriteria()
    : solveCriteria(), numRhs(1)
    {}
  BlockSolveCriteria( const SolveCriteria<Scalar> &_solveCriteria, int _numRhs )
    : solveCriteria(_solveCriteria), numRhs(_numRhs)
    {}
};

class CatastrophicSolveFailure : public std::runtime_error
{public: CatastrophicSolveFailure(const std::string& what_arg) : std::runtime_error(what_arg) {}};

enum ESolveStatus {
  SOLVE_STATUS_CONVERGED        
  ,SOLVE_STATUS_UNCONVERGED     
  ,SOLVE_STATUS_UNKNOWN         
};

inline
const char* toString(const ESolveStatus solveStatus)
{
  switch(solveStatus) {
    case SOLVE_STATUS_CONVERGED:    return "SOLVE_STATUS_CONVERGED";
    case SOLVE_STATUS_UNCONVERGED:  return "SOLVE_STATUS_UNCONVERGED";
    case SOLVE_STATUS_UNKNOWN:      return "SOLVE_STATUS_UNKNOWN";
    default: TEST_FOR_EXCEPT(true);
  }
  return ""; // Never be called!
}

template <class Scalar>
struct SolveStatus {
  typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType ScalarMag;
  static ScalarMag unknownTolerance() { return ScalarMag(-1); }
  ESolveStatus solveStatus;
  ScalarMag achievedTol;
  int iterations;
  std::string message;
  SolveStatus()
    :solveStatus(SOLVE_STATUS_UNKNOWN), achievedTol(unknownTolerance())
     ,iterations(1)
    {}

  static std::string achievedTolToString( const ScalarMag &achievedTol )
    {
      if(achievedTol==unknownTolerance()) return "unknownTolerance()";
      std::ostringstream oss; oss << achievedTol; return oss.str();
    }
};

enum EPreconditionerInputType {
  PRECONDITIONER_INPUT_TYPE_AS_OPERATOR  
  ,PRECONDITIONER_INPUT_TYPE_AS_MATRIX   
};

template <class Scalar>
void accumulateSolveStatus(
  const SolveCriteria<Scalar>    &overallSolveCriteria
  ,const SolveStatus<Scalar>     &solveStatus
  ,SolveStatus<Scalar>           *overallSolveStatus
  )
{
#ifdef _DEBUG
  TEST_FOR_EXCEPT(overallSolveStatus==NULL);
#endif
  switch(solveStatus.solveStatus) {
    case SOLVE_STATUS_UNCONVERGED:
    {
      // First, if we see any unconverged solve status, then the entire block is
      // unconverged!
      overallSolveStatus->solveStatus = SOLVE_STATUS_UNCONVERGED;
      overallSolveStatus->message = solveStatus.message;
      break;
    }
    case SOLVE_STATUS_UNKNOWN:
    {
      // Next, if any solve status is unknown, then if the overall solve status
      // says converged, then we have to mark it as unknown.  Note than unknown
      // could mean that the system is actually converged!
      switch(overallSolveStatus->solveStatus) {
        case SOLVE_STATUS_CONVERGED:
          overallSolveStatus->solveStatus = SOLVE_STATUS_UNKNOWN;
          break;
        case SOLVE_STATUS_UNCONVERGED:
        case SOLVE_STATUS_UNKNOWN:
          // If we get here then the overall solve status is either unknown
          // already or says unconverged and this will not change here!
          overallSolveStatus->message = solveStatus.message;
          break;
        default:
          TEST_FOR_EXCEPT(true); // Corrupted enum?
      }
      break;
    }
    case SOLVE_STATUS_CONVERGED:
    {
      // If we get here then the overall solve status is either unknown,
      // unconverged, or converged and this will not change here!
      if(overallSolveStatus->message == "")
        overallSolveStatus->message = solveStatus.message;
      break;
    }
    default:
      TEST_FOR_EXCEPT(true); // Corrupted enum?
  }
  // Update the achieved tolerence to the maximum returned
  if( solveStatus.achievedTol > overallSolveStatus->achievedTol ) {
    overallSolveStatus->achievedTol = solveStatus.achievedTol;
  }
  // Accumulate the total number of iterations
  overallSolveStatus->iterations += solveStatus.iterations;
  // Set a message if non is set
  if(overallSolveStatus->message == "")
    overallSolveStatus->message = solveStatus.message;
}

} // namespace Thyra

#endif // THYRA_SOLVE_SUPPORT_TYPES_HPP

---