anchor.diff:    real(x_precision)                    :: delta_S_max, delta_T_max, S_rel_diff
     real(x_precision)                    :: t
     real(x_precision)                    :: dt_nu, dt_T, pf
  +  integer :: p
  +  real(x_precision), dimension(n_cell) :: mean
     logical                              :: T_converged
     logical                              :: unstable, wasunstable
   
     real(x_precision), dimension(n_cell) :: dist_crit
   
     call getarg(1, arg)
  -
  +  mean       = 0._x_precision
     !----------------------------------------------
     ! Convergence criteria for S and T
     !----------------------------------------------
  
positive.code:         ! Mdot kick
          ! Increase Mdot at r_max if S is stalled
          !
          ! the convergence parameter is the maximum relative variation of
          ! s. We give a kick when it's below a constant.
          ! Since it depends linearly on the time, we take it into account
          !----------------------------------------------
          S_rel_diff = maxval(dabs(prev_S - s%S) / s%S / dt_nu)
          if (S_rel_diff < delta_S_max) then
             s%Mdot(n_cell) = s%Mdot(n_cell) * 2._x_precision
             print*, 'Mdot kick! YOLOOOOO', s%Mdot(n_cell)
          else
             ! print*, S_rel_diff
          end if
  
       end if
  
       ! Recompute variables when the system is stable
  
negative.code:         end if
  
       else
          S_steps = S_steps + 1 
          if (wasunstable) then
             print*, 'Switched back to implicit mode!'
             s%Mdot(n_cell) = s%Mdot(n_cell) / 4._x_precision
  