Assemble and solve the equation system.
326{
329
331 auto& minus_delta_x =
334
335 bool error_norms_met = false;
336
337
338
340
342
343 int iteration = 1;
345 {
347 double time_dirichlet = 0.0;
348
350 time_iteration.
start();
351
352 timer_dirichlet.
start();
353 sys.computeKnownSolutions(*x[process_id], process_id);
354 sys.applyKnownSolutions(*x[process_id]);
355 time_dirichlet += timer_dirichlet.
elapsed();
356
357 sys.preIteration(iteration, *x[process_id]);
358
360 time_assembly.
start();
361 try
362 {
363 sys.assemble(x, x_prev, process_id);
364 }
365 catch (AssemblyException const& e)
366 {
367 ERR(
"Abort nonlinear iteration. Repeating timestep. Reason: {:s}",
368 e.what());
369 error_norms_met = false;
371 break;
372 }
373 sys.getResidual(*x[process_id], *x_prev[process_id], res);
374 sys.getJacobian(J);
375 INFO(
"[time] Assembly took {:g} s.", time_assembly.
elapsed());
376
377
380
381 minus_delta_x.setZero();
382
383 timer_dirichlet.
start();
384 sys.applyKnownSolutionsNewton(J, res, minus_delta_x);
385 time_dirichlet += timer_dirichlet.
elapsed();
386 INFO(
"[time] Applying Dirichlet BCs took {:g} s.", time_dirichlet);
387
389 {
391 }
392
394 time_linear_solver.
start();
396 INFO(
"[time] Linear solver took {:g} s.", time_linear_solver.
elapsed());
397
398 if (!iteration_succeeded)
399 {
400 ERR(
"Newton: The linear solver failed.");
401 }
402 else
403 {
404
405
406
407
408
409 std::vector<GlobalVector*> x_new{x};
410 x_new[process_id] =
414
415 if (postIterationCallback)
416 {
417 postIterationCallback(iteration, x_new);
418 }
419
420 switch (sys.postIteration(*x_new[process_id]))
421 {
423 break;
425 ERR(
"Newton: The postIteration() hook reported a "
426 "non-recoverable error.");
427 iteration_succeeded = false;
428 break;
431 "Newton: The postIteration() hook decided that this "
432 "iteration has to be repeated.");
433
435 *x_new[process_id]);
436 continue;
437 }
438
440 *x[process_id]);
442 *x_new[process_id]);
443 }
444
445 if (!iteration_succeeded)
446 {
447
448 error_norms_met = false;
449 break;
450 }
451
452 if (sys.isLinear())
453 {
454 error_norms_met = true;
455 }
456 else
457 {
459 {
460
462 *x[process_id]);
463 }
464
466 }
467
468 INFO(
"[time] Iteration #{:d} took {:g} s.", iteration,
470
471 if (error_norms_met)
472 {
473 break;
474 }
475
476
477
479 {
480 break;
481 }
482 }
483
485 {
486 ERR(
"Newton: Could not solve the given nonlinear system within {:d} "
487 "iterations",
489 }
490
494
495 return {error_norms_met, iteration};
496}
void ERR(fmt::format_string< Args... > fmt, Args &&... args)
double elapsed() const
Get the elapsed time in seconds.
void start()
Start the timer.
bool solve(EigenMatrix &A, EigenVector &b, EigenVector &x)
virtual void preFirstIteration()
virtual void checkResidual(GlobalVector const &residual)=0
Check if the residual satisfies the convergence criterion.
virtual bool hasResidualCheck() const =0
virtual bool isSatisfied() const
Tell if the convergence criterion is satisfied.
virtual void checkDeltaX(GlobalVector const &minus_delta_x, GlobalVector const &x)=0
virtual bool hasDeltaXCheck() const =0
virtual void releaseMatrix(GlobalMatrix const &A)=0
virtual GlobalMatrix & getMatrix(std::size_t &id)=0
Get an uninitialized matrix with the given id.
std::size_t _J_id
ID of the Jacobian matrix.
std::size_t _x_new_id
ID of the vector storing .
std::size_t _res_id
ID of the residual vector.
std::size_t _minus_delta_x_id
ID of the vector.
void copy(PETScVector const &x, PETScVector &y)
void axpy(PETScVector &y, PetscScalar const a, PETScVector const &x)
static NUMLIB_EXPORT MatrixProvider & provider