Belos  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
BelosStatusTestUserOutput.hpp
Go to the documentation of this file.
1 //@HEADER
2 // ************************************************************************
3 //
4 // Belos: Block Linear Solvers Package
5 // Copyright 2004 Sandia Corporation
6 //
7 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8 // the U.S. Government retains certain rights in this software.
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 Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ************************************************************************
40 //@HEADER
41 //
42 
43 #ifndef BELOS_STATUS_TEST_USER_OUTPUT_HPP
44 #define BELOS_STATUS_TEST_USER_OUTPUT_HPP
45 
51 #include <vector>
52 #include "BelosConfigDefs.hpp"
53 #include "BelosTypes.hpp"
54 #include "BelosIteration.hpp"
55 
56 #include "BelosStatusTest.hpp"
57 #include "BelosStatusTestCombo.hpp"
61 
62 namespace Belos {
63 
69 template <class ScalarType, class MV, class OP>
70 class StatusTestUserOutput : public StatusTestOutput<ScalarType,MV,OP> {
71 
76 
77  public:
79 
80 
98  Teuchos::RCP<std::map<std::string,Teuchos::RCP<StatusTest<ScalarType,MV,OP> > > > taggedTests,
99  int mod = 1,
100  int printStates = Passed)
101  : printer_(printer),
102  taggedTests_(taggedTests),
103  state_(Undefined),
104  headerPrinted_(false),
105  stateTest_(printStates),
106  modTest_(mod),
107  lastNumIters_(-1),
108  comboType_(0),
109  blockSize_(1),
110  currNumRHS_(0),
111  currLSNum_(0),
112  numLSDgts_(1),
113  numIterDgts_(1)
114  {
115  // Set the input test.
116  setChild(test);
117  }
118 
120  virtual ~StatusTestUserOutput() {};
122 
124 
125 
143  {
144  TEUCHOS_TEST_FOR_EXCEPTION(iterTest_ == Teuchos::null,StatusTestError,"StatusTestUserOutput::checkStatus(): iteration test pointer is null.");
145  TEUCHOS_TEST_FOR_EXCEPTION(resTestVec_.size() == 0,StatusTestError,"StatusTestUserOutput::checkStatus(): residual test pointer is null.");
146  state_ = test_->checkStatus(solver);
147 
148  // Update some information for the header, if it has not printed or the linear system has changed.
149  LinearProblem<ScalarType,MV,OP> currProb = solver->getProblem();
150  //if (!headerPrinted_ || currLSNum_ != currProb.getLSNumber()) {
151  if (currLSNum_ != currProb.getLSNumber()) {
152  currLSNum_ = currProb.getLSNumber();
153  blockSize_ = solver->getBlockSize();
154  currIdx_ = currProb.getLSIndex();
155  currNumRHS_ = currIdx_.size();
156  numLSDgts_ = (int)std::floor((double)MVT::GetNumberVecs(*(currProb.getRHS())))+1;
157  numIterDgts_ = (int)std::floor(std::log10((double)iterTest_->getMaxIters()))+1;
158  }
159  // Print out current iteration information if it hasn't already been printed, or the status has changed
160  if (((iterTest_->getNumIters() % modTest_ == 0) && (iterTest_->getNumIters()!=lastNumIters_)) || (state_ == Passed)) {
161  lastNumIters_ = iterTest_->getNumIters();
162  if ( (state_ & stateTest_) == state_) {
163  if ( printer_->isVerbosity(StatusTestDetails) ) {
164  print( printer_->stream(StatusTestDetails) );
165  }
166  else if ( printer_->isVerbosity(Debug) ) {
167  print( printer_->stream(Debug) );
168  }
169  }
170  }
171 
172  return state_;
173  }
174 
177  return state_;
178  }
180 
181 
183 
184 
187  void setOutputManager(const Teuchos::RCP<OutputManager<ScalarType> > &printer) { printer_ = printer; }
188 
191  void setOutputFrequency(int mod) { modTest_ = mod; }
192 
198 
199  // First check to see if this test is a combination test
200  Teuchos::RCP<StatusTestCombo_t> comboTest = Teuchos::rcp_dynamic_cast<StatusTestCombo_t>(test);
201  TEUCHOS_TEST_FOR_EXCEPTION(comboTest == Teuchos::null,StatusTestError,"StatusTestUserOutput::setChild: The parameter \"test\" must be a Belos::StatusTestCombo.");
202  std::vector<Teuchos::RCP<StatusTest<ScalarType,MV,OP> > > tmpVec = comboTest->getStatusTests();
203 
204  // Get the number of tests.
205  int numTests = tmpVec.size();
206 
207  // Find the maximum iteration and residual tests
208  for (int i=0; i<numTests; ++i) {
209 
210  // Check if this is a maximum iteration test.
211  Teuchos::RCP<StatusTestMaxIters_t> tmpItrTest = Teuchos::rcp_dynamic_cast<StatusTestMaxIters_t>(tmpVec[i]);
212  if (tmpItrTest != Teuchos::null) {
213  iterTest_ = tmpItrTest;
214  continue;
215  }
216 
217  // Check if this is a single residual test
218  // this should only be the case if there are no user-specified tests when using the hard-coded solver-
219  // specific tests only
220  Teuchos::RCP<StatusTestResNorm_t> tmpResTest = Teuchos::rcp_dynamic_cast<StatusTestResNorm_t>(tmpVec[i]);
221  // If the residual status test is a single test, put in the vector
222  if (tmpResTest != Teuchos::null) {
223  resTestVec_.resize( 1 );
224  resTestVec_[0] = tmpResTest;
225  resTestNamesVec_.resize( 1 );
226  resTestNamesVec_[0] = "IMPLICIT RES";
227  continue;
228  }
229 
230  // Check if the residual test is a combination of several StatusTestResNorm objects.
231  // This should be the standard: we have a combo of the solver-specific tests and user-specific tests
232  // The user-specific tests are probably a combo again.
233  Teuchos::RCP<StatusTestCombo_t> tmpComboTest = Teuchos::rcp_dynamic_cast<StatusTestCombo_t>(tmpVec[i]);
234  TEUCHOS_TEST_FOR_EXCEPTION(tmpComboTest == Teuchos::null,StatusTestError,"StatusTestUserOutput(): test must be Belos::StatusTest[MaxIters|ResNorm|Combo].");
235  tmpVec = tmpComboTest->getStatusTests();
236  comboType_ = tmpComboTest->getComboType();
237 
238  // Add only status tests which are not in the user-specified list of tagged status tests
239  // More specifically: we want to add the implicit residual test
240  typename std::map<std::string,Teuchos::RCP<StatusTest<ScalarType,MV,OP> > >::iterator it;
241  for (size_t j=0; j<tmpVec.size(); ++j) {
242  tmpResTest = Teuchos::rcp_dynamic_cast<StatusTestResNorm_t>(tmpVec[j]);
243 
244  if(tmpResTest == Teuchos::null) continue;
245 
246  bool bFound = false;
247  for(it = taggedTests_->begin(); it != taggedTests_->end(); ++it) {
248  if(tmpVec[j] == it->second) { bFound = true; break; }
249  }
250  if(!bFound) {
251  resTestVec_.push_back(tmpResTest);
252  resTestNamesVec_.push_back("IMPLICIT RES");
253  }
254  }
255 
256  // add all tagged tests (the ordering is by the Tag names in alphabetical ordering)
257  for(it = taggedTests_->begin(); it != taggedTests_->end(); ++it) {
258  resTestVec_.push_back(it->second);
259  resTestNamesVec_.push_back(it->first);
260  }
261  }
262 
263  // Keep the pointer to the new test and reset the state to Undefined.
264  test_ = test;
265  state_ = Undefined;
266 
267  }
268 
271  return test_;
272  }
273 
276  void setSolverDesc(const std::string& solverDesc) { solverDesc_ = solverDesc; }
277 
280  void setPrecondDesc(const std::string& precondDesc) { precondDesc_ = precondDesc; }
282 
283 
285 
286 
291  void reset() {
292  state_ = Undefined;
293  test_->reset();
294  lastNumIters_ = -1;
295  headerPrinted_ = false;
296  }
297 
299 
302  void resetNumCalls() {}
303 
305 
307 
308 
310  void print(std::ostream& os, int indent = 0) const {
311  std::string ind(indent,' ');
312  std::string starLine(55,'*');
313  std::string starFront(5,'*');
314 
315  std::ios_base::fmtflags osFlags(os.flags());
316 
317  os.setf(std::ios::scientific, std::ios::floatfield);
318  os.precision(6);
319 
320  // Print header if this is the first call to this output status test.
321  if (!headerPrinted_) {
322  os << std::endl << ind << starLine << std::endl;
323  os << ind << starFront << " Belos Iterative Solver: " << solverDesc_ << std::endl;
324  if (precondDesc_ != "")
325  os << ind << starFront << " Preconditioner: " << precondDesc_ << std::endl;
326  os << ind << starFront << " Maximum Iterations: " << iterTest_->getMaxIters() << std::endl;
327  os << ind << starFront << " Block Size: " << blockSize_ << std::endl;
328  os << ind << starFront << " Status tests: " << std::endl;
329  test_->print(os,indent + 3);
330  os << ind << starLine << std::endl;
331  os.setf(std::ios_base::right, std::ios_base::adjustfield);
332  std::string indheader( 7 + numIterDgts_, ' ' );
333  os << indheader;
334  for (int i=0; i<currNumRHS_; ++i) {
335  if ( i > 0 && currIdx_[i]!=-1 ) {
336  // Put in space where 'Iter :' is in the previous lines
337  os << ind << indheader;
338  }
339  os << "[" << std::setw(numLSDgts_) << currIdx_[i]+1 << "] : ";
340  for (size_t j=0; j<resTestVec_.size(); ++j) {
341  os << std::setw(15) << resTestNamesVec_[j];
342  }
343  os << std::endl;
344  }
345  headerPrinted_ = true;
346  }
347 
348  // Print out residuals for each residual test.
349  os.setf(std::ios_base::right, std::ios_base::adjustfield);
350  std::string ind2( 7 + numIterDgts_, ' ' );
351  os << ind << "Iter " << std::setw(numIterDgts_) << iterTest_->getNumIters() << ", ";
352  for (int i=0; i<currNumRHS_; ++i) {
353  if ( i > 0 && currIdx_[i]!=-1 ) {
354  // Put in space where 'Iter :' is in the previous lines
355  os << ind << ind2;
356  }
357  os << "[" << std::setw(numLSDgts_) << currIdx_[i]+1 << "] : ";
358  for (size_t j=0; j<resTestVec_.size(); ++j) {
359  if ( resTestVec_[j]->getStatus() != Undefined && currIdx_[i]!=-1 ) {
360  // distinguish between ResNormTest derived and others
361  Teuchos::RCP<StatusTestResNorm_t> tempResTest = Teuchos::rcp_dynamic_cast<StatusTestResNorm_t>(resTestVec_[j]);
362  if(tempResTest != Teuchos::null)
363  os << std::setw(15) << (*tempResTest->getTestValue())[currIdx_[i]];
364  else {
365  if(resTestVec_[j]->getStatus() == Belos::Passed)
366  os << std::setw(15) << "Passed";
367  else if(resTestVec_[j]->getStatus() == Belos::Failed)
368  os << std::setw(15) << "Failed";
369  else os << std::setw(15) << "Undefined";
370  }
371  } else {
372  os << std::setw(15) << "---";
373  }
374  }
375  os << std::endl;
376  }
377  // reset os format
378  os.flags(osFlags);
379  }
380 
382 
383  private:
384  // Output manager.
386 
387  // Overall status test.
389 
390  // tagged tests
392 
393  // Iteration test (as passed in).
395 
397  std::vector<Teuchos::RCP<StatusTest<ScalarType,MV,OP> > > resTestVec_;
398 
400  std::vector<std::string> resTestNamesVec_;
401 
402  std::string solverDesc_;
403  std::string precondDesc_;
404  std::vector<int> currIdx_;
405  StatusType state_;
406  mutable bool headerPrinted_;
407  int stateTest_, modTest_;
408  int lastNumIters_, comboType_;
409  int blockSize_;
410  int currNumRHS_, currLSNum_;
411  int numLSDgts_, numIterDgts_;
412 };
413 
414 } // end of Belos namespace
415 
416 #endif /* BELOS_STATUS_TEST_USER_OUTPUT_HPP */
Collection of types and exceptions used within the Belos solvers.
Belos&#39;s basic output manager for sending information of select verbosity levels to the appropriate ou...
void reset()
Informs the status test that it should reset its internal configuration to the uninitialized state...
virtual int getBlockSize() const =0
Get the blocksize to be used by the iterative solver in solving this linear problem.
Exception thrown to signal error in a status test during Belos::StatusTest::checkStatus().
void print(std::ostream &os, int indent=0) const
Output formatted description of stopping test to output stream.
A special StatusTest for printing other status tests in a simple format.
Virtual base class for StatusTest that printing status tests.
Pure virtual base class for defining the status testing capabilities of Belos.
ComboType getComboType() const
Return the type of combination (OR, AND, or SEQ).
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
An abstract class of StatusTest for stopping criteria using residual norms.
void setSolverDesc(const std::string &solverDesc)
Set a short solver description for output clarity.
Belos::StatusTest class for specifying a maximum number of iterations.
Pure virtual base class which describes the basic interface to the linear solver iteration.
Teuchos::RCP< const MV > getRHS() const
A pointer to the right-hand side B.
A pure virtual class for defining the status tests for the Belos iterative solvers.
StatusType
Whether the StatusTest wants iteration to stop.
Definition: BelosTypes.hpp:189
static int GetNumberVecs(const MV &mv)
Obtain the number of vectors in mv.
Traits class which defines basic operations on multivectors.
Belos::StatusTest for logically combining several status tests.
A Belos::StatusTest class for specifying a maximum number of iterations.
const std::vector< int > & getLSIndex() const
(Zero-based) indices of the linear system(s) currently being solved.
void resetNumCalls()
Informs the outputting status test that it should reset the number of calls to zero.
A linear system to solve, and its associated information.
int getLSNumber() const
The number of linear systems that have been set.
Teuchos::RCP< StatusTest< ScalarType, MV, OP > > getChild() const
Get child test.
StatusTestUserOutput(const Teuchos::RCP< OutputManager< ScalarType > > &printer, Teuchos::RCP< StatusTest< ScalarType, MV, OP > > test, Teuchos::RCP< std::map< std::string, Teuchos::RCP< StatusTest< ScalarType, MV, OP > > > > taggedTests, int mod=1, int printStates=Passed)
Constructor.
virtual const LinearProblem< ScalarType, MV, OP > & getProblem() const =0
Get a constant reference to the linear problem.
virtual const std::vector< MagnitudeType > * getTestValue() const =0
Returns the test value, , computed in most recent call to CheckStatus.
Belos::StatusTest abstract class for specifying a residual norm stopping criteria.
void setPrecondDesc(const std::string &precondDesc)
Set a short preconditioner description for output clarity.
void setChild(Teuchos::RCP< StatusTest< ScalarType, MV, OP > > test)
Set child test, which must be a combination of a Belos::StatusTestMaxIters AND a combination of Belos...
void setOutputManager(const Teuchos::RCP< OutputManager< ScalarType > > &printer)
Set the output manager.
StatusType getStatus() const
Return the result of the most recent checkStatus call, or undefined if it has not been run...
StatusType checkStatus(Iteration< ScalarType, MV, OP > *solver)
A class for extending the status testing capabilities of Belos via logical combinations.
st_vector getStatusTests()
Return the vector of status tests.
Belos header file which uses auto-configuration information to include necessary C++ headers...
A virtual base class for StatusTest that print other status tests.
void setOutputFrequency(int mod)
Set how often the child test is printed.

Generated on Thu Mar 28 2024 09:27:41 for Belos by doxygen 1.8.5