#!/usr/bin/env python3
"""The difficulty controller against a wide DAG: the coupling between red share, chain-step spacing and the estimator
(Horizon network lane, 6 October 2026).

The rule as the fork runs it (vendor/igneum-node release-0.3.14-node, consensus/src/processes/difficulty.rs
`igneum_difficulty_bits`, consensus/core/src/igneum.rs `difficulty`): the lanes walk the SELECTED CHAIN; every step
carries `work = blue_work(b) - blue_work(selected parent)` (the mergeset BLUES' work only, red work is not in blue
work) and `solvetime = min(c(b) - c(p), CAP_BLOCKS x T)` with CAP_BLOCKS = 20 (2 s at T = 100 ms); the target moves at
most HARDEN_PERCENT 3 (difficulty up) or EASE_PERCENT 10 (difficulty down) per chain block (spec 02 2.3).
Kaspa's rule (vendor/rusty-kaspa consensus/src/processes/window.rs `push_mergeset`, difficulty.rs
`calculate_difficulty_bits`): the DAA window holds every mergeset block, blue AND red, above the window's blue-score
floor, and the new target is the window's average target x measured span / (window size x T): the whole-DAG rate
over the real span.

The abstract DAG regime (no block-level simulation; the inputs are the propagation lane's):
  spacing = max(1 / lambda, d)    seconds between selected-chain blocks (every block a chain block when the DAG is
                                  narrow; one chain block per propagation delay when it is wide)
  m = lambda x spacing            mergeset per chain block
  blues = min(m, k + 1)           the k-cluster admits at most k blues beside the selected parent per chain step
                                  (protocol.rs); the rest of the mergeset is red: r = 1 - blues / m
  (the PHANTOM tail P(Poisson(2 d lambda) > k) that k is derived from, bps.rs `calculate_ghostdag_k`, says when
  reds appear at all; it is printed beside r)
Estimators of the network hash H (work per block = D, the difficulty):
  igneum   = blues x D / min(spacing, cap)                    (blue work over the capped chain step)
  whole    = m x D / spacing                                  (every mergeset block's work over the real span) = H
  blue-rc  = blues x D / spacing / (1 - r_obs)                (blue rate with the observed red share corrected) = H
The controller sets the next D toward estimate x T through the 3% / 10% clamps once per chain block.

Run: python3 sim/horizon/network/controller.py [--bps 10 --k 124 --delays 0.3,1,2,5,10 --minutes 30 --lambda0 26]
"""
import argparse, math

def pois_tail_gt(k, x):
    """P(Poisson(x) > k)"""
    if x <= 0: return 0.0
    # sum P(N <= k) in log space
    lp = -x; s = math.exp(lp)
    for i in range(1, k + 1):
        lp += math.log(x / i); s += math.exp(lp)
        if lp < -745: break
    return max(0.0, min(1.0, 1.0 - s))

def run(est, bps, k, d, minutes, lam0, H=1e9):
    T = 1.0 / bps; cap = 20 * T
    D = H / lam0  # difficulty giving the initial production rate
    t = 0.0; rows = []
    while t < minutes * 60:
        lam = H / D
        x = 2 * d * lam
        spacing = max(1.0 / lam, d)
        m = lam * spacing
        blues = min(m, k + 1)
        r = 1 - blues / m
        if est == "igneum":
            h = blues * D / min(spacing, cap)
        elif est == "whole":
            h = m * D / spacing
        else:
            h = blues * D / spacing / max(1e-9, 1 - r)
        want = h * T
        D_new = min(D * 1.03, max(D * 0.90, want))  # harden at most 3%, ease at most 10%
        rows.append((t, lam, r, m, blues, h / H, D, pois_tail_gt(k, x)))
        D = D_new; t += spacing
    return rows

def main():
    ap = argparse.ArgumentParser()
    ap.add_argument("--bps", type=float, default=10); ap.add_argument("--k", type=int, default=124)
    ap.add_argument("--delays", default="0.3,1,2,5,10"); ap.add_argument("--minutes", type=float, default=30)
    ap.add_argument("--lambda0", type=float, default=26.0, help="initial production rate, blocks/s (run A's genesis burst was 26)")
    ap.add_argument("--out", default=None)
    a = ap.parse_args()
    lines = []
    def emit(s): print(s); lines.append(s)
    emit(f"# controller.py: bps {a.bps:g}, k {a.k}, cap {20/a.bps:g} s, start at {a.lambda0:g} blocks/s, {a.minutes:g} min; settled = mean of the last 5 minutes")
    emit("")
    emit("| estimator | delay d, s | DAG rate settled, blocks/s | blue rate, blocks/s | red share | P(anticone > k) | mergeset per chain block | estimate / true hash | difficulty vs correct | time to within 10% of target DAG rate, min |")
    emit("|---|---|---|---|---|---|---|---|---|---|")
    for est in ("igneum", "whole", "blue-rc"):
        for d in [float(v) for v in a.delays.split(",")]:
            rows = run(est, a.bps, a.k, d, a.minutes, a.lambda0)
            tail = [r for r in rows if r[0] >= (a.minutes - 5) * 60] or rows[-5:]
            lam = sum(r[1] for r in tail) / len(tail); red = sum(r[2] for r in tail) / len(tail)
            m = sum(r[3] for r in tail) / len(tail); ratio = sum(r[5] for r in tail) / len(tail)
            Dcorrect = 1e9 / a.bps; Dratio = sum(r[6] for r in tail) / len(tail) / Dcorrect
            within = next((r[0] / 60 for r in rows if abs(r[1] - a.bps) <= 0.1 * a.bps), None)
            wtxt = f"{within:.1f}" if within is not None else "never"
            tail_p = sum(r[7] for r in tail) / len(tail)
            emit(f"| {est} | {d:g} | {lam:.2f} | {lam*(1-red):.2f} | {100*red:.1f}% | {tail_p:.2f} | {m:.1f} | {ratio:.2f}x | {Dratio:.2f}x | {wtxt} |")
    if a.out:
        with open(a.out, "w") as f: f.write("\n".join(lines) + "\n")

if __name__ == "__main__":
    main()