gassmoeller
commented
2 weeks ago Thanks for the review and the testing! I worked a bit on reproducing the results of the nonlinear channel flow benchmark from your paper (and of course found and fixed a bunch of things that changed in the run and plotting scripts). But the good news is I get quite similar results to the paper. I noticed too late that I am now running with the GMG instead of the AMG solver, but I guess it is good news that there are only some minor differences between them.
Summarizing the differences:
-
With traction boundary conditions we now converge to much smaller nonlinear tolerances in all cases without A block stabilization.
-
With velocity boundary conditions we now converge in way more models. This is really good news!
-
The newton convergence with stabilized A block across all models is slower and less consistent than in the paper (but still faster than Picard). I would guess this is related to the changes of the stabilization we did (or to the switch to GMG).
In summary I am happy with this PR and would suggest we merge this if you agree with my changes to the benchmark scripts.
There is still #6139 to fix. I included a proposed fix in this PR, but will rebase the other PR after this one was merged to make sure the tests capture the parallel case.
For comparison the figures:
Fig.1 from Fraters et al. 2019:
And my reproduction with GMG and ASPECT this branch:
Fig. 2 from Fraters et al. 2019:
And my reproduction with GMG and ASPECT this branch:
MFraters
This is an alternative to #6116 that would be my suggestion for a fix for #6046. The PR consists of the following parts:
compute_pressure_scalingis no longer called as part ofassemble_stokes, instead the calling functions call it (or not). This allows us to always call it for the normal solve, and only call it once at the beginning of each newton iteration. This way it is no longer called when the velocity block is set to 0 and we compute the correct pressure scaling. This alone leads to a significant reduction in linear solver iterations for the Newton solver (at the cost of slightly more nonlinear iterations).solutionalways contains the actual solution, and only for the shortest possible time contains the update of the solution. Immediately after the call tosolve_stokesthe newton solver takes care to restore the solution vector to its previous value and instead stores the update in the new local vectorsearch_directiondcr.residual_old) is updated meant that the Eisenstat-Walker method saw almost no difference betweenresidualandresidual_old. I moved the place whereresidual_oldis updated so that the Eisenstat-Walker method actually sees the difference between the residual of the last iteration and the current one. The placement is a bit tricky, because unlike the normal Stokes solver the Newton solver actually computes the new residual after solving the equation. So storing the old residual at the end of the last iteration is almost identical to computing it again at the beginning of the next iteration. The important part is to update it before the new residual is computed, but after the Eisenstat-Walker computation has been done. I am not sure if this fixes #5146 as well, at least it seems related.This is based on #6133, I will rebase and update tests when that one is merged.