Stratimikos Package Browser (Single Doxygen Collection)  Version of the Day
test_stratimikos_belos_nan_handler.cpp
Go to the documentation of this file.
1 
2 #include "Thyra_BelosLinearOpWithSolveFactory.hpp"
3 #include "Thyra_LinearOpWithSolveFactoryHelpers.hpp"
4 #include "Thyra_MultiVectorStdOps.hpp"
5 #include "Thyra_VectorBase.hpp"
6 #include "Thyra_VectorStdOps.hpp"
7 #include "Thyra_EpetraLinearOp.hpp"
8 #include "EpetraExt_readEpetraLinearSystem.h"
9 #include "Epetra_SerialComm.h"
10 
11 #include "Teuchos_UnitTestHarness.hpp"
12 
13 namespace Thyra {
14 
15  /* This test is to make sure that stratimikos wrappers to belos
16  correctly catch belos exceptions thrown and report back a failed
17  solve. If belos detects a nan in the linear solve, it will throw
18  an exception. This was causing the termination of the
19  simulation. We have codes where cutting the time step allows
20  recovery so it makes sense that stratimikos should catch the
21  thrown exception and report that the solve failed. Then the time
22  integrator can retake a time step. This test makes sure the
23  stratimikos interface catches the exception and reports a failed
24  solve.
25  */
26  TEUCHOS_UNIT_TEST(belos, nan_handling)
27  {
28 
29  Epetra_SerialComm comm;
30  Teuchos::RCP<Epetra_CrsMatrix> epetra_A;
31  std::string matrixFile = "FourByFour.mtx";
32  EpetraExt::readEpetraLinearSystem( matrixFile, comm, &epetra_A );
33 
34  Teuchos::RCP<const LinearOpBase<double> > A = epetraLinearOp(epetra_A);
35 
36  Teuchos::RCP<LinearOpWithSolveFactoryBase<double> >
37  lowsFactory;
38  {
39  Teuchos::RCP<BelosLinearOpWithSolveFactory<double> >
40  belosLowsFactory = Teuchos::rcp(new BelosLinearOpWithSolveFactory<double>());
41  lowsFactory = belosLowsFactory;
42  }
43 
44  Teuchos::ParameterList belosLOWSFPL;
45  {
46  belosLOWSFPL.set("Solver Type","Block GMRES");
47 
48  Teuchos::ParameterList& belosLOWSFPL_solver =
49  belosLOWSFPL.sublist("Solver Types");
50 
51  Teuchos::ParameterList& belosLOWSFPL_gmres =
52  belosLOWSFPL_solver.sublist("Block GMRES");
53 
54  belosLOWSFPL_gmres.set("Maximum Iterations",int(4));
55  belosLOWSFPL_gmres.set("Convergence Tolerance",double(1.0e-4));
56  belosLOWSFPL_gmres.set("Maximum Restarts",int(0));
57  belosLOWSFPL_gmres.set("Block Size",int(1));
58  belosLOWSFPL_gmres.set("Num Blocks",int(4));
59  belosLOWSFPL_gmres.set("Output Frequency",int(1));
60  belosLOWSFPL_gmres.set("Show Maximum Residual Norm Only",bool(false));
61 
62  lowsFactory->setParameterList(Teuchos::rcp(&belosLOWSFPL,false));
63  }
64 
65  Teuchos::RCP<LinearOpWithSolveBase<double> > nsA = lowsFactory->createOp();
66  Thyra::initializeOp<double>(*lowsFactory, A, nsA.ptr());
67 
68  Teuchos::RCP<Thyra::VectorBase<double> > x = Thyra::createMember(A->domain());
69  Teuchos::RCP<Thyra::VectorBase<double> > f = Thyra::createMember(A->range());
70 
71  Thyra::put_scalar(0.0,x.ptr());
72  Thyra::put_scalar(1.0,f.ptr());
73 
74  // Insert a nan
75  Thyra::set_ele(0,Teuchos::ScalarTraits<double>::nan(),x.ptr());
76 
77  Thyra::SolveStatus<double> status = Thyra::solve<double>(*nsA,Thyra::NOTRANS,*f,x.ptr());
78 
79  TEST_EQUALITY(status.solveStatus, Thyra::SOLVE_STATUS_UNCONVERGED);
80  }
81 
82 }
TEUCHOS_UNIT_TEST(belos, nan_handling)
LinearOpWithSolveFactoryBase subclass implemented in terms of Belos.