test_single_stratimikos_solver.cpp

Go to the documentation of this file.
00001 // @HEADER
00002 // ***********************************************************************
00003 // 
00004 //         Stratimikos: Thyra-based strategies for linear solvers
00005 //                Copyright (2006) Sandia Corporation
00006 // 
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 // 
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //  
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //  
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Roscoe A. Bartlett (rabartl@sandia.gov) 
00025 // 
00026 // ***********************************************************************
00027 // @HEADER
00028 
00029 #include "test_single_stratimikos_solver.hpp"
00030 #include "Stratimikos_DefaultLinearSolverBuilder.hpp"
00031 #include "Thyra_EpetraLinearOp.hpp"
00032 #include "Thyra_LinearOpTester.hpp"
00033 #include "Thyra_LinearOpWithSolveTester.hpp"
00034 #include "Thyra_LinearOpWithSolveFactoryExamples.hpp"
00035 #include "Thyra_LinearOpWithSolveFactoryHelpers.hpp"
00036 #include "EpetraExt_readEpetraLinearSystem.h"
00037 #include "Teuchos_ParameterList.hpp"
00038 
00039 #ifdef HAVE_MPI
00040 #  include "Epetra_MpiComm.h"
00041 #else
00042 #  include "Epetra_SerialComm.h"
00043 #endif
00044 
00045 bool Thyra::test_single_stratimikos_solver(
00046   Teuchos::ParameterList *paramList_inout
00047   ,const bool dumpAll
00048   ,Teuchos::FancyOStream *out
00049   )
00050 {
00051 
00052   using Teuchos::rcp;
00053   using Teuchos::RCP;
00054   using Teuchos::OSTab;
00055   using Teuchos::ParameterList;
00056   using Teuchos::getParameter;
00057   typedef double Scalar;
00058 
00059   bool success = true, result = false;
00060 
00061   try {
00062 
00063     TEST_FOR_EXCEPT(!paramList_inout);
00064 
00065     RCP<ParameterList>
00066       paramList = rcp(paramList_inout,false);
00067 
00068     if(out) {
00069       *out << "\nEchoing input parameters ...\n";
00070       paramList->print(*out,1,true,false);
00071     }
00072 
00073     // Create list of valid parameter sublists
00074     Teuchos::ParameterList validParamList("test_single_stratimikos_solver");
00075     validParamList.set("Matrix File","fileName");
00076     validParamList.set("Solve Adjoint",false);
00077     validParamList.sublist("Linear Solver Builder").disableRecursiveValidation();
00078     validParamList.sublist("LinearOpWithSolveTester").disableRecursiveValidation();
00079     
00080     if(out) *out << "\nValidating top-level input parameters ...\n";
00081     paramList->validateParametersAndSetDefaults(validParamList);
00082 
00083     const std::string
00084       &matrixFile = getParameter<std::string>(*paramList,"Matrix File");
00085     const bool
00086       solveAdjoint = getParameter<bool>(*paramList,"Solve Adjoint");
00087     RCP<ParameterList>
00088       solverBuilderSL  = sublist(paramList,"Linear Solver Builder",true),
00089       lowsTesterSL     = sublist(paramList,"LinearOpWithSolveTester",true);
00090 
00091     if(out) *out << "\nReading in an epetra matrix A from the file \'"<<matrixFile<<"\' ...\n";
00092   
00093 #ifdef HAVE_MPI
00094     Epetra_MpiComm comm(MPI_COMM_WORLD);
00095 #else
00096     Epetra_SerialComm comm;
00097 #endif
00098     RCP<Epetra_CrsMatrix> epetra_A;
00099     EpetraExt::readEpetraLinearSystem( matrixFile, comm, &epetra_A );
00100 
00101     RCP<const LinearOpBase<double> >
00102       A = Thyra::epetraLinearOp(epetra_A);
00103 
00104     if(out) *out << "\nCreating a Stratimikos::DefaultLinearSolverBuilder object ...\n";
00105     
00106     RCP<Thyra::LinearSolverBuilderBase<double> >
00107       linearSolverBuilder = rcp(new Stratimikos::DefaultLinearSolverBuilder);
00108 
00109     if(out) {
00110       *out << "\nValid parameters for DefaultLinearSolverBuilder ...\n";
00111       linearSolverBuilder->getValidParameters()->print(*out,1,true,false);
00112     }
00113 
00114     linearSolverBuilder->setParameterList(solverBuilderSL);
00115 
00116     if(out) *out << "\nCreating the LinearOpWithSolveFactoryBase object lowsFactory ...\n";
00117     RCP<LinearOpWithSolveFactoryBase<double> >
00118       lowsFactory = createLinearSolveStrategy(*linearSolverBuilder);
00119     if(out) *out << "\nlowsFactory described as:\n" << describe(*lowsFactory,Teuchos::VERB_MEDIUM) << std::endl;
00120 
00121     if(out) *out << "\nRunning example use cases for not externally preconditioned ...\n";
00122     nonExternallyPreconditionedLinearSolveUseCases(
00123       *A, *lowsFactory, solveAdjoint, *out
00124       );
00125 
00126     Thyra::LinearOpWithSolveTester<Scalar> linearOpWithSolveTester;
00127     linearOpWithSolveTester.setParameterList(lowsTesterSL);
00128     linearOpWithSolveTester.turn_off_all_tests();
00129     linearOpWithSolveTester.check_forward_default(true);
00130     linearOpWithSolveTester.check_forward_residual(true);
00131     if (solveAdjoint) {
00132       linearOpWithSolveTester.check_adjoint_default(true);
00133       linearOpWithSolveTester.check_adjoint_residual(true);
00134     }
00135     // ToDo: Use parameter lists for the above
00136 
00137     if(out) *out << "\nChecking the LOWSB interface ...\n";
00138     RCP<Thyra::LinearOpWithSolveBase<Scalar> >
00139       lowsA = Thyra::linearOpWithSolve<Scalar>(*lowsFactory, A);
00140     result = linearOpWithSolveTester.check(*lowsA, out);
00141     if (!result) success = false;
00142     
00143     if(out) {
00144       *out << "\nPrinting the parameter list (showing what was used) ...\n";
00145       paramList->print(*out,1,true,true);
00146     }
00147     
00148   }
00149   catch( const std::exception &excpt ) {
00150     std::cerr << "*** Caught standard exception : " << excpt.what() << std::endl;
00151     success = false;
00152   }
00153   
00154   return success;
00155   
00156 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends
Generated on Wed Apr 13 10:19:51 2011 for Stratimikos Package Browser (Single Doxygen Collection) by  doxygen 1.6.3