Damping for Picard fixed point iterations
Part of the !5714 (merged)
Test case
The test case is a heat conduction problem:
Governing PDE (transient heat conduction, 1D, x \in [0,1]) reduces to:
\frac{\partial T}{\partial t} - \frac{\partial^2 T}{\partial x^2} = 0
Initial condition T(x,0) = 0. Boundary conditions:
- Left, x=0 — Dirichlet: T(0,t) = 0
- Right, x=1 — nonlinear Robin flux: q(T) = \alpha(T),(T - T_\infty), \qquad \alpha(T) = \frac{h_0}{1 + T^2} with T_\infty = 20, h_0 = 100.
Damping
Picard treats the nonlinear coefficient with a lag: at iterate x_k it freezes \alpha at the old temperature, solves one linear system, gets raw output g(x_k). Fixpoint iteration:
x_{k+1}=g(x_k)
converges only if the map g is a contraction near the solution, i.e. spectral radius \rho(g'(x^*))<1.
Here the coefficient is violently nonlinear:
\alpha(0)=h_0=100,\qquad \alpha(T_\infty)=\frac{100}{1+400}\approx 0.25
Two orders of magnitude swing over the operating range, plus large driving gap (T_\infty=20 vs IC 0). The frozen-coefficient lag over-corrects each step: guess low T gives huge \alpha, huge flux, overshoots high; high T gives tiny \alpha, undershoots back. So g'(x^*) is negative with magnitude >1 - iteration oscillates and diverges. Plain Picard never converges (that is why undamped test caps at 200 iters yet still fails).
Under-relaxation blends old and raw:
x_{k+1}=(1-\beta)\,x_k+\beta\,g(x_k),\qquad \beta\in(0,1]
Effective map derivative becomes \;1-\beta+\beta\,g'\;. Picking \beta small enough pulls the magnitude below 1:
|1-\beta(1-g')|<1
For oscillatory divergence (g'<-1) this needs \beta<\tfrac{2}{1-g'}. \beta=0.39 sits inside that band; \beta=1 does not. Damping trades speed for a guaranteed contraction.
- Tests covering your feature were added?
- Any new feature or behaviour change was documented?
- Add this merge request to the changelog