Bug in Simulator<dim>::assemble_and_solve_stokes() when computing the nonlinear residual?

[bangerth]

In trying to figure out why #5291 does not test cleanly, I've come across an issue in the Simulator<dim>::assemble_and_solve_stokes() function. Here is how it looks, in redux on the current main branch:

template <int dim>
  double Simulator<dim>::assemble_and_solve_stokes (const bool compute_initial_residual,
                                                    double *initial_nonlinear_residual)
  {
[...]

    if (compute_initial_residual)
      {
        Assert(initial_nonlinear_residual != nullptr, ExcInternalError());
        *initial_nonlinear_residual = compute_initial_stokes_residual();
      }
[...]
    if ((initial_nonlinear_residual != nullptr) && (*initial_nonlinear_residual > 0))
      return current_nonlinear_residual / *initial_nonlinear_residual;
    else
      return 0.0;
  }

The bug (?) is that the computation of the nonlinear residual (and putting it into *initial_nonlinear_residual only happens if compute_initial_residual==true. But later on, we use that value.

There are two ways to see this issue:

• Either it's a bug where we use a value we haven't computed.
• Or it's on purpose -- I'm leaning towards that. The idea then is that compute_initial_residual is only true in the very first iteration, and in all other iterations we scale the computed residual by the value computed the very first time around. I think that's what the variable name represents.

I'm going to fix up #5291 based on that second assumption by separating the initial_nonlinear_residual into two variables: One we're asking to compute, and one that is passed in by value.

@MFraters Does that match your interpretation?

[bangerth]

Yea, that second interpretation is right. I was just confused by the naming of the function arguments.