|
| 1 | +#!/usr/bin/env python3 |
| 2 | +"""Optimize parameters for string tension, screening length, and spectral gap.""" |
| 3 | +from __future__ import annotations |
| 4 | + |
| 5 | +import math |
| 6 | +import time |
| 7 | + |
| 8 | +import torch |
| 9 | + |
| 10 | +from fragile.physics.app.coupling_diagnostics import ( |
| 11 | + CouplingDiagnosticsConfig, |
| 12 | + compute_coupling_diagnostics, |
| 13 | +) |
| 14 | +from fragile.physics.fractal_gas.cloning import CloneOperator |
| 15 | +from fragile.physics.fractal_gas.euclidean_gas import EuclideanGas |
| 16 | +from fragile.physics.fractal_gas.fitness import FitnessOperator |
| 17 | +from fragile.physics.fractal_gas.kinetic_operator import KineticOperator |
| 18 | + |
| 19 | + |
| 20 | +def run_and_diagnose(params: dict, seed: int = 42) -> dict: |
| 21 | + p = params |
| 22 | + t0 = time.monotonic() |
| 23 | + |
| 24 | + kinetic_op = KineticOperator( |
| 25 | + gamma=float(p["gamma"]), beta=float(p["beta"]), |
| 26 | + delta_t=float(p["delta_t"]), temperature=float(p["temperature"]), |
| 27 | + nu=float(p["nu"]), use_viscous_coupling=True, |
| 28 | + viscous_length_scale=float(p["viscous_length_scale"]), |
| 29 | + viscous_neighbor_weighting=str(p["viscous_neighbor_weighting"]), |
| 30 | + beta_curl=float(p["beta_curl"]), |
| 31 | + ) |
| 32 | + if p.get("auto_thermostat", True): |
| 33 | + kinetic_op.auto_thermostat = True |
| 34 | + |
| 35 | + cloning = CloneOperator( |
| 36 | + p_max=float(p["p_max"]), epsilon_clone=float(p["epsilon_clone"]), |
| 37 | + sigma_x=float(p["sigma_x"]), alpha_restitution=float(p["alpha_restitution"]), |
| 38 | + ) |
| 39 | + fitness_op = FitnessOperator( |
| 40 | + alpha=float(p["fitness_alpha"]), beta=float(p["fitness_beta"]), |
| 41 | + eta=float(p["eta"]), sigma_min=float(p["sigma_min"]), A=float(p["A"]), |
| 42 | + ) |
| 43 | + |
| 44 | + gas = EuclideanGas( |
| 45 | + N=int(p["N"]), d=int(p["d"]), |
| 46 | + kinetic_op=kinetic_op, cloning=cloning, fitness_op=fitness_op, |
| 47 | + device=torch.device("cpu"), dtype="float32", |
| 48 | + clone_every=int(p["clone_every"]), |
| 49 | + neighbor_graph_update_every=int(p["neighbor_graph_update_every"]), |
| 50 | + neighbor_weight_modes=list(p["neighbor_weight_modes"]), |
| 51 | + tessellation_timing="after_cloning", |
| 52 | + ) |
| 53 | + |
| 54 | + N, d = int(p["N"]), int(p["d"]) |
| 55 | + x_init = torch.randn(N, d) * float(p["init_spread"]) + float(p["init_offset"]) |
| 56 | + v_init = torch.randn(N, d) * float(p["init_velocity_scale"]) |
| 57 | + |
| 58 | + history = gas.run( |
| 59 | + int(p["n_steps"]), x_init=x_init, v_init=v_init, |
| 60 | + record_every=int(p["record_every"]), seed=seed, show_progress=False, |
| 61 | + ) |
| 62 | + |
| 63 | + config = CouplingDiagnosticsConfig(warmup_fraction=0.15) |
| 64 | + output = compute_coupling_diagnostics(history, config=config) |
| 65 | + duration = time.monotonic() - t0 |
| 66 | + |
| 67 | + summary = {} |
| 68 | + for k, v in output.summary.items(): |
| 69 | + if isinstance(v, float) and (math.isnan(v) or math.isinf(v)): |
| 70 | + summary[k] = None |
| 71 | + else: |
| 72 | + summary[k] = v |
| 73 | + return {"summary": summary, "duration": round(duration, 1)} |
| 74 | + |
| 75 | + |
| 76 | +def defaults(n_steps=1000, N=500): |
| 77 | + return { |
| 78 | + "N": N, "d": 3, "n_steps": n_steps, "record_every": 1, |
| 79 | + "init_offset": 0.0, "init_spread": 0.0, "init_velocity_scale": 0.0, |
| 80 | + "gamma": 1.0, "beta": 1.0, "auto_thermostat": True, |
| 81 | + "delta_t": 0.01, "temperature": 0.5, "nu": 1.0, |
| 82 | + "use_viscous_coupling": True, "viscous_length_scale": 1.0, |
| 83 | + "viscous_neighbor_weighting": "riemannian_kernel_volume", |
| 84 | + "beta_curl": 1.0, |
| 85 | + "p_max": 1.0, "epsilon_clone": 1e-6, "sigma_x": 0.01, |
| 86 | + "alpha_restitution": 1.0, |
| 87 | + "fitness_alpha": 1.0, "fitness_beta": 1.0, |
| 88 | + "eta": 0.0, "sigma_min": 0.0, "A": 2.0, |
| 89 | + "neighbor_graph_update_every": 1, |
| 90 | + "neighbor_weight_modes": ["inverse_riemannian_distance", "kernel", "riemannian_kernel_volume"], |
| 91 | + "clone_every": 1, "dtype": "float32", |
| 92 | + } |
| 93 | + |
| 94 | + |
| 95 | +def fmt(v, w=8): |
| 96 | + if v is None: |
| 97 | + return " " * (w - 4) + "None" |
| 98 | + return f"{v:{w}.4f}" |
| 99 | + |
| 100 | + |
| 101 | +def report(name: str, result: dict): |
| 102 | + s = result["summary"] |
| 103 | + sigma = s.get("string_tension_sigma") |
| 104 | + xi = s.get("screening_length_xi") |
| 105 | + asym = s.get("re_im_asymmetry_mean") |
| 106 | + score = s.get("regime_score") |
| 107 | + poly = s.get("polyakov_abs") |
| 108 | + sg_fiedler = s.get("spectral_gap_fiedler") |
| 109 | + sg_ac = s.get("spectral_gap_autocorrelation") |
| 110 | + sg_tm = s.get("spectral_gap_transfer_matrix") |
| 111 | + dur = result["duration"] |
| 112 | + |
| 113 | + print( |
| 114 | + f" {name:40s} | σ={fmt(sigma)} | ξ={fmt(xi)} | asym={fmt(asym)} | " |
| 115 | + f"poly={fmt(poly)} | fiedler={fmt(sg_fiedler)} | ac_gap={fmt(sg_ac)} | " |
| 116 | + f"tm_gap={fmt(sg_tm)} | score={fmt(score)} | {dur:.0f}s", |
| 117 | + flush=True, |
| 118 | + ) |
| 119 | + |
| 120 | + |
| 121 | +if __name__ == "__main__": |
| 122 | + import warnings |
| 123 | + warnings.filterwarnings("ignore") |
| 124 | + |
| 125 | + header = ( |
| 126 | + f" {'name':40s} | {'σ':>8s} | {'ξ':>8s} | {'asym':>8s} | " |
| 127 | + f"{'poly':>8s} | {'fiedler':>8s} | {'ac_gap':>8s} | " |
| 128 | + f"{'tm_gap':>8s} | {'score':>8s} | time" |
| 129 | + ) |
| 130 | + |
| 131 | + # Focused runs with 500 walkers, 1000 steps |
| 132 | + print("=" * 160) |
| 133 | + print("PARAMETER OPTIMIZATION (500 walkers, 1000 steps)") |
| 134 | + print("=" * 160) |
| 135 | + print(header) |
| 136 | + print("-" * 160) |
| 137 | + |
| 138 | + combos = [ |
| 139 | + # Baseline |
| 140 | + ("default", {}), |
| 141 | + # Temperature sweep (moderate range, keeping asym healthy) |
| 142 | + ("T=0.5", {"temperature": 0.5}), |
| 143 | + ("T=0.7", {"temperature": 0.7}), |
| 144 | + ("T=1.0", {"temperature": 1.0}), |
| 145 | + # Best individual effects combined |
| 146 | + ("T=0.5_nu=2", {"temperature": 0.5, "nu": 2.0}), |
| 147 | + ("T=0.7_nu=2", {"temperature": 0.7, "nu": 2.0}), |
| 148 | + ("T=1.0_nu=2", {"temperature": 1.0, "nu": 2.0}), |
| 149 | + ("T=0.5_fit_a=2", {"temperature": 0.5, "fitness_alpha": 2.0}), |
| 150 | + ("T=0.7_fit_a=2", {"temperature": 0.7, "fitness_alpha": 2.0}), |
| 151 | + # viscous_length + nu |
| 152 | + ("T=0.5_nu=2_vl=2", {"temperature": 0.5, "nu": 2.0, "viscous_length_scale": 2.0}), |
| 153 | + ("T=0.7_nu=2_vl=0.1", {"temperature": 0.7, "nu": 2.0, "viscous_length_scale": 0.1}), |
| 154 | + # beta_curl combos |
| 155 | + ("T=0.5_bc=2_nu=2", {"temperature": 0.5, "beta_curl": 2.0, "nu": 2.0}), |
| 156 | + ("T=0.7_bc=0.5_nu=2", {"temperature": 0.7, "beta_curl": 0.5, "nu": 2.0}), |
| 157 | + # Moderate push on multiple fronts |
| 158 | + ("T=0.5_nu=2_fit_a=2", {"temperature": 0.5, "nu": 2.0, "fitness_alpha": 2.0}), |
| 159 | + ("T=0.7_nu=2_fit_a=2", {"temperature": 0.7, "nu": 2.0, "fitness_alpha": 2.0}), |
| 160 | + ("T=0.5_nu=2_fit_a=2_vl=2", {"temperature": 0.5, "nu": 2.0, "fitness_alpha": 2.0, "viscous_length_scale": 2.0}), |
| 161 | + # Higher nu with moderate T |
| 162 | + ("T=0.5_nu=3", {"temperature": 0.5, "nu": 3.0}), |
| 163 | + ("T=0.7_nu=3", {"temperature": 0.7, "nu": 3.0}), |
| 164 | + ] |
| 165 | + |
| 166 | + for name, overrides in combos: |
| 167 | + p = defaults() |
| 168 | + p.update(overrides) |
| 169 | + try: |
| 170 | + r = run_and_diagnose(p) |
| 171 | + report(name, r) |
| 172 | + except Exception as e: |
| 173 | + print(f" {name:40s} | ERROR: {e}", flush=True) |
0 commit comments