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