44 #include "EpetraModelEval2DSim.hpp" 45 #include "EpetraModelEval4DOpt.hpp" 47 #include "Thyra_DefaultModelEvaluatorWithSolveFactory.hpp" 48 #include "Stratimikos_DefaultLinearSolverBuilder.hpp" 49 #include "Thyra_DampenedNewtonNonlinearSolver.hpp" 50 #include "Teuchos_VerboseObject.hpp" 51 #include "Teuchos_CommandLineProcessor.hpp" 52 #include "Teuchos_StandardCatchMacros.hpp" 53 #include "Teuchos_VerbosityLevelCommandLineProcessorHelpers.hpp" 59 const Teuchos::RCP<const Epetra_Vector>
60 createScalingVec(
const double &scale,
const Epetra_Map &map)
62 Teuchos::RCP<Epetra_Vector> scalingVec = Teuchos::rcp(
new Epetra_Vector(map));
63 scalingVec->PutScalar(scale);
68 void scaleEpetraModelEvaluator(
const double &s_x,
const double &s_f,
69 const Teuchos::Ptr<Thyra::EpetraModelEvaluator> &model
73 model->setStateVariableScalingVec(
74 createScalingVec(s_x, *model->getEpetraModel()->get_x_map())
78 model->setStateFunctionScalingVec(
79 createScalingVec(s_f, *model->getEpetraModel()->get_f_map())
88 int main(
int argc,
char* argv[] )
93 using Teuchos::CommandLineProcessor;
94 using Teuchos::outArg;
95 typedef RCP<Thyra::VectorBase<double> > VectorPtr;
105 CommandLineProcessor clp;
106 clp.throwExceptions(
false);
107 clp.addOutputSetupOptions(
true);
110 "This example program solves a simple 2 x 2 set of nonlinear equations using a simple\n" 111 "dampened Newton method.\n\n" 113 "The equations that are solved are:\n\n" 115 " f[0] = x[0] + x[1]*x[1] - p[0];\n" 116 " f[1] = d * ( x[0]*x[0] - x[1] - p[1] );\n\n" 118 "The Jacobian for these equations is nonsingular for every point except x=(-0.5,0.5)\n" 119 "and x=(0.5,-0.5) You can cause the Jacobian to be singular at the solution by setting\n" 120 "p[0]=x[0]+x[1]*x[1] and p[1] = x[0]*x[0]-x[1] for these values of x.\n\n" 122 "The equations are solved using a simple dampended Newton method that uses a Armijo\n" 123 "line search which is implemented in the general class Thyra::DampenedNewtonNonlinearsolver\n" 124 "You can get different levels of detail about the Newton method by adjustingthe command-line\n" 125 "option \"verb-level\" (see above)\n" 129 clp.setOption(
"d", &d,
"Model constant d" );
131 clp.setOption(
"p0", &p0,
"Model constant p[0]" );
133 clp.setOption(
"p1", &p1,
"Model constant p[1]" );
135 clp.setOption(
"x00", &x00,
"Initial guess for x[0]" );
137 clp.setOption(
"x01", &x01,
"Initial guess for x[1]" );
138 Teuchos::EVerbosityLevel verbLevel = Teuchos::VERB_DEFAULT;
139 setVerbosityLevelOption(
"verb-level", &verbLevel,
"Verbosity level.", &clp );
141 clp.setOption(
"tol", &tol,
"Nonlinear solve tolerance" );
143 clp.setOption(
"max-iters", &maxIters,
"Maximum number of nonlinear iterations" );
144 bool use4DOpt =
false;
145 clp.setOption(
"use-4D-opt",
"use-2D-sim", &use4DOpt,
146 "Determines if the EpetraModelEval4DOpt or EpetraModelEval2DSim subclasses are used" );
147 bool externalFactory =
false;
148 clp.setOption(
"external-lowsf",
"internal-lowsf", &externalFactory,
149 "Determines of the Thyra::LinearOpWithSolveFactory is used externally or internally to the Thyra::EpetraModelEvaluator object" );
150 bool showSetInvalidArg =
false;
151 clp.setOption(
"show-set-invalid-arg",
"no-show-set-invalid-arg", &showSetInvalidArg,
152 "Determines if an attempt is made to set an invalid/unsupported ModelEvaluator input argument" );
153 bool showGetInvalidArg =
false;
154 clp.setOption(
"show-get-invalid-arg",
"no-show-get-invalid-arg", &showGetInvalidArg,
155 "Determines if an attempt is made to get an invalid/unsupported ModelEvaluator output argument (2DSim only)" );
157 clp.setOption(
"x-scale", &s_x,
"State variables scaling." );
159 clp.setOption(
"f-scale", &s_f,
"State function scaling." );
161 CommandLineProcessor::EParseCommandLineReturn
162 parse_return = clp.parse(argc,argv,&std::cerr);
164 if( parse_return != CommandLineProcessor::PARSE_SUCCESSFUL )
167 RCP<Teuchos::FancyOStream>
168 out = Teuchos::VerboseObjectBase::getDefaultOStream();
170 *out <<
"\nCreating the nonlinear equations object ...\n";
172 RCP<EpetraExt::ModelEvaluator> epetraModel;
174 epetraModel = rcp(
new EpetraModelEval4DOpt(0.0,0.0,p0,p1,d,x00,x01,p0,p1));
177 epetraModel = rcp(
new EpetraModelEval2DSim(d,p0,p1,x00,x01,showGetInvalidArg));
180 *out <<
"\nCreating linear solver strategy ...\n";
182 Stratimikos::DefaultLinearSolverBuilder linearSolverBuilder;
183 linearSolverBuilder.setParameterList(Teuchos::parameterList());
184 RCP<Thyra::LinearOpWithSolveFactoryBase<double> >
185 lowsFactory = linearSolverBuilder.createLinearSolveStrategy(
"Amesos");
191 RCP<Thyra::EpetraModelEvaluator>
194 RCP<Thyra::ModelEvaluator<double> > thyraModel;
195 if(externalFactory) {
196 epetraThyraModel->initialize(epetraModel, Teuchos::null);
198 new Thyra::DefaultModelEvaluatorWithSolveFactory<double>(
199 epetraThyraModel, lowsFactory
204 epetraThyraModel->initialize(epetraModel, lowsFactory);
205 thyraModel = epetraThyraModel;
208 scaleEpetraModelEvaluator( s_x, s_f, epetraThyraModel.ptr() );
210 if( showSetInvalidArg ) {
211 *out <<
"\nAttempting to set an invalid input argument that throws an exception ...\n\n";
212 Thyra::ModelEvaluatorBase::InArgs<double> inArgs = thyraModel->createInArgs();
213 inArgs.set_x_dot(createMember(thyraModel->get_x_space()));
216 *out <<
"\nCreating the nonlinear solver and solving the equations ...\n\n";
218 Thyra::DampenedNewtonNonlinearSolver<double> newtonSolver;
219 newtonSolver.setVerbLevel(verbLevel);
221 VectorPtr x = createMember(thyraModel->get_x_space());
222 V_V( &*x, *thyraModel->getNominalValues().get_x() );
224 Thyra::SolveCriteria<double> solveCriteria;
225 solveCriteria.solveMeasureType.set(Thyra::SOLVE_MEASURE_NORM_RESIDUAL,Thyra::SOLVE_MEASURE_NORM_RHS);
226 solveCriteria.requestedTol = tol;
227 solveCriteria.extraParameters = Teuchos::parameterList(
"Nonlinear Solve");
228 solveCriteria.extraParameters->set(
"Max Iters",
int(maxIters));
230 newtonSolver.setModel(thyraModel);
231 Thyra::SolveStatus<double>
232 solveStatus = Thyra::solve( newtonSolver, &*x, &solveCriteria );
234 *out <<
"\nNonlinear solver return status:\n";
236 Teuchos::OSTab tab(out);
239 *out <<
"\nFinal solution for x=\n" << *x;
242 TEUCHOS_STANDARD_CATCH_STATEMENTS(
true,std::cerr,success)
244 return success ? 0 : 1;
Concrete Adapter subclass that takes an EpetraExt::ModelEvaluator object and wraps it as a Thyra::Mod...
int main(int argc, char *argv[])