1import json
 2import numpy as np
 3from anytime_pd import solve_full, hierarchy, stop_level
 4
 5
 6def main():
 7    rng = np.random.default_rng(2930)
 8    n = 6
 9    # Box constraints plus two coupling inequalities create a nontrivial LP.
10    I = np.eye(n)
11    A = np.vstack([I, -I, rng.normal(size=(2, n))])
12    b = np.concatenate([np.ones(n), np.ones(n), np.array([0.85, 0.65])])
13    c = np.array([0.80, 0.55, 0.35, 0.25, 0.15, 0.10])
14
15    Lstar, Ustar, vstar, ystar = solve_full(A, b, c)
16    # Coordinate subspaces are nested and reach the full primal space.
17    bases = [np.eye(n)[:, :k] for k in range(1, n + 1)]
18    # ystar is already a feasible dual certificate. Adding nonnegative
19    # coordinate directions gives nested dual feasible sets.
20    dual_bases = [np.column_stack([ystar, np.eye(A.shape[0])[:, :k]])
21                  for k in range(1, A.shape[0] + 1)]
22    records = hierarchy(A, b, c, bases, dual_bases, full_opt=Lstar)
23
24    lower = [r.lower for r in records]
25    upper = [r.upper for r in records]
26    widths = [r.width for r in records]
27    checks = {
28        "all_cover_full": all(r.covers_full for r in records),
29        "all_primal_feasible": all(r.primal_residual <= 1e-7 for r in records),
30        "all_dual_feasible": all(r.dual_eq_residual <= 1e-7 and r.dual_nonnegative for r in records),
31        "lower_monotone": all(x <= y + 1e-8 for x, y in zip(lower, lower[1:])),
32        "upper_monotone": all(x >= y - 1e-8 for x, y in zip(upper, upper[1:])),
33        "width_monotone": all(x >= y - 1e-8 for x, y in zip(widths, widths[1:])),
34        "exact_final_primal": abs(lower[-1] - Lstar) <= 1e-7,
35        "strong_duality": abs(Lstar - Ustar) <= 1e-7,
36    }
37    result = {
38        "full_primal": Lstar,
39        "full_dual": Ustar,
40        "levels": [r.__dict__ for r in records],
41        "stop_level_2pct": stop_level(records, tau=0.02),
42        "checks": checks,
43    }
44    print(json.dumps(result, indent=2))
45    assert all(checks.values())
46
47
48if __name__ == "__main__":
49    main()