From: Nils Schween Date: Thu, 28 Sep 2023 14:05:26 +0000 (+0200) Subject: Exchange Richardson with GMRES, Fine tune GMRES params and adapt doc X-Git-Tag: relicensing~404^2~4 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=5af1624abb791c97a52319714aca5023420fec36;p=dealii.git Exchange Richardson with GMRES, Fine tune GMRES params and adapt doc --- diff --git a/examples/step-12/step-12.cc b/examples/step-12/step-12.cc index 453284723c..f7662a6cc6 100644 --- a/examples/step-12/step-12.cc +++ b/examples/step-12/step-12.cc @@ -458,8 +458,9 @@ namespace Step12 // @sect3{All the rest} // - // For this simple problem we use the simplest possible solver, called - // Richardson iteration, that represents a simple defect correction. This, in + // For this simple problem we use a standard iterative solver, called GMRES, + // that creates approximate solutions minimising the residual in each + // iterations by adding a new basis vector to the Krylov subspace. This, in // combination with a block SSOR preconditioner, that uses the special block // matrix structure of system matrices arising from DG discretizations. The // size of these blocks are the number of DoFs per cell. Here, we use a SSOR @@ -470,8 +471,17 @@ namespace Step12 template void AdvectionProblem::solve() { - SolverControl solver_control(1000, 1e-12); - SolverRichardson> solver(solver_control); + SolverControl solver_control(1000, 1e-12); + // We create an additional data object for the GMRES solver to increase the + // maximum number of basis vectors of the Krylov subspace. When this number + // is reached the GMRES algorithm is restarted using the solution of the + // previous iteration as the starting approximation. The choice of the + // number of basis vectors is a trade- off between memory consumption and + // convergence speed, since a longer basis means minimization over a larger + // space. + SolverGMRES>::AdditionalData additional_data; + additional_data.max_n_tmp_vectors = 100; + SolverGMRES> solver(solver_control, additional_data); // Here we create the preconditioner, PreconditionBlockSSOR> preconditioner;