Thyra  Version of the Day
Thyra_SolveSupportTypes.hpp
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Thyra: Interfaces and Support for Abstract Numerical Algorithms
5 // Copyright (2004) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Roscoe A. Bartlett (bartlettra@ornl.gov)
38 //
39 // ***********************************************************************
40 // @HEADER
41 
42 #ifndef THYRA_SOLVE_SUPPORT_TYPES_HPP
43 #define THYRA_SOLVE_SUPPORT_TYPES_HPP
44 
45 #include "Thyra_OperatorVectorTypes.hpp"
46 #include "Teuchos_ParameterList.hpp"
47 #include "Teuchos_FancyOStream.hpp"
48 #include "Teuchos_Describable.hpp"
49 
50 
51 namespace Thyra {
52 
53 
71 };
72 
73 
78 inline
79 const std::string toString(const ESolveMeasureNormType solveMeasureNormType)
80 {
81  switch(solveMeasureNormType) {
82  case SOLVE_MEASURE_ONE:
83  return "SOLVE_MEASURE_ONE";
85  return "SOLVE_MEASURE_NORM_RESIDUAL";
87  return "SOLVE_MEASURE_NORM_SOLUTION";
89  return "SOLVE_MEASURE_NORM_INIT_RESIDUAL";
91  return "SOLVE_MEASURE_NORM_RHS";
92  default:
93  TEUCHOS_TEST_FOR_EXCEPT(true);
94  }
95  return NULL; // Never be called!
96 }
97 
98 
120  {}
123  :numerator(_numerator),denominator(_denominator)
124  {}
126  void set(ESolveMeasureNormType _numerator, ESolveMeasureNormType _denominator)
127  { numerator = _numerator; denominator = _denominator; }
131  bool useDefault() const
135  ESolveMeasureNormType denominator_in
136  ) const
137  { return ( numerator==numerator_in && denominator==denominator_in ); }
139  bool contains(ESolveMeasureNormType measure) const
140  { return ( numerator==measure || denominator==measure ); }
141 };
142 
143 
148 inline
149 std::ostream& operator<<(std::ostream &out, const SolveMeasureType &solveMeasureType)
150 {
151  out << "("<<toString(solveMeasureType.numerator)
152  << "/"<<toString(solveMeasureType.denominator)<<")";
153  return out;
154 }
155 
156 
160 template<class Scalar>
161 class ReductionFunctional : public Teuchos::Describable {
162 public:
163 
166 
175  typename ScalarTraits<Scalar>::magnitudeType
176  reduce( const VectorBase<Scalar> &v ) const
177  {
178 #ifdef THYRA_DEBUG
179  TEUCHOS_TEST_FOR_EXCEPTION(!isCompatible(v), Exceptions::IncompatibleVectorSpaces,
180  "Error, the vector v="<<v.description()<<" is not compatiable with"
181  " *this="<<this->description()<<"!");
182 #endif
183  return reduceImpl(v);
184  }
185 
189  bool isCompatible( const VectorBase<Scalar> &v ) const
190  { return isCompatibleImpl(v); }
191 
193 
194 protected:
195 
198 
200  virtual typename ScalarTraits<Scalar>::magnitudeType
201  reduceImpl( const VectorBase<Scalar> &v ) const = 0;
202 
204  virtual bool isCompatibleImpl( const VectorBase<Scalar> &v ) const = 0;
205 
207 
208 };
209 
210 
311 template <class Scalar>
314  typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType ScalarMag;
316  static ScalarMag unspecifiedTolerance() { return ScalarMag(-1.0); }
328  RCP<ParameterList> extraParameters;
331  RCP<const ReductionFunctional<Scalar> > numeratorReductionFunc;
334  RCP<const ReductionFunctional<Scalar> > denominatorReductionFunc;
338  {}
341  SolveMeasureType solveMeasureType_in,
342  ScalarMag requestedTol_in,
343  const RCP<ParameterList> &extraParameters_in = Teuchos::null,
344  const RCP<ReductionFunctional<Scalar> > &numeratorReductionFunc_in = Teuchos::null,
345  const RCP<ReductionFunctional<Scalar> > &denominatorReductionFunc_in = Teuchos::null
346  )
347  : solveMeasureType(solveMeasureType_in),
348  requestedTol(requestedTol_in),
349  extraParameters(extraParameters_in),
350  numeratorReductionFunc(numeratorReductionFunc_in),
351  denominatorReductionFunc(denominatorReductionFunc_in)
352  {}
353 };
354 
355 
360 template<class Scalar>
361 std::ostream& operator<<(std::ostream &out, const SolveCriteria<Scalar> &solveCriteria)
362 {
363  out << typeName(solveCriteria) << "{";
364  out << "solveMeasureType="<<solveCriteria.solveMeasureType;
365  out << ", requestedTol="<<solveCriteria.requestedTol;
366  if (nonnull(solveCriteria.extraParameters)) {
367  out << ", extraParameters="<<solveCriteria.extraParameters;
368  }
369  if (nonnull(solveCriteria.numeratorReductionFunc)) {
370  out << ", numeratorReductionFunc="<<solveCriteria.numeratorReductionFunc->description();
371  }
372  if (nonnull(solveCriteria.denominatorReductionFunc)) {
373  out << ", denominatorReductionFunc="<<solveCriteria.denominatorReductionFunc->description();
374  }
375  out << "}";
376  return out;
377 }
378 
379 
384 class CatastrophicSolveFailure : public std::runtime_error
385 {public: CatastrophicSolveFailure(const std::string& what_arg) : std::runtime_error(what_arg) {}};
386 
387 
396 };
397 
398 
403 inline
404 const std::string toString(const ESolveStatus solveStatus)
405 {
406  switch(solveStatus) {
407  case SOLVE_STATUS_CONVERGED: return "SOLVE_STATUS_CONVERGED";
408  case SOLVE_STATUS_UNCONVERGED: return "SOLVE_STATUS_UNCONVERGED";
409  case SOLVE_STATUS_UNKNOWN: return "SOLVE_STATUS_UNKNOWN";
410  default: TEUCHOS_TEST_FOR_EXCEPT(true);
411  }
412  return ""; // Never be called!
413 }
414 
415 
422 template <class Scalar>
423 struct SolveStatus {
425  typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType ScalarMag;
427  static ScalarMag unknownTolerance() { return ScalarMag(-1); }
435  std::string message;
438  RCP<ParameterList> extraParameters;
442  {}
445  static std::string achievedTolToString( const ScalarMag &achievedTol )
446  {
447  if(achievedTol==unknownTolerance()) return "unknownTolerance()";
448  std::ostringstream oss; oss << achievedTol; return oss.str();
449  }
450 };
451 
452 
457 template <class Scalar>
458 std::ostream& operator<<( std::ostream& out_arg, const SolveStatus<Scalar> &solveStatus )
459 {
460  RCP<Teuchos::FancyOStream>
461  out = Teuchos::getFancyOStream(Teuchos::rcp(&out_arg,false));
462  Teuchos::OSTab tab(out);
463  *out
464  << "solveStatus = " << toString(solveStatus.solveStatus) << std::endl
465  << "achievedTol = " << SolveStatus<Scalar>::achievedTolToString(solveStatus.achievedTol) << std::endl;
466  *out << "message:";
467  if (solveStatus.message.length()) {
468  Teuchos::OSTab tab2(out);
469  *out << "\n" << solveStatus.message << "\n";
470  }
471  *out << "extraParameters:";
472  if(solveStatus.extraParameters.get()) {
473  *out << "\n";
474  Teuchos::OSTab tab3(out);
475  solveStatus.extraParameters->print(*out, 10, true);
476  }
477  else {
478  *out << " NONE\n";
479  }
480  return out_arg;
481 }
482 
483 
494 };
495 
496 
504 };
505 
506 
511 template <class Scalar>
513  const Ptr<SolveStatus<Scalar> > &overallSolveStatus
514  )
515 {
516  overallSolveStatus->solveStatus = SOLVE_STATUS_CONVERGED;
517 }
518 
519 
536 template <class Scalar>
538  const SolveCriteria<Scalar>, // ToDo: Never used, need to take this out!
539  const SolveStatus<Scalar> &solveStatus,
540  const Ptr<SolveStatus<Scalar> > &overallSolveStatus
541  )
542 {
543  switch(solveStatus.solveStatus) {
545  {
546  // First, if we see any unconverged solve status, then the entire block is
547  // unconverged!
548  overallSolveStatus->solveStatus = SOLVE_STATUS_UNCONVERGED;
549  overallSolveStatus->message = solveStatus.message;
550  overallSolveStatus->extraParameters = solveStatus.extraParameters;
551  break;
552  }
554  {
555  // Next, if any solve status is unknown, then if the overall solve
556  // status says converged, then we have to mark it as unknown. Note that
557  // unknown could mean that the system is actually converged!
558  switch(overallSolveStatus->solveStatus) {
560  overallSolveStatus->solveStatus = SOLVE_STATUS_UNKNOWN;
561  break;
564  // If we get here then the overall solve status is either unknown
565  // already or says unconverged and this will not change here!
566  overallSolveStatus->message = solveStatus.message;
567  overallSolveStatus->extraParameters = solveStatus.extraParameters;
568  break;
569  default:
570  TEUCHOS_TEST_FOR_EXCEPT(true); // Corrupted enum?
571  }
572  break;
573  }
575  {
576  // If we get here then the overall solve status is either unknown,
577  // unconverged, or converged and this will not change here!
578  if(overallSolveStatus->message == "")
579  overallSolveStatus->message = solveStatus.message;
580  break;
581  }
582  default:
583  TEUCHOS_TEST_FOR_EXCEPT(true); // Corrupted enum?
584  }
585  // Update the achieved tolerence to the maximum returned
586  if( solveStatus.achievedTol > overallSolveStatus->achievedTol ) {
587  overallSolveStatus->achievedTol = solveStatus.achievedTol;
588  }
589  // Set a message if none is set
590  if(overallSolveStatus->message == "")
591  overallSolveStatus->message = solveStatus.message;
592  // Set the extra parameters if none is set
593  if(overallSolveStatus->extraParameters.get()==NULL)
594  overallSolveStatus->extraParameters = solveStatus.extraParameters;
595 }
596 
597 
598 } // namespace Thyra
599 
600 
601 #endif // THYRA_SOLVE_SUPPORT_TYPES_HPP
const std::string toString(const ESolveStatus solveStatus)
bool useDefault() const
Return if this is a default solve measure (default constructed).
How the output LOWSB object will be useded for solves in unspecified.
SolveCriteria()
Default construction to use default solve criteria.
Norm of the current solution vector (i.e. ||x||)
ScalarMag requestedTol
The requested solve tolerance (what the client would like to see). Only significant if !this->solveMe...
virtual ScalarTraits< Scalar >::magnitudeType reduceImpl(const VectorBase< Scalar > &v) const =0
Thrown if vector spaces are incompatible.
The output LOWSB object will only be used for forward solves.
std::string message
A simple one-line message (i.e. no newlines) returned from the solver.
The output LOWSB object will used for forward and transpose solves.
The final solution status is unknown but he solve did not totally fail.
bool operator()(ESolveMeasureNormType numerator_in, ESolveMeasureNormType denominator_in) const
Return if (numerator,denominataor) matches this.
void accumulateSolveStatus(const SolveCriteria< Scalar >, const SolveStatus< Scalar > &solveStatus, const Ptr< SolveStatus< Scalar > > &overallSolveStatus)
Accumulate solve status objects for solving a block of RHSs is smaller sub-blocks.
ScalarTraits< Scalar >::magnitudeType reduce(const VectorBase< Scalar > &v) const
Compute the reduction over a vector.
void accumulateSolveStatusInit(const Ptr< SolveStatus< Scalar > > &overallSolveStatus)
Initial overallSolveStatus before calling accumulateSolveStatus().
SolveCriteria(SolveMeasureType solveMeasureType_in, ScalarMag requestedTol_in, const RCP< ParameterList > &extraParameters_in=Teuchos::null, const RCP< ReductionFunctional< Scalar > > &numeratorReductionFunc_in=Teuchos::null, const RCP< ReductionFunctional< Scalar > > &denominatorReductionFunc_in=Teuchos::null)
Construct with a specified solve criteria.
SolveMeasureType(ESolveMeasureNormType _numerator, ESolveMeasureNormType _denominator)
std::ostream & operator<<(std::ostream &out, const SolveMeasureType &solveMeasureType)
Output operator.
ESolveMeasureNormType
Type of solve measure norm.
The output LOWSB object will only be used for transpose solves.
bool contains(ESolveMeasureNormType measure) const
bool isCompatible(const VectorBase< Scalar > &v) const
Returns true if v is compatible with *this.
ESolveStatus
Solution status.
static ScalarMag unknownTolerance()
Abstract interface for finite-dimensional dense vectors.
Simple struct for the return status from a solve.
Norm of the right-hand side (i.e. ||b||)
ESupportSolveUse
Enum that specifies how a LinearOpWithSolveBase object will be used for solves after it is constructe...
ESolveStatus solveStatus
The return status of the solve.
Teuchos::ScalarTraits< Scalar >::magnitudeType ScalarMag
static ScalarMag unspecifiedTolerance()
No solve measure (i.e. same as 1.0)
ScalarMag achievedTol
The maximum final tolerance actually achieved by the (block) linear solve. A value of unknownToleranc...
static std::string achievedTolToString(const ScalarMag &achievedTol)
Output the achieveTol field.
RCP< ParameterList > extraParameters
Any extra control parameters (e.g. max iterations).
virtual bool isCompatibleImpl(const VectorBase< Scalar > &v) const =0
RCP< const ReductionFunctional< Scalar > > numeratorReductionFunc
Reduction function to be used in place of the natural norm of the numerator.
The input preconditioner should just be applied as an operator.
A general reduction functional to be used in specialized solve convergence criteria.
RCP< const ReductionFunctional< Scalar > > denominatorReductionFunc
Reduction function to be used in place of the natural norm of the numerator.
EPreconditionerInputType
Enum defining the status of a preconditioner object.
RCP< ParameterList > extraParameters
Any extra status parameters. Note that the contents of this parameter list is totally undefined...
SolveMeasureType solveMeasureType
The type of solve tolerance requested as given in this->requestedTol.
The requested solution criteria has likely been achieved.
Exception type thrown on an catastrophic solve failure.
Teuchos::ScalarTraits< Scalar >::magnitudeType ScalarMag
The requested solution criteria has likely not been achieved.
Simple struct that defines the requested solution criteria for a solve.
The input preconditioner should viewed as a matrix to be factored then backsolved as a preconditioner...
ESolveMeasureNormType denominator
Norm of the initial residual vector given a non-zero guess (i.e. ||A*xo-b||)
Norm of the current residual vector (i.e. ||A*x-b||)