import json, math, random, time from pathlib import Path import numpy as np SEED = 3070 np.random.seed(SEED); random.seed(SEED) def math_checks(): # Piecewise-constant temperature field: the stated integral reduces to a # harmonic mean over the uphill half-period. T = np.array([0.5, 1.0, 2.0, 4.0]) Tact = len(T) / np.sum(1.0 / T) numerical = 2.0 / np.mean(1.0 / T) / 2.0 # same harmonic mean, explicit average # T <- T r^k gives log(r_next)=(1-k)log(r) when U is locally fixed. rs = np.array([0.05, 0.2, 3.0, 20.0]) k = 0.08 next_rs = rs ** (1-k) contraction = np.max(np.abs(np.log(next_rs))/np.abs(np.log(rs))) # Arrhenius claim: log escape probability is linear in -U/T. ratios = np.linspace(0.3, 3.0, 12) p = np.exp(-ratios) slope = np.polyfit(ratios, np.log(p), 1)[0] return {"Tact_harmonic": float(Tact), "Tact_formula_error": float(abs(Tact-numerical)), "controller_log_contraction": float(contraction), "arrhenius_log_slope": float(slope)} def double_well_check(): # V(x)=(x^2-1)^2 has minima +/-1 and barrier U=1 at x=0. # Parallel trajectories reduce rare-transition sampling noise. rng = np.random.default_rng(SEED) dt, steps, burn, ntraj = 0.002, 100000, 5000, 256 temps = np.array([0.18, 0.24, 0.32, 0.42, 0.55]) rates = [] for temp in temps: x = np.full(ntraj, -1.0) side = np.full(ntraj, -1) transitions = 0 for i in range(steps): grad = 4*x*(x*x - 1) x += -grad*dt + math.sqrt(2*temp*dt)*rng.normal(size=ntraj) if i >= burn: newside = np.where(x > 0, 1, -1) transitions += int(np.count_nonzero(newside != side)) side = newside rates.append(transitions/(ntraj*(steps-burn))) rates = np.maximum(np.asarray(rates), 1e-12) fit = np.polyfit(1/temps, np.log(rates), 1) return {"temperatures": temps.tolist(), "transition_rates": rates.tolist(), "log_rate_vs_inverse_temperature_slope": float(fit[0]), "expected_negative_barrier_sign": bool(fit[0] < 0)} def make_data(n=512): rng=np.random.default_rng(SEED) x=rng.normal(size=(n,2)).astype(np.float32) y=((x[:,0]**2 + x[:,1]**2 + .35*x[:,0]) > 1.0).astype(np.int64) # introduce a nontrivial but tiny classification task return x,y def run_mlp(controller, steps=500): import torch torch.manual_seed(SEED); np.random.seed(SEED) device = "cuda" if torch.cuda.is_available() else "cpu" try: dev=torch.device(device) X,y=make_data(); X=torch.tensor(X,device=dev); y=torch.tensor(y,device=dev) model=torch.nn.Sequential(torch.nn.Linear(2,16),torch.nn.Tanh(),torch.nn.Linear(16,2)).to(dev) lossfn=torch.nn.CrossEntropyLoss() lr=0.035; T=0.0025; batch=64; d=sum(p.numel() for p in model.parameters()) losses=[]; controls=[]; transitions=0; last_pred=None for step in range(steps): ix=torch.randint(0,len(X),(batch,),device=dev) logits=model(X[ix]); loss=lossfn(logits,y[ix]) grads=torch.autograd.grad(loss,tuple(model.parameters()),create_graph=False) with torch.no_grad(): for p,g in zip(model.parameters(),grads): # Langevinized SGD. T is the injected update temperature. p.add_(-lr*g + math.sqrt(2*max(T,1e-12))*lr*torch.randn_like(p)) losses.append(float(loss.detach().cpu())) if step % 20 == 0: with torch.no_grad(): pred=model(X).argmax(1) basin=int((pred==y).sum().item() > len(y)*.5) if last_pred is not None and basin != last_pred: transitions += 1 last_pred=basin if controller and step % 50 == 0: # Gradient covariance from 8 independent minibatches; P=lr I. gs=[] for _ in range(8): j=torch.randint(0,len(X),(batch,),device=dev) ll=lossfn(model(X[j]),y[j]) gg=torch.autograd.grad(ll,tuple(model.parameters())) gs.append(torch.cat([z.detach().reshape(-1) for z in gg])) G=torch.stack(gs); cov_trace=((G-G.mean(0))**2).sum(1).mean() That=float((lr*lr*cov_trace/(2*d)).detach().cpu()) # Two perturbed endpoints and a linear interpolation barrier. base=[p.detach().clone() for p in model.parameters()] eps=[torch.randn_like(p)*math.sqrt(2*max(T,1e-12))*0.7 for p in model.parameters()] vals=[] for alpha in torch.linspace(0,1,7,device=dev): with torch.no_grad(): for p,b,e in zip(model.parameters(),base,eps): p.copy_(b+(2*alpha-1)*e) vals.append(float(lossfn(model(X),y).detach().cpu())) with torch.no_grad(): for p,b in zip(model.parameters(),base): p.copy_(b) U=max(0.0,max(vals)-vals[3]) effective=max(That, T) r=U/max(effective,1e-10) old=T T=float(np.clip(T*math.exp(0.08*np.clip(math.log(max(r,1e-8)),-4,4)),1e-5,0.08)) controls.append((step,U,That,r,old,T)) with torch.no_grad(): final_loss=float(lossfn(model(X),y).cpu()); acc=float((model(X).argmax(1)==y).float().mean().cpu()) return {"final_loss":final_loss,"accuracy":acc,"mean_last100_loss":float(np.mean(losses[-100:])), "transitions":transitions,"final_temperature":T,"controls":controls,"device":str(dev)} except Exception as e: # Explicit CPU fallback as required for shared-GPU robustness. if device == "cuda": torch.cuda.empty_cache() torch.cuda.is_available=lambda: False return run_mlp(controller,steps) raise def main(): t=time.time() out={"seed":SEED,"math":math_checks(),"double_well":double_well_check(), "mlp_fixed":run_mlp(False),"mlp_barrier_matching":run_mlp(True),"seconds":time.time()-t} # Keep the report compact while preserving all controller observations. Path("results.json").write_text(json.dumps(out,indent=2)) print(json.dumps({k:v for k,v in out.items() if k not in ("mlp_fixed","mlp_barrier_matching")},indent=2)) for k in ("mlp_fixed","mlp_barrier_matching"): v=out[k]; print(k, {z:v[z] for z in ("final_loss","accuracy","mean_last100_loss","transitions","final_temperature","device")}) if __name__ == "__main__": main()