Assemble and solve the equation system.
266{
269
271 auto& minus_delta_x =
274
275 bool error_norms_met = false;
276
277
278
280
282
283 int iteration = 1;
285 {
287 double time_dirichlet = 0.0;
288
290 time_iteration.
start();
291
292 timer_dirichlet.
start();
293 sys.computeKnownSolutions(*x[process_id], process_id);
294 sys.applyKnownSolutions(*x[process_id]);
295 time_dirichlet += timer_dirichlet.
elapsed();
296
297 sys.preIteration(iteration, *x[process_id]);
298
300 time_assembly.
start();
301 try
302 {
303 sys.assemble(x, x_prev, process_id);
304 }
305 catch (AssemblyException const& e)
306 {
307 ERR(
"Abort nonlinear iteration. Repeating timestep. Reason: {:s}",
308 e.what());
309 error_norms_met = false;
311 break;
312 }
313 sys.getResidual(*x[process_id], *x_prev[process_id], res);
314 sys.getJacobian(J);
315 INFO(
"[time] Assembly took {:g} s.", time_assembly.
elapsed());
316
317
320
321 minus_delta_x.setZero();
322
323 timer_dirichlet.
start();
324 sys.applyKnownSolutionsNewton(J, res, minus_delta_x);
325 time_dirichlet += timer_dirichlet.
elapsed();
326 INFO(
"[time] Applying Dirichlet BCs took {:g} s.", time_dirichlet);
327
329 {
331 }
332
334 time_linear_solver.
start();
336 INFO(
"[time] Linear solver took {:g} s.", time_linear_solver.
elapsed());
337
338 if (!iteration_succeeded)
339 {
340 ERR(
"Newton: The linear solver failed.");
341 }
342 else
343 {
344
345
346
347
348
349 std::vector<GlobalVector*> x_new{x};
350 x_new[process_id] =
354
355 if (postIterationCallback)
356 {
357 postIterationCallback(iteration, x_new);
358 }
359
360 switch (sys.postIteration(*x_new[process_id]))
361 {
363 break;
365 ERR(
"Newton: The postIteration() hook reported a "
366 "non-recoverable error.");
367 iteration_succeeded = false;
368 break;
371 "Newton: The postIteration() hook decided that this "
372 "iteration has to be repeated.");
373
375 *x_new[process_id]);
376 continue;
377 }
378
380 *x[process_id]);
382 *x_new[process_id]);
383 }
384
385 if (!iteration_succeeded)
386 {
387
388 error_norms_met = false;
389 break;
390 }
391
392 if (sys.isLinear())
393 {
394 error_norms_met = true;
395 }
396 else
397 {
399 {
400
402 *x[process_id]);
403 }
404
406 }
407
408 INFO(
"[time] Iteration #{:d} took {:g} s.", iteration,
410
411 if (error_norms_met)
412 {
413 break;
414 }
415
416
417
419 {
420 break;
421 }
422 }
423
425 {
426 ERR(
"Newton: Could not solve the given nonlinear system within {:d} "
427 "iterations",
429 }
430
434
435 return {error_norms_met, iteration};
436}
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