"""Check the staircase lower bound (Proposition in Section 4.4).

f(x) = sum_i hub_{s_i}(x_i), hub_s' (t) = clip(t, -s, s)  (convex, 1-smooth, minimizer 0),
x_0 = (a, ..., a), s_i = a / 2^i.  GD with step 1 moves coordinate i linearly to 0, arriving at k = 2^i.
Claim: with eps = a/2 the sequence has at least d eps-fluctuations, while ||x_0 - 0|| = a sqrt(d) = 2 sqrt(d) eps.
We count the maximal number of eps-fluctuations exactly by the greedy rule
(earliest possible right end point), which is optimal for this interval-scheduling structure.
"""
import numpy as np


def trajectory(d, a=1.0):
    s = a / 2.0 ** np.arange(1, d + 1)
    x = np.full(d, a)
    xs = [x.copy()]
    while np.any(x != 0):
        x = x - np.clip(x, -s, s)
        xs.append(x.copy())
    return np.array(xs)


def max_fluctuations(X, eps):
    """Greedy: scan j, keep the set of candidate left points since the last cut; cut at the first j with
    max_i ||x_i - x_j|| >= eps. Optimal for chains i1<j1<=i2<j2<=..."""
    count, start, n = 0, 0, len(X)
    j = start + 1
    while j < n:
        dmax = np.max(np.linalg.norm(X[start:j] - X[j], axis=1))
        if dmax >= eps - 1e-12:
            count += 1
            start = j
        j += 1
    return count


if __name__ == "__main__":
    for d in range(1, 15):
        X = trajectory(d)
        eps = 0.5
        D = np.linalg.norm(X[0])
        print(f"d={d:2d} steps={len(X)-1:6d} D/eps={D/eps:6.3f} (2sqrt d={2*np.sqrt(d):6.3f}) "
              f"fluct={max_fluctuations(X, eps):3d} pathlen/D={np.sum(np.linalg.norm(np.diff(X,axis=0),axis=1))/D:6.3f}")