Finite Element Domain Decomposition Library
FEDDLib
Loading...
Searching...
No Matches
NonLinearSolver_def.hpp
1#ifndef NONLINEARSOLVER_DEF_hpp
2#define NONLINEARSOLVER_DEF_hpp
3#include "NonLinearSolver_decl.hpp"
12
13namespace FEDD {
14
15
16template<class SC,class LO,class GO,class NO>
17NonLinearSolver<SC,LO,GO,NO>::NonLinearSolver():
18type_("")
19{}
20
21
22template<class SC,class LO,class GO,class NO>
23NonLinearSolver<SC,LO,GO,NO>::NonLinearSolver(std::string type):
24type_(type)
25{}
26
27template<class SC,class LO,class GO,class NO>
28NonLinearSolver<SC,LO,GO,NO>::~NonLinearSolver(){
29
30}
31
32template<class SC,class LO,class GO,class NO>
33void NonLinearSolver<SC,LO,GO,NO>::solve(NonLinearProblem_Type &problem,vec_dbl_ptr_Type valuesForExport){
34
35 if (!type_.compare("FixedPoint")) {
36 solveFixedPoint(problem,valuesForExport);
37 }
38 else if(!type_.compare("Newton")){
39 solveNewton(problem,valuesForExport);
40 }
41 else if(!type_.compare("NOX")){
42#ifdef FEDD_HAVE_NOX
43 solveNOX(problem,valuesForExport);
44#endif
45 }
46
47}
48
49template<class SC,class LO,class GO,class NO>
50void NonLinearSolver<SC,LO,GO,NO>::solve(TimeProblem_Type &problem, double time, vec_dbl_ptr_Type valuesForExport){
51
52 if (!type_.compare("FixedPoint")) {
53 solveFixedPoint(problem,time);
54 }
55 else if(!type_.compare("Newton")){
56 solveNewton(problem,time, valuesForExport);
57 }
58 else if(!type_.compare("NOX")){
59#ifdef FEDD_HAVE_NOX
60 solveNOX(problem, valuesForExport);
61#endif
62 }
63 else if(!type_.compare("Extrapolation")){
64 solveExtrapolation(problem, time);
65 }
66}
67
68#ifdef FEDD_HAVE_NOX
69template<class SC,class LO,class GO,class NO>
70void NonLinearSolver<SC,LO,GO,NO>::solveNOX(NonLinearProblem_Type &problem,vec_dbl_ptr_Type valuesForExport){
71
72 bool verbose = problem.getVerbose();
73 Teuchos::RCP<NonLinearProblem<SC,LO,GO,NO> > problemPtr = Teuchos::rcpFromRef(problem);
74 Teuchos::RCP<Teuchos::ParameterList> p = sublist(problemPtr->getParameterList(),"ThyraSolver");
75// p->set("Preconditioner Type", "None"); // CH 16.04.19: preconditioner will be built seperately
76// sublist( sublist(p, "Linear Solver Types") , "Belos")->set("Left Preconditioner If Unspecified",true);
77 problemPtr->getLinearSolverBuilder()->setParameterList(p);
78
79 Teuchos::RCP<Thyra::LinearOpWithSolveFactoryBase<SC> > lowsFactory = problemPtr->getLinearSolverBuilder()->createLinearSolveStrategy("");
80
81 //problemPtr->set_W_factory(lowsFactory);
82
83 // Create the initial guess
84 Teuchos::RCP<Thyra::VectorBase<SC> > initial_guess = problemPtr->getNominalValues().get_x()->clone_v();
85 Thyra::V_S(initial_guess.ptr(),Teuchos::ScalarTraits<SC>::zero());
86
87 Teuchos::RCP<Thyra::LinearOpBase<SC> > W_op = problemPtr->create_W_op();
88 Teuchos::RCP<Thyra::PreconditionerBase<SC> > W_prec = problemPtr->create_W_prec();
89
90 // cout << " NonLinearSolver<SC,LO,GO,NO>::solveNOX " << endl;
91
92 // problemPtr->create_W_op();
93
94 Teuchos::RCP<NOX::Thyra::Group> nox_group(new NOX::Thyra::Group(initial_guess,
95 problemPtr.getConst(),
96 W_op,
97 lowsFactory.getConst(),
98 W_prec,
99 Teuchos::null,
100 Teuchos::null,
101 Teuchos::null));
102
103 nox_group->computeF();
104
105 // Create the NOX status tests and the solver
106 // Create the convergence tests
107
108 Teuchos::RCP<NOX::StatusTest::NormUpdate> updateTol =
109 Teuchos::rcp(new NOX::StatusTest::NormUpdate( problemPtr->getParameterList()->sublist("Parameter").get("updateTol",1.0e-6) ) );
110
111 Teuchos::RCP<NOX::StatusTest::RelativeNormF> relresid =
112 Teuchos::rcp(new NOX::StatusTest::RelativeNormF( problemPtr->getParameterList()->sublist("Parameter").get("relNonLinTol",1.0e-4) ) );
113
114 Teuchos::RCP<NOX::StatusTest::NormWRMS> wrms =
115 Teuchos::rcp(new NOX::StatusTest::NormWRMS(problemPtr->getParameterList()->sublist("Parameter").get("relNonLinTol",1.0e-4), problemPtr->getParameterList()->sublist("Parameter").get("absNonLinTol",1.0e-6)));
116
117 Teuchos::RCP<NOX::StatusTest::NormF> absRes =
118 Teuchos::rcp(new NOX::StatusTest::NormF( problemPtr->getParameterList()->sublist("Parameter").get("absNonLinTol",1.0e-6) ) );
119
120 Teuchos::RCP<NOX::StatusTest::Combo> converged;
121
122 if ( !problemPtr->getParameterList()->sublist("Parameter").get("Combo","AND").compare("AND") )
123 converged = Teuchos::rcp(new NOX::StatusTest::Combo(NOX::StatusTest::Combo::AND));
124 else if (!problemPtr->getParameterList()->sublist("Parameter").get("Combo","AND").compare("OR") )
125 converged = Teuchos::rcp(new NOX::StatusTest::Combo(NOX::StatusTest::Combo::OR));
126
127 if ( problemPtr->getParameterList()->sublist("Parameter").get("Use rel tol",true) )
128 converged->addStatusTest(relresid);
129 if ( problemPtr->getParameterList()->sublist("Parameter").get("Use update tol",false) )
130 converged->addStatusTest(updateTol);
131 if (problemPtr->getParameterList()->sublist("Parameter").get("Use WRMS",false))
132 converged->addStatusTest(wrms);
133 if (problemPtr->getParameterList()->sublist("Parameter").get("Use abs tol",false))
134 converged->addStatusTest(absRes);
135
136 Teuchos::RCP<NOX::StatusTest::MaxIters> maxiters =
137 Teuchos::rcp(new NOX::StatusTest::MaxIters(problemPtr->getParameterList()->sublist("Parameter").get("MaxNonLinIts",10)));
138 Teuchos::RCP<NOX::StatusTest::FiniteValue> fv =
139 Teuchos::rcp(new NOX::StatusTest::FiniteValue);
140 Teuchos::RCP<NOX::StatusTest::Combo> combo =
141 Teuchos::rcp(new NOX::StatusTest::Combo(NOX::StatusTest::Combo::OR));
142 combo->addStatusTest(fv);
143 combo->addStatusTest(converged);
144 combo->addStatusTest(maxiters);
145
146 // Create nox parameter list
147 Teuchos::RCP<Teuchos::ParameterList> nl_params = sublist(problemPtr->getParameterList(),"NOXSolver");
148
149 // Create the solver
150 Teuchos::RCP<NOX::Solver::Generic> solver = NOX::Solver::buildSolver(nox_group, combo, nl_params);
151 NOX::StatusTest::StatusType solveStatus = solver->solve();
152 double nonLinearIts = solver->getSolverStatistics()->linearSolve.allNonlinearSolves_NumLinearSolves;
153 double linearIts = solver->getSolverStatistics()->linearSolve.allNonlinearSolves_NumLinearIterations;
154
155 linearIts/=nonLinearIts;
156 if (verbose){
157 std::cout << "############################################################" << std::endl;
158 std::cout << "### Total nonlinear iterations : " << nonLinearIts << " with an average of " << linearIts << " linear iterations ###" << std::endl;
159 std::cout << "############################################################" << std::endl;
160 }
161
162 if ( problemPtr->getParameterList()->sublist("Parameter").get("Cancel MaxNonLinIts",false) ) {
163 TEUCHOS_TEST_FOR_EXCEPTION((int)nonLinearIts == problemPtr->getParameterList()->sublist("Parameter").get("MaxNonLinIts",10) ,std::runtime_error,"Maximum nonlinear Iterations reached. Problem might have converged in the last step. Still we cancel here.");
164 }
165 nonLinearIts_ = nonLinearIts;
166
167}
168
169template<class SC,class LO,class GO,class NO>
170void NonLinearSolver<SC,LO,GO,NO>::solveNOX(TimeProblem_Type &problem, vec_dbl_ptr_Type valuesForExport){
171
172 bool verbose = problem.getVerbose();
173 Teuchos::RCP<TimeProblem_Type> problemPtr = Teuchos::rcpFromRef(problem);
174 Teuchos::RCP<Teuchos::ParameterList> p = sublist(problemPtr->getParameterList(),"ThyraSolver");
175
176 problemPtr->getLinearSolverBuilder()->setParameterList(p);
177
178 Teuchos::RCP<Thyra::LinearOpWithSolveFactoryBase<SC> >
179 lowsFactory = problemPtr->getLinearSolverBuilder()->createLinearSolveStrategy("");
180
181 TEUCHOS_TEST_FOR_EXCEPTION(problemPtr->getSolution()->getNumVectors()>1, std::runtime_error, "With the current implementation NOX can only be used with 1 MultiVector column.");
182 // Create the initial guess and fill with last solution
183 Teuchos::RCP<Thyra::VectorBase<SC> > initialGuess = problemPtr->getNominalValues().get_x()->clone_v();
184 // Try to convert to a ProductVB. If resulting pointer is not null we need to use the ProductMV below, otherwise it is a monolithic vector.
185 Teuchos::RCP<Thyra::ProductVectorBase<SC> > initialGuessProd = Teuchos::rcp_dynamic_cast<Thyra::ProductVectorBase<SC> >(initialGuess);
186 Teuchos::RCP<Thyra::MultiVectorBase<SC> > solMV;
187 if (!initialGuessProd.is_null())
188 solMV = problemPtr->getSolution()->getProdThyraMultiVector();
189 else
190 solMV = problemPtr->getSolution()->getThyraMultiVector();
191
192 Thyra::assign(initialGuess.ptr(), *solMV->col(0));
193 //Thyra::V_S(initialGuess.ptr(),Teuchos::ScalarTraits<SC>::zero());
194 Teuchos::RCP<Thyra::LinearOpBase<SC> > W_op = problemPtr->create_W_op();
195 Teuchos::RCP<Thyra::PreconditionerBase<SC> > W_prec = problemPtr->create_W_prec();
196 // problemPtr->create_W_op();
197
198 Teuchos::RCP<NOX::Thyra::Group> nox_group(new NOX::Thyra::Group(initialGuess,
199 problemPtr.getConst(),
200 W_op,
201 lowsFactory.getConst(),
202 W_prec,
203 Teuchos::null,
204 Teuchos::null,
205 Teuchos::null));
206
207 nox_group->computeF();
208
209 // Create the NOX status tests and the solver
210 // Create the convergence tests
211 Teuchos::RCP<NOX::StatusTest::NormUpdate> updateTol =
212 Teuchos::rcp(new NOX::StatusTest::NormUpdate( problemPtr->getParameterList()->sublist("Parameter").get("updateTol",1.0e-6) ) );
213
214 Teuchos::RCP<NOX::StatusTest::RelativeNormF> relresid =
215 Teuchos::rcp(new NOX::StatusTest::RelativeNormF( problemPtr->getParameterList()->sublist("Parameter").get("relNonLinTol",1.0e-4) ) );
216
217 Teuchos::RCP<NOX::StatusTest::NormWRMS> wrms =
218 Teuchos::rcp(new NOX::StatusTest::NormWRMS(problemPtr->getParameterList()->sublist("Parameter").get("relNonLinTol",1.0e-4), problemPtr->getParameterList()->sublist("Parameter").get("absNonLinTol",1.0e-6)));
219
220 Teuchos::RCP<NOX::StatusTest::NormF> absRes =
221 Teuchos::rcp(new NOX::StatusTest::NormF( problemPtr->getParameterList()->sublist("Parameter").get("absNonLinTol",1.0e-6) ) );
222
223 Teuchos::RCP<NOX::StatusTest::Combo> converged;
224
225 if ( !problemPtr->getParameterList()->sublist("Parameter").get("Combo","AND").compare("AND") )
226 converged = Teuchos::rcp(new NOX::StatusTest::Combo(NOX::StatusTest::Combo::AND));
227 else if (!problemPtr->getParameterList()->sublist("Parameter").get("Combo","AND").compare("OR") )
228 converged = Teuchos::rcp(new NOX::StatusTest::Combo(NOX::StatusTest::Combo::OR));
229
230 if ( problemPtr->getParameterList()->sublist("Parameter").get("Use rel tol",true) )
231 converged->addStatusTest(relresid);
232 if ( problemPtr->getParameterList()->sublist("Parameter").get("Use update tol",false) )
233 converged->addStatusTest(updateTol);
234 if (problemPtr->getParameterList()->sublist("Parameter").get("Use WRMS",false))
235 converged->addStatusTest(wrms);
236 if (problemPtr->getParameterList()->sublist("Parameter").get("Use abs tol",false))
237 converged->addStatusTest(absRes);
238
239 Teuchos::RCP<NOX::StatusTest::MaxIters> maxiters =
240 Teuchos::rcp(new NOX::StatusTest::MaxIters(problemPtr->getParameterList()->sublist("Parameter").get("MaxNonLinIts",10)));
241 Teuchos::RCP<NOX::StatusTest::FiniteValue> fv =
242 Teuchos::rcp(new NOX::StatusTest::FiniteValue);
243 Teuchos::RCP<NOX::StatusTest::Combo> combo =
244 Teuchos::rcp(new NOX::StatusTest::Combo(NOX::StatusTest::Combo::OR));
245 combo->addStatusTest(fv);
246 combo->addStatusTest(converged);
247 combo->addStatusTest(maxiters);
248
249 // Create nox parameter list
250 Teuchos::RCP<Teuchos::ParameterList> nl_params = sublist(problemPtr->getParameterList(),"NOXSolver");
251
252 // Create the solver
253 Teuchos::RCP<NOX::Solver::Generic> solver =
254 NOX::Solver::buildSolver(nox_group, combo, nl_params);
255 NOX::StatusTest::StatusType solveStatus = solver->solve();
256
257 double nonLinearIts = solver->getSolverStatistics()->linearSolve.allNonlinearSolves_NumLinearSolves;
258 double linearIts = solver->getSolverStatistics()->linearSolve.allNonlinearSolves_NumLinearIterations;
259
260 linearIts/=nonLinearIts;
261 if (verbose){
262 std::cout << "############################################################" << std::endl;
263 std::cout << "### Total nonlinear iterations : " << nonLinearIts << " with an average of " << linearIts << " linear iterations ###" << std::endl;
264 std::cout << "############################################################" << std::endl;
265 }
266
267 if ( problemPtr->getParameterList()->sublist("Parameter").get("Cancel MaxNonLinIts",false) ) {
268 TEUCHOS_TEST_FOR_EXCEPTION((int)nonLinearIts == problemPtr->getParameterList()->sublist("Parameter").get("MaxNonLinIts",10) ,std::runtime_error,"Maximum nonlinear Iterations reached. Problem might have converged in the last step. Still we cancel here.");
269 }
270
271 if (!valuesForExport.is_null()) {
272 if (valuesForExport->size() == 2){
273 (*valuesForExport)[0] = linearIts;
274 (*valuesForExport)[1] = nonLinearIts;
275 }
276 }
277}
278#endif
279
280template<class SC,class LO,class GO,class NO>
281void NonLinearSolver<SC,LO,GO,NO>::solveFixedPoint(NonLinearProblem_Type &problem,vec_dbl_ptr_Type valuesForExport){
282
283 bool verbose = problem.getVerbose();
284 TEUCHOS_TEST_FOR_EXCEPTION(problem.getRhs()->getNumVectors()!=1,std::logic_error,"We need to change the code for numVectors>1.");
285 // -------
286 // fix point iteration
287 // -------
288 double gmresIts = 0.;
289 double residual0 = 1.;
290 double residual = 1.;
291
292 double tol = problem.getParameterList()->sublist("Parameter").get("relNonLinTol",1.0e-6);
293 int maxNonLinIts = problem.getParameterList()->sublist("Parameter").get("MaxNonLinIts",10);
294 int nlIts=0;
295
296 double criterionValue = 1.;
297 std::string criterion = problem.getParameterList()->sublist("Parameter").get("Criterion","Residual");
298
299 while ( nlIts < maxNonLinIts ) {
300
301 problem.calculateNonLinResidualVec("reverse");
302
303 problem.setBoundariesSystem();
304
305 if (criterion=="Residual")
306 residual = problem.calculateResidualNorm();
307
308 if (nlIts==0)
309 residual0 = residual;
310
311 if (criterion=="Residual"){
312 criterionValue = residual/residual0;
313 if (verbose)
314 std::cout << "### Fixed Point iteration : " << nlIts << " relative nonlinear residual : " << criterionValue << std::endl;
315 if ( criterionValue < tol )
316 break;
317 }
318
319
320 gmresIts += problem.solveAndUpdate( criterion, criterionValue );
321 nlIts++;
322 if(criterion=="Update"){
323 if (verbose)
324 std::cout << "### Fixed Point iteration : " << nlIts << " residual of update : " << criterionValue << std::endl;
325 if ( criterionValue < tol )
326 break;
327 }
328 // ####### end FPI #######
329 }
330
331 gmresIts/=nlIts;
332 if (verbose)
333 std::cout << "### Total FPI : " << nlIts << " with average gmres its : " << gmresIts << std::endl;
334 if ( problem.getParameterList()->sublist("Parameter").get("Cancel MaxNonLinIts",false) ) {
335 TEUCHOS_TEST_FOR_EXCEPTION( nlIts == maxNonLinIts ,std::runtime_error,"Maximum nonlinear Iterations reached. Problem might have converged in the last step. Still we cancel here.");
336 }
337
338}
339
340template<class SC,class LO,class GO,class NO>
341void NonLinearSolver<SC,LO,GO,NO>::solveNewton( NonLinearProblem_Type &problem, vec_dbl_ptr_Type valuesForExport ){
342
343 bool verbose = problem.getVerbose();
344
345 TEUCHOS_TEST_FOR_EXCEPTION(problem.getRhs()->getNumVectors()!=1,std::logic_error,"We need to change the code for numVectors>1.")
346 // -------
347 // Newton
348 // -------
349 double gmresIts = 0.;
350 double residual0 = 1.;
351 double residual = 1.;
352 double tol = problem.getParameterList()->sublist("Parameter").get("relNonLinTol",1.0e-6);
353 int nlIts=0;
354 int maxNonLinIts = problem.getParameterList()->sublist("Parameter").get("MaxNonLinIts",10);
355 double criterionValue = 1.;
356 std::string criterion = problem.getParameterList()->sublist("Parameter").get("Criterion","Residual");
357
358 while ( nlIts < maxNonLinIts ) {
359 //this makes only sense for Navier-Stokes/Stokes, for other problems, e.g., non linear elasticity, it should do nothing.
360
361 problem.calculateNonLinResidualVec("reverse");
362
363 if (criterion=="Residual")
364 residual = problem.calculateResidualNorm();
365
366 problem.assemble("Newton");
367
368 problem.setBoundariesSystem();
369
370 if (nlIts==0)
371 residual0 = residual;
372
373 if (criterion=="Residual"){
374 criterionValue = residual/residual0;
375 if (verbose)
376 std::cout << "### Newton iteration : " << nlIts << " relative nonlinear residual : " << criterionValue << std::endl;
377 if ( criterionValue < tol )
378 break;
379 }
380 // PRINT INFOS
381 gmresIts += problem.solveAndUpdate( criterion, criterionValue );
382 nlIts++;
383 if(criterion=="Update"){
384 if (verbose)
385 std::cout << "### Newton iteration : " << nlIts << " residual of update : " << criterionValue << std::endl;
386 if ( criterionValue < tol )
387 break;
388 }
389
390 // ####### end FPI #######
391 }
392
393 gmresIts/=nlIts;
394 if (verbose)
395 std::cout << "### Total Newton iterations : " << nlIts << " with average gmres its : " << gmresIts << std::endl;
396 if ( problem.getParameterList()->sublist("Parameter").get("Cancel MaxNonLinIts",false) ) {
397 TEUCHOS_TEST_FOR_EXCEPTION(nlIts == maxNonLinIts ,std::runtime_error,"Maximum nonlinear Iterations reached. Problem might have converged in the last step. Still we cancel here.");
398 }
399
400 if (!valuesForExport.is_null()) {
401 if (valuesForExport->size() == 2){
402 (*valuesForExport)[0] = gmresIts;
403 (*valuesForExport)[1] = nlIts;
404 }
405 }
406}
407
408template<class SC,class LO,class GO,class NO>
409void NonLinearSolver<SC,LO,GO,NO>::solveFixedPoint(TimeProblem_Type &problem, double time){
410
411 bool verbose = problem.getVerbose();
412 problem.setBoundariesRHS(time);
413 TEUCHOS_TEST_FOR_EXCEPTION(problem.getRhs()->getNumVectors()!=1,std::logic_error,"We need to change the code for numVectors>1.")
414
415 // -------
416 // fix point iteration
417 // -------
418 double gmresIts = 0.;
419 double residual0 = 1.;
420 double residual = 1.;
421 double tol = problem.getParameterList()->sublist("Parameter").get("relNonLinTol",1.0e-6);
422 int nlIts=0;
423 int maxNonLinIts = problem.getParameterList()->sublist("Parameter").get("MaxNonLinIts",10);
424 double criterionValue = 1.;
425 std::string criterion = problem.getParameterList()->sublist("Parameter").get("Criterion","Residual");
426
427 while ( nlIts < maxNonLinIts ) {
428
429 problem.calculateNonLinResidualVec("reverse", time);
430
431 if (criterion=="Residual")
432 residual = problem.calculateResidualNorm();
433
434 if (nlIts==0)
435 residual0 = residual;
436
437 // Linearization of system matrix is done in calculateNonLinResidualVec
438 // Now we need to combine it with the mass matrix
439 problem.combineSystems();
440
441 problem.setBoundariesSystem();
442
443 if (criterion=="Residual"){
444 criterionValue = residual/residual0;
445 if (verbose)
446 std::cout << "### Fixed Point iteration : " << nlIts << " relative nonlinear residual : " << criterionValue << std::endl;
447 if ( criterionValue < tol )
448 break;
449 }
450
451 gmresIts += problem.solveAndUpdate( criterion, criterionValue );
452
453 nlIts++;
454 if(criterion=="Update"){
455 if (verbose)
456 std::cout << "### Fixed Point iteration : " << nlIts << " residual of update : " << criterionValue << std::endl;
457 if ( criterionValue < tol )
458 break;
459 }
460 // ####### end FPI #######
461 }
462
463 gmresIts/=nlIts;
464 if (verbose)
465 std::cout << "### Total FPI : " << nlIts << " with average gmres its : " << gmresIts << std::endl;
466 if ( problem.getParameterList()->sublist("Parameter").get("Cancel MaxNonLinIts",false) ) {
467 TEUCHOS_TEST_FOR_EXCEPTION( nlIts == maxNonLinIts ,std::runtime_error,"Maximum nonlinear Iterations reached. Problem might have converged in the last step. Still we cancel here.");
468 }
469}
470
471
472
473template<class SC,class LO,class GO,class NO>
474void NonLinearSolver<SC,LO,GO,NO>::solveNewton(TimeProblem_Type &problem, double time, vec_dbl_ptr_Type valuesForExport ){
475
476 bool verbose = problem.getVerbose();
477 problem.setBoundariesRHS(time);
478
479
480 TEUCHOS_TEST_FOR_EXCEPTION(problem.getRhs()->getNumVectors()!=1,std::logic_error,"We need to change the code for numVectors>1.")
481
482 // -------
483 // Newton iteration
484 // -------
485 double gmresIts = 0.;
486 double residual0 = 1.;
487 double residual = 1.;
488 double tol = problem.getParameterList()->sublist("Parameter").get("relNonLinTol",1.0e-6);
489 int nlIts=0;
490 int maxNonLinIts = problem.getParameterList()->sublist("Parameter").get("MaxNonLinIts",10);
491 double criterionValue = 1.;
492 std::string criterion = problem.getParameterList()->sublist("Parameter").get("Criterion","Residual");
493 std::string timestepping = problem.getParameterList()->sublist("Timestepping Parameter").get("Class","Singlestep");
494
495 while ( nlIts < maxNonLinIts ) {
496 if (timestepping == "External")
497 problem.calculateNonLinResidualVec("external", time);
498 else
499 problem.calculateNonLinResidualVec("reverse", time);
500 if (criterion=="Residual")
501 residual = problem.calculateResidualNorm();
502
503 if (nlIts==0)
504 residual0 = residual;
505
506 if (criterion=="Residual"){
507 criterionValue = residual/residual0;
508// exporterTxt->exportData( criterionValue );
509 if (verbose)
510 std::cout << "### Newton iteration : " << nlIts << " relative nonlinear residual : " << criterionValue << std::endl;
511 if ( criterionValue < tol )
512 break;
513 }
514
515 // Systems are combined in timeProblem.assemble("Newton") and then combined
516 problem.assemble("Newton");
517
518 problem.setBoundariesSystem();
519
520 problem.getSystem()->writeMM("Assembled");
521
522
523 if (timestepping == "External"){//AceGen
524 gmresIts += problem.solveAndUpdate( "ResidualAceGen", criterionValue );
525 // exporterTxt->exportData( criterionValue );
526
527 //problem.assembleExternal( "OnlyUpdate" );// update AceGEN internal variables
528 }
529 else
530 gmresIts += problem.solveAndUpdate( criterion, criterionValue );
531
532 nlIts++;
533
534 //problem.getSolution()->getBlock(0)->print();
535 if(criterion=="Update"){
536 if (verbose)
537 std::cout << "### Newton iteration : " << nlIts << " residual of update : " << criterionValue << std::endl;
538 if ( criterionValue < tol )
539 break;
540 }
541
542 // ####### end FPI #######
543 }
544
545 gmresIts/=nlIts;
546 if (verbose)
547 std::cout << "### Total Newton iteration : " << nlIts << " with average gmres its : " << gmresIts << std::endl;
548 if ( problem.getParameterList()->sublist("Parameter").get("Cancel MaxNonLinIts",false) ) {
549 TEUCHOS_TEST_FOR_EXCEPTION(nlIts == maxNonLinIts ,std::runtime_error,"Maximum nonlinear Iterations reached. Problem might have converged in the last step. Still we cancel here.");
550 }
551 if (!valuesForExport.is_null()) {
552 if (valuesForExport->size() == 2){
553 (*valuesForExport)[0] = gmresIts;
554 (*valuesForExport)[1] = nlIts;
555 }
556
557
558 }
559
560}
561
562
563template<class SC,class LO,class GO,class NO>
564void NonLinearSolver<SC,LO,GO,NO>::solveExtrapolation(TimeProblem<SC,LO,GO,NO> &problem, double time){
565
566 bool verbose = problem.getVerbose();
567
568 problem.assemble("Extrapolation");
569
570 problem.setBoundaries(time); // Setting boundaries to system rhs. The rest of the rhs (e.g. M*u_t) must/should be implemented in DAESolver
571
572 int gmresIts = problem.solve( );
573
574 if (verbose) {
575 std::cout << "### GMRES Its : " << gmresIts << std::endl;
576 }
577}
578}
579#endif
FEDD::NonLinearSolver
Definition NonLinearSolver_decl.hpp:32
FEDD::NonLinearSolver::solve
void solve(NonLinearProblem_Type &problem, vec_dbl_ptr_Type valuesForExport=Teuchos::null)
Call for solving a nonlinear problem, depending on linearization solveNOX/solveFixedPoint/solveNewton...
Definition NonLinearSolver_def.hpp:33
FEDD
Adaptive Mesh Refinement.
Definition AdaptiveMeshRefinement_decl.hpp:33