47 #include "Tpetra_Core.hpp"
48 #include "Tpetra_CrsMatrix.hpp"
49 #include "MatrixMarket_Tpetra.hpp"
53 #include "Thyra_LinearOpBase.hpp"
54 #include "Thyra_LinearOpWithSolveFactoryExamples.hpp"
55 #include "Thyra_LinearOpWithSolveFactoryHelpers.hpp"
59 int main(
int argc,
char* argv[])
66 using Teuchos::getParameter;
67 typedef double Scalar;
68 typedef Tpetra::Map<>::local_ordinal_type
LO;
69 typedef Tpetra::Map<>::global_ordinal_type
GO;
70 typedef Tpetra::Map<>::node_type NO;
79 out = Teuchos::VerboseObjectBase::getDefaultOStream();
82 Tpetra::ScopeGuard tpetraScope(&argc,&argv);
89 std::string inputFile =
"";
90 std::string extraParams =
"";
91 bool printUnusedParams =
false;
92 RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
93 int MyPID = rank(*comm);
94 verbose = ( out && (MyPID==0) );
96 CommandLineProcessor clp(
false);
97 clp.addOutputSetupOptions(
true);
98 clp.setOption(
"input-file", &inputFile,
"Input file [Required].",
true );
99 clp.setOption(
"extra-params", &extraParams,
"Extra parameters overriding the parameters read in from --input-file");
100 clp.setOption(
"print-unused-params",
"no-print-unused-params", &printUnusedParams,
"Whether to print all unused parameters at the end.");
102 "Testing program for Trilinos Tpetra linear solvers access through Thyra."
105 CommandLineProcessor::EParseCommandLineReturn parse_return = clp.parse(argc,argv);
106 if( parse_return != CommandLineProcessor::PARSE_SUCCESSFUL )
return parse_return;
109 if(verbose) *out <<
"\nReading parameters from XML file \""<<inputFile<<
"\" ...\n";
110 Teuchos::updateParametersFromXmlFile(inputFile, Teuchos::inOutArg(*paramList));
111 if(extraParams.length()) {
112 if(verbose) *out <<
"\nAppending extra parameters from the XML string \""<<extraParams<<
"\" ...\n";
113 Teuchos::updateParametersFromXmlString(extraParams, Teuchos::inOutArg(*paramList));
117 *out <<
"\nEchoing input parameters ...\n";
118 paramList->print(*out,1,
true,
false);
123 validParamList.
set(
"Matrix File",
"fileName");
124 validParamList.
sublist(
"Linear Solver Builder").disableRecursiveValidation();
126 if(verbose) *out <<
"\nValidating top-level input parameters ...\n";
127 paramList->validateParametersAndSetDefaults(validParamList);
130 &matrixFile = getParameter<std::string>(*paramList,
"Matrix File");
132 solverBuilderSL = sublist(paramList,
"Linear Solver Builder",
true);
134 if(verbose) *out <<
"\nReading in an tpetra matrix A from the file \'"<<matrixFile<<
"\' ...\n";
136 RCP<const Tpetra::CrsMatrix<Scalar, LO, GO, NO>> tpetra_A = Tpetra::MatrixMarket::Reader<Tpetra::CrsMatrix<Scalar, LO, GO, NO>>::readSparseFile( matrixFile, comm );
138 RCP<Tpetra::Vector<Scalar,LO,GO,NO>> tpetra_x0 =
rcp(
new Tpetra::Vector<Scalar,LO,GO,NO>(tpetra_A->getDomainMap()));
139 RCP<Tpetra::Vector<Scalar,LO,GO,NO>> tpetra_b =
rcp(
new Tpetra::Vector<Scalar,LO,GO,NO>(tpetra_A->getRangeMap()));
140 tpetra_x0->putScalar(1.0);
141 tpetra_A->apply(*tpetra_x0, *tpetra_b);
142 tpetra_x0->putScalar(0.0);
144 RCP<const Thyra::LinearOpBase<double> >
147 RCP<Thyra::VectorBase<double> >
149 RCP<const Thyra::VectorBase<double> >
152 if(verbose) *out <<
"\nCreating a Stratimikos::DefaultLinearSolverBuilder object ...\n";
158 Stratimikos::enableBelosPrecTpetra<Tpetra::CrsMatrix<Scalar,LO,GO,NO>>(linearSolverBuilder);
162 if(verbose) *out <<
"\nCreating the LinearOpWithSolveFactoryBase object lowsFactory ...\n";
163 RCP<Thyra::LinearOpWithSolveFactoryBase<double> >
164 lowsFactory = createLinearSolveStrategy(linearSolverBuilder);
166 if(verbose) *out <<
"\nChecking the LOWSB interface ...\n";
167 RCP<Thyra::LinearOpWithSolveBase<Scalar> >
168 lowsA = Thyra::linearOpWithSolve<Scalar>(*lowsFactory, A);
171 Thyra::SolveStatus<double> status;
172 status = Thyra::solve<double>(*lowsA, Thyra::NOTRANS, *b, x0.ptr());
174 success = (status.solveStatus == Thyra::SOLVE_STATUS_CONVERGED ?
true :
false);
176 if(verbose && printUnusedParams) {
177 *out <<
"\nPrinting the parameter list (showing what was used) ...\n";
178 paramList->print(*out,1,
true,
true);
183 catch(
const std::exception &excpt ) {
184 std::cerr <<
"*** Caught standard exception : " << excpt.what() << std::endl;
189 if(success) *out <<
"\nCongratulations! The test passed!\n";
190 else *out <<
"\nOh no! At least one of the solves failed!\n";
193 return ( success ? EXIT_SUCCESS : EXIT_FAILURE );
int main(int argc, char *argv[])
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
RCP< ParameterList > sublist(const RCP< ParameterList > ¶mList, const std::string &name, bool mustAlreadyExist=false, const std::string &docString="")
basic_OSTab< char > OSTab
RCP< const LinearOpBase< Scalar > > createConstLinearOp(const RCP< const Tpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraOperator, const RCP< const VectorSpaceBase< Scalar > > rangeSpace=Teuchos::null, const RCP< const VectorSpaceBase< Scalar > > domainSpace=Teuchos::null)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
RCP< VectorBase< Scalar > > createVector(const RCP< Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraVector, const RCP< const VectorSpaceBase< Scalar > > space=Teuchos::null)
void setParameterList(RCP< ParameterList > const ¶mList)
Concrete subclass of Thyra::LinearSolverBuilderBase for creating Thyra::LinearOpWithSolveFactoryBase ...
map_type::local_ordinal_type LO
map_type::global_ordinal_type GO