1import json,time
 2import numpy as np
 3L,W=4,2
 4
 5def make(seed,K):
 6 r=np.random.default_rng(seed); E=[]
 7 for t in range(L):
 8  src=range(t*W,(t+1)*W) if t else [-1]; dst=range((t+1)*W,(t+2)*W) if t<L-1 else [-2]
 9  E += [(a,b) for a in src for b in dst]
10 F=[]
11 for _ in E:
12  lo=np.maximum(.03,r.normal(1.2,.25,K)); hi=lo+r.uniform(.05,.75,K); m=r.dirichlet(np.ones(K)); F.append((lo,hi,m))
13 return E,F
14
15def allroutes(E):
16 ix={p:i for i,p in enumerate(E)}; R=[]
17 for q in np.ndindex(*(W,)*(L-1)):
18  ns=[-1]+[t*W+q[t-1] for t in range(1,L)]+[-2]; x=np.zeros(len(E),dtype=np.int8)
19  for a,b in zip(ns[:-1],ns[1:]): x[ix[a,b]]=1
20  R.append(x)
21 return np.array(R)
22
23def base(F): return np.array([np.dot(m,lo) for lo,hi,m in F])
24def blockcost(lo,m,g):
25 # singleton is exact; merged block is conservatively raised to max lower.
26 mass=m[g].sum(); raw=np.dot(m[g],lo[g])
27 return float(raw if len(g)==1 else max(raw,mass*max(lo[g])))
28
29def compress(F,X,target,mode):
30 G=[[[i] for i in range(len(F[e][0]))] for e in range(len(F))]
31 c=base(F); x0=X[np.argmin(X@c)]; total=sum(map(len,G[0]))*len(G)
32 while sum(len(g) for z in G for g in z)>target*len(G):
33  best=(1e99,None,None)
34  for e,gs in enumerate(G):
35   lo,hi,m=F[e]
36   if len(gs)<2: continue
37   old={id(g):blockcost(lo,m,g) for g in gs}
38   for i in range(len(gs)):
39    for j in range(i+1,len(gs)):
40     a,b=gs[i],gs[j]; merged=a+b; d=blockcost(lo,m,merged)-old[id(a)]-old[id(b)]
41     if mode=='regret': s=d*(1 if x0[e] else 0)+.05*d
42     elif mode=='random': s=float((e*997+i*31+j)%1000)
43     elif mode=='mass': s=-(m[a].sum()+m[b].sum())
44     else:
45      span=max(max(hi[a]),max(hi[b]))-min(min(lo[a]),min(lo[b])); overlap=max(0,min(max(hi[a]),max(hi[b]))-max(min(lo[a]),min(lo[b])))
46      s=-overlap/max(span,1e-9)
47     if s<best[0]: best=(s,e,(i,j,d))
48  if best[1] is None: break
49  _,e,(i,j,_)=best; G[e]=G[e][:i]+[G[e][i]+G[e][j]]+G[e][i+1:j]+G[e][j+1:]
50 return np.array([sum(blockcost(lo,m,g) for g in G[e]) for e,(lo,hi,m) in enumerate(F)])
51
52def run():
53 result={'theorem':{},'sweep':{},'baselines':{}}
54 maxviol=0.; ratios=[]
55 for seed in range(10):
56  E,F=make(seed,8); X=allroutes(E); c=base(F); ch=compress(F,X,3,'regret'); delta=X@(ch-c); s=np.argmin(X@c); h=np.argmin(X@ch); reg=(X[h]-X[s])@c
57  maxviol=max(maxviol,float(reg-delta[s])); ratios.append(reg/(delta[s]+1e-12))
58 result['theorem']={'max_bound_violation':float(maxviol),'mean_regret_over_bound':float(np.mean(ratios))}
59 for K in (4,8,12):
60  for target in (1,2,3):
61   a=[]
62   for seed in range(3):
63    E,F=make(100+seed,K); X=allroutes(E); c=base(F); s=np.argmin(X@c); ch=compress(F,X,target,'regret'); h=np.argmin(X@ch); d=X@(ch-c); a.append(((X[h]-X[s])@c,h!=s,d[s],sum(ch-c)))
64   result['sweep'][f'K{K}_to{target}']={'regret':float(np.mean([z[0] for z in a])),'flip':float(np.mean([z[1] for z in a])),'bound':float(np.mean([z[2] for z in a])),'inflation':float(np.mean([z[3] for z in a]))}
65 for mode in ('regret','jaccard','random'):
66  a=[]
67  for seed in range(5):
68   E,F=make(900+seed,8); X=allroutes(E); c=base(F); s=np.argmin(X@c); ch=compress(F,X,2,mode); h=np.argmin(X@ch); a.append(((X[h]-X[s])@c,h!=s))
69  result['baselines'][mode]={'regret':float(np.mean([z[0] for z in a])),'flip':float(np.mean([z[1] for z in a]))}
70 print(json.dumps(result,indent=2))
71if __name__=='__main__': run()