docs: def divmod_neg(a, b):
      """Return divmod with closest result to zero"""
      q, r = divmod(a, b)
      # make sure r is close to zero
      sr = np.sign(r)
      if np.abs(r) > b/2:
          q += sr
          r -= b * sr
      return q, r