Add an assertion to the computation of the SPD factor.

[bangerth]

This is a follow-up to #5671 and codifies an assumption in the design of the algorithm by @YiminJin .

@naliboff FYI.

[gassmoeller]

It looks like you got a merge conflict and a bunch of newton solver tests fail.

[bangerth]

@YiminJin I could use your help here. This patch adds an assertion that the computation of the SPD factor is really only called when eps is parallel to d eta/d eps. It fails in a bunch of our nonlinear tests, such as nonlinear_channel_flow_tractions_Newton_Stokes_no_deviator, but I don't understand why. In that test, for example, we use the benchmarks/newton_solver_benchmark_set/nonlinear_channel_flow/simple_nonlinear.cc material model that computes the viscosity as a function of only e_dot_ii = eps:eps/2. The derivative of such a quantity should always be parallel to eps.

We compute the derivative with Utilities::derivative_of_weighted_p_norm_average() from the array composition_viscosities_derivatives. I don't see anything that is particularly wrong here, but perhaps I'm missing something. Or does the averaging just not preserve the property we need?

[YiminJin]

@bangerth This is because I forced the assembler to use the deviatoric strain rate (line 179 and 404 in source/simulator/assemblers/newton_stokes.cc). The exception handler will certainly be triggered when a material model calculates the effective viscosity with the full strain rate. I made this modification because the newton solver was meant to serve the only two rheology models that calculate the viscosity derivatives (ViscoPlastic and CompositeViscoPlastic), and both of them use the deviatoric strain rate when calculating the effective viscosity (line 135 in source/material_model/rheology/visco_plastic.cc and line 133 in source/material_model/rheology/composite_visco_plastic.cc).

Personally I think it is better to use the deviatoric strain rate when calculating the viscosity, since the viscosity is only related with the deviatoric parts of strain and stress in most cases. But if we are really in need of calculating the viscosity with the full strain rate, we may distinguish such cases by passing a boolean to the assembler (somewhere in Scratch), or add an entry in the parameter list. What do you think?

[YiminJin]

I have opened a PR for the SPD factor in general cases (#5754 ). This can circumvent the problem, but I doubt if it is a good idea, because the convergence rate of the Newton solver is very sensitive to the descending direction, and a small deviation in the descending direction may significantly increase the number of iterations required for convergence.