Stratimikos Package Browser (Single Doxygen Collection) Version of the Day
test_single_belos_thyra_solver.cpp
Go to the documentation of this file.
00001 #include "test_single_belos_thyra_solver.hpp"
00002 
00003 #ifndef SUN_CXX
00004 
00005 #include "Thyra_BelosLinearOpWithSolveFactory.hpp"
00006 #include "Thyra_LinearOpWithSolveFactoryHelpers.hpp"
00007 #include "Thyra_EpetraLinearOp.hpp"
00008 #include "Thyra_LinearOpTester.hpp"
00009 #include "Thyra_LinearOpWithSolveBase.hpp"
00010 #include "Thyra_LinearOpWithSolveTester.hpp"
00011 #include "Thyra_MultiVectorStdOps.hpp"
00012 #include "Thyra_VectorStdOps.hpp"
00013 #include "EpetraExt_readEpetraLinearSystem.h"
00014 #include "Epetra_SerialComm.h"
00015 #include "Teuchos_ParameterList.hpp"
00016 #ifdef HAVE_BELOS_IFPACK
00017 #  include "Thyra_IfpackPreconditionerFactory.hpp"
00018 #endif
00019 
00020 #endif // SUN_CXX
00021 
00022 bool Thyra::test_single_belos_thyra_solver(
00023   const std::string                       matrixFile
00024   ,const bool                             testTranspose
00025   ,const bool                             usePreconditioner
00026   ,const int                              numRhs
00027   ,const int                              numRandomVectors
00028   ,const double                           maxFwdError
00029   ,const double                           maxResid
00030   ,const double                           maxSolutionError
00031   ,const bool                             showAllTests
00032   ,const bool                             dumpAll
00033   ,Teuchos::ParameterList                 *belosLOWSFPL
00034   ,Teuchos::ParameterList                 *precPL
00035   ,Teuchos::FancyOStream                  *out_arg
00036   )
00037 {
00038   using Teuchos::rcp;
00039   using Teuchos::OSTab;
00040   bool result, success = true;
00041 
00042   Teuchos::RCP<Teuchos::FancyOStream> out = Teuchos::rcp(out_arg,false);
00043 
00044   try {
00045 
00046 #ifndef SUN_CXX
00047 
00048     if(out.get()) {
00049       *out << "\n***"
00050            << "\n*** Testing Thyra::BelosLinearOpWithSolveFactory (and Thyra::BelosLinearOpWithSolve)"
00051            << "\n***\n"
00052            << "\nEchoing input options:"
00053            << "\n  matrixFile             = " << matrixFile
00054            << "\n  testTranspose          = " << testTranspose
00055            << "\n  usePreconditioner      = " << usePreconditioner
00056            << "\n  numRhs                 = " << numRhs
00057            << "\n  numRandomVectors       = " << numRandomVectors
00058            << "\n  maxFwdError            = " << maxFwdError
00059            << "\n  maxResid               = " << maxResid
00060            << "\n  showAllTests           = " << showAllTests
00061            << "\n  dumpAll                = " << dumpAll
00062            << std::endl;
00063     }
00064 
00065     if(out.get()) *out << "\nA) Reading in an epetra matrix A from the file \'"<<matrixFile<<"\' ...\n";
00066   
00067     Epetra_SerialComm comm;
00068     Teuchos::RCP<Epetra_CrsMatrix> epetra_A;
00069     EpetraExt::readEpetraLinearSystem( matrixFile, comm, &epetra_A );
00070 
00071     Teuchos::RCP<const LinearOpBase<double> > A = epetraLinearOp(epetra_A);
00072 
00073     if(out.get() && dumpAll) *out << "\ndescribe(A) =\n" << describe(*A,Teuchos::VERB_EXTREME);
00074 
00075     if(out.get()) *out << "\nB) Creating a BelosLinearOpWithSolveFactory object opFactory ...\n";
00076 
00077     Teuchos::RCP<LinearOpWithSolveFactoryBase<double> >
00078       lowsFactory;
00079     {
00080       Teuchos::RCP<BelosLinearOpWithSolveFactory<double> >
00081         belosLowsFactory = Teuchos::rcp(new BelosLinearOpWithSolveFactory<double>());
00082       lowsFactory = belosLowsFactory;
00083     }
00084 
00085     if(out.get()) {
00086       *out << "\nlowsFactory.getValidParameters() before setting preconditioner factory:\n";
00087       lowsFactory->getValidParameters()->print(OSTab(out).o(),0,true,false);
00088     }
00089 
00090     if(usePreconditioner) {
00091 #ifdef HAVE_BELOS_IFPACK
00092       if(out.get()) {
00093         *out << "\nSetting an Ifpack preconditioner factory ...\n";
00094       }
00095       RCP<PreconditionerFactoryBase<double> >
00096         precFactory = Teuchos::rcp(new IfpackPreconditionerFactory());
00097       if (precPL)
00098         precFactory->setParameterList(rcp(precPL,false));
00099       lowsFactory->setPreconditionerFactory(precFactory,"Ifpack");
00100 #else
00101       TEUCHOS_TEST_FOR_EXCEPT(usePreconditioner);
00102 #endif
00103     }
00104     
00105     if(out.get()) {
00106       *out << "\nlowsFactory.getValidParameters() after setting preconditioner factory:\n";
00107       lowsFactory->getValidParameters()->print(OSTab(out).o(),0,true,false);
00108       *out << "\nbelosLOWSFPL before setting parameters:\n";
00109       belosLOWSFPL->print(OSTab(out).o(),0,true);
00110     }
00111 
00112     lowsFactory->setParameterList(Teuchos::rcp(belosLOWSFPL,false));
00113 
00114     if(out.get()) {
00115       *out << "\nbelosLOWSFPL after setting parameters:\n";
00116       belosLOWSFPL->print(OSTab(out).o(),0,true);
00117     }
00118 
00119     if(out.get()) *out << "\nC) Creating a BelosLinearOpWithSolve object nsA from A ...\n";
00120 
00121     Teuchos::RCP<LinearOpWithSolveBase<double> > nsA = lowsFactory->createOp();
00122     Thyra::initializeOp<double>(*lowsFactory,  A, nsA.ptr());
00123 
00124     if(out.get()) *out << "\nD) Testing the LinearOpBase interface of nsA ...\n";
00125 
00126     LinearOpTester<double> linearOpTester;
00127     linearOpTester.check_adjoint(testTranspose);
00128     linearOpTester.num_rhs(numRhs);
00129     linearOpTester.num_random_vectors(numRandomVectors);
00130     linearOpTester.set_all_error_tol(maxFwdError);
00131     linearOpTester.set_all_warning_tol(1e-2*maxFwdError);
00132     linearOpTester.show_all_tests(showAllTests);
00133     linearOpTester.dump_all(dumpAll);
00134     Thyra::seed_randomize<double>(0);
00135     result = linearOpTester.check(*nsA,Teuchos::Ptr<Teuchos::FancyOStream>(out.get()));
00136     if(!result) success = false;
00137 
00138     if(out.get()) *out << "\nE) Testing the LinearOpWithSolveBase interface of nsA ...\n";
00139     
00140     LinearOpWithSolveTester<double> linearOpWithSolveTester;
00141     linearOpWithSolveTester.num_rhs(numRhs);
00142     linearOpWithSolveTester.turn_off_all_tests();
00143     linearOpWithSolveTester.check_forward_default(true);
00144     linearOpWithSolveTester.check_forward_residual(true);
00145     if(testTranspose) {
00146       linearOpWithSolveTester.check_adjoint_default(true);
00147       linearOpWithSolveTester.check_adjoint_residual(true);
00148     }
00149     else {
00150       linearOpWithSolveTester.check_adjoint_default(false);
00151       linearOpWithSolveTester.check_adjoint_residual(false);
00152     }
00153     linearOpWithSolveTester.set_all_solve_tol(maxResid);
00154     linearOpWithSolveTester.set_all_slack_error_tol(maxResid);
00155     linearOpWithSolveTester.set_all_slack_warning_tol(1e+1*maxResid);
00156     linearOpWithSolveTester.forward_default_residual_error_tol(2*maxResid);
00157     linearOpWithSolveTester.forward_default_solution_error_error_tol(maxSolutionError);
00158     linearOpWithSolveTester.adjoint_default_residual_error_tol(2*maxResid);
00159     linearOpWithSolveTester.adjoint_default_solution_error_error_tol(maxSolutionError);
00160     linearOpWithSolveTester.show_all_tests(showAllTests);
00161     linearOpWithSolveTester.dump_all(dumpAll);
00162     Thyra::seed_randomize<double>(0);
00163     result = linearOpWithSolveTester.check(*nsA,out.get());
00164     if(!result) success = false;
00165 
00166     if(out.get()) {
00167       *out << "\nbelosLOWSFPL after solving:\n";
00168       belosLOWSFPL->print(OSTab(out).o(),0,true);
00169     }
00170     
00171 #else // SUN_CXX
00172     
00173     if(out.get()) *out << "\nTest failed since is was not even compiled since SUN_CXX was defined!\n";
00174     success = false;
00175 
00176 #endif // SUN_CXX
00177 
00178   }
00179   catch( const std::exception &excpt ) {
00180     if(out.get()) *out << std::flush;
00181     std::cerr << "*** Caught standard exception : " << excpt.what() << std::endl;
00182     success = false;
00183   }
00184    
00185   return success;
00186     
00187 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines