Line search in Integrator for 3D-0D coupling

Current situation

In structural mechanics simulations (and electromechanics), implicit coupling of 3D and 0D models might face convergence issues in the nonlinear solver in presence of valves in the 0D model (be them piecewise or tanh valves).

The convergence issue typically manifests as stagnation in the nonlinear solver. The residual does not decrease with iterations, but instead bounces between two different values until hitting the maximum iteration number.

This behavior can be understood in terms of the discontinuous (or nearly discontinuous) tangent of the valve model. Near to the switching point, the linear approximation inherent to Newton's method is a poor approximation of the valve relation, so that Newton iterations keep switching between open and close state, never converging to the solution.

This is also discussed by @dseyler in this svZeroDSolver issue. The proposed solution there (only implemented for piecewise valves, but extended to tanh valves here), is to use explicit discretization for valves (also called "freezing"): the valve resistance is evaluated using the old pressure jump, and so for a fixed time the valve tangent is fixed.

Valve freezing has two potential issues as far as I can see, one minor and one more critical.

1. The update lag may create artificial overshoots at the instants when the valve switches; this is a minor issue, and is expected to vanish as $\Delta t \to 0$.
2. The explicit discretization may give rise to instabilities. This is particularly evident for tanh valves (see the attached picture below). These instabilities, in my experiments, get slightly worse as $\Delta t$ is reduced.

Pressure and volume traces with tanh valves and freezing

Proposed solution

I propose to address this by implementing backtracking/line search in the nonlinear solver, at least for 3D-0D coupled structural problems. This is expected to fix the issue because it prevents the inaccurate tangent for the valve model from overshooting the Newton iteration.

The gist of it is that, after solving one iteration of the nonlinear system, the residual is evaluated. If it hasn't reduced enough compared to previous iteration, then the increment is rescaled by a factor $\alpha < 1$.

I have experimentally implemented this in this branch, and it seems to solve the issue (see the attached figure below, and notice the absence of oscillations). Line search only slightly increases the nonlinear iteration count (for most of the simulation no step is rejected).

Pressure and volume traces with tanh valves and line search (without freezing)

The code linked above is Claude-generated, and it is very struct-specific. Since it is in the Integrator class, which should be equation-agnostic, I will try to make it a bit more general and overall cleaner before proceeding to a PR.

FYI @dseyler @aabrown100-git @Eleven7825