bench/
bench.rs

1//! Ising sweep timing; prints ms/sweep and a state checksum that behavior-neutral
2//! changes must leave unchanged. Configure with environment variables (defaults):
3//! `PEAPODS_L` (128), `PEAPODS_DIM` (2), `PEAPODS_TEMPS` (16), `PEAPODS_REPLICAS` (2),
4//! `PEAPODS_SWEEPS` (50), `PEAPODS_NREAL` (100), `PEAPODS_TMIN` (0.1), `PEAPODS_TMAX`
5//! (5.0, geometric ladder), `PEAPODS_MODE` (cmr | metropolis | pt | sw | wolff | sw_pt |
6//! houdayer | jorg | any overlap build mode string, which runs with PT),
7//! `PEAPODS_SWEEP` (metropolis | gibbs), `PEAPODS_COUPLINGS` (bimodal | ferro | gaussian),
8//! and the flags `PEAPODS_SEQUENTIAL`, `PEAPODS_OVERLAP_WOLFF` (Wolff overlap clusters
9//! instead of SW) and `PEAPODS_GENERIC_LATTICE` (generic neighbor table).
10
11use std::env;
12use std::hash::{DefaultHasher, Hash, Hasher};
13use std::sync::atomic::AtomicBool;
14use std::time::Instant;
15
16use rand::{Rng, RngCore, SeedableRng};
17use rand_xoshiro::Xoshiro256StarStar;
18use spin_sim::config::*;
19use spin_sim::geometry::hypercubic;
20use spin_sim::{run_sweep_parallel, Lattice, Realization};
21
22fn env_usize(name: &str, default: usize) -> usize {
23    env::var(name).map_or(default, |value| {
24        value
25            .parse()
26            .unwrap_or_else(|_| panic!("{name} must be an integer"))
27    })
28}
29
30fn main() {
31    let l = env_usize("PEAPODS_L", 128);
32    let n_dims = env_usize("PEAPODS_DIM", 2);
33    let n_temps = env_usize("PEAPODS_TEMPS", 16);
34    let n_replicas = env_usize("PEAPODS_REPLICAS", 2);
35    let n_sweeps = env_usize("PEAPODS_SWEEPS", 50);
36    let n_realizations = env_usize("PEAPODS_NREAL", 100);
37    let sequential = env::var_os("PEAPODS_SEQUENTIAL").is_some();
38    let overlap_cluster_mode = if env::var_os("PEAPODS_OVERLAP_WOLFF").is_some() {
39        ClusterMode::Wolff
40    } else {
41        ClusterMode::Sw
42    };
43    let generic_lattice = env::var_os("PEAPODS_GENERIC_LATTICE").is_some();
44    let mode = env::var("PEAPODS_MODE").unwrap_or_else(|_| "cmr".to_string());
45    let couplings_kind = env::var("PEAPODS_COUPLINGS").unwrap_or_else(|_| "bimodal".to_string());
46    let sweep_mode = match env::var("PEAPODS_SWEEP").as_deref().unwrap_or("metropolis") {
47        "metropolis" => SweepMode::Metropolis,
48        "gibbs" => SweepMode::Gibbs,
49        other => panic!("unknown PEAPODS_SWEEP '{other}'"),
50    };
51    let shape = vec![l; n_dims];
52    let lattice = if generic_lattice {
53        Lattice::with_offsets(shape, hypercubic(n_dims))
54    } else {
55        Lattice::new(shape)
56    };
57    let n_spins = lattice.n_spins;
58    let n_neighbors = lattice.n_neighbors;
59
60    let t_min: f32 = env::var("PEAPODS_TMIN").map_or(0.1, |v| v.parse().unwrap());
61    let t_max: f32 = env::var("PEAPODS_TMAX").map_or(5.0, |v| v.parse().unwrap());
62    let temps: Vec<f32> = (0..n_temps)
63        .map(|i| t_min * (t_max / t_min).powf(i as f32 / (n_temps.max(2) - 1) as f32))
64        .collect();
65
66    let n_pairs = n_replicas / 2;
67    let n_systems = n_replicas * n_temps;
68    let rngs_per_real = n_systems + n_temps * n_pairs;
69
70    let mut rng = Xoshiro256StarStar::seed_from_u64(0x5eed);
71    let mut realizations = Vec::with_capacity(n_realizations);
72    for r in 0..n_realizations {
73        let couplings: Vec<f32> = (0..n_spins * n_neighbors)
74            .map(|_| match couplings_kind.as_str() {
75                "bimodal" => {
76                    if rng.gen::<bool>() {
77                        1.0
78                    } else {
79                        -1.0
80                    }
81                }
82                "ferro" => 1.0,
83                // Box-Muller keeps the bench free of extra distribution crates.
84                "gaussian" => {
85                    let u1: f64 = 1.0 - rng.gen::<f64>();
86                    let u2: f64 = rng.gen();
87                    ((-2.0 * u1.ln()).sqrt() * (std::f64::consts::TAU * u2).cos()) as f32
88                }
89                other => panic!("unknown PEAPODS_COUPLINGS '{other}'"),
90            })
91            .collect();
92        let base_seed = 42 + (r * rngs_per_real) as u64;
93        realizations.push(Realization::new(
94            &lattice, couplings, &temps, n_replicas, base_seed,
95        ));
96    }
97
98    let interrupted = AtomicBool::new(false);
99    let (cluster_update, pt_interval, overlap_cluster) = match mode.as_str() {
100        "metropolis" => (None, None, None),
101        "pt" => (None, Some(1), None),
102        "sw" => (
103            Some(ClusterConfig {
104                interval: 1,
105                mode: ClusterMode::Sw,
106                action: ClusterAction::Update,
107                collect_stats: false,
108            }),
109            None,
110            None,
111        ),
112        "wolff" => (
113            Some(ClusterConfig {
114                interval: 1,
115                mode: ClusterMode::Wolff,
116                action: ClusterAction::Update,
117                collect_stats: false,
118            }),
119            None,
120            None,
121        ),
122        "cmr" => (
123            None,
124            Some(1),
125            Some(OverlapClusterConfig {
126                interval: 1,
127                modes: vec![OverlapClusterBuildMode::Cmr],
128                cluster_mode: overlap_cluster_mode,
129                action: ClusterAction::Update,
130                collect_stats: false,
131                snapshot_interval: None,
132                max_temperature: None,
133            }),
134        ),
135        "sw_pt" => (
136            Some(ClusterConfig {
137                interval: 1,
138                mode: ClusterMode::Sw,
139                action: ClusterAction::Update,
140                collect_stats: false,
141            }),
142            Some(1),
143            None,
144        ),
145        "houdayer" | "jorg" => (
146            None,
147            Some(1),
148            Some(OverlapClusterConfig {
149                interval: 1,
150                modes: vec![if mode == "jorg" {
151                    OverlapClusterBuildMode::Jorg(2)
152                } else {
153                    OverlapClusterBuildMode::Houdayer(2)
154                }],
155                cluster_mode: overlap_cluster_mode,
156                action: ClusterAction::Update,
157                collect_stats: false,
158                snapshot_interval: None,
159                max_temperature: None,
160            }),
161        ),
162        // Any other overlap build mode string, e.g. "pair4", "jorg4" or "rmc", with PT.
163        other => (
164            None,
165            Some(1),
166            Some(OverlapClusterConfig {
167                interval: 1,
168                modes: parse_overlap_modes(other)
169                    .unwrap_or_else(|err| panic!("unknown PEAPODS_MODE '{other}': {err}")),
170                cluster_mode: overlap_cluster_mode,
171                action: ClusterAction::Update,
172                collect_stats: false,
173                snapshot_interval: None,
174                max_temperature: None,
175            }),
176        ),
177    };
178
179    let config = SimConfig {
180        n_sweeps,
181        warmup_sweeps: 0,
182        sweep_mode,
183        cluster_update,
184        pt_interval,
185        pt_schedule: PtSchedule::SingleRandomEdge,
186        overlap_cluster,
187        autocorrelation_max_lag: None,
188        autocorrelation_backend: AutocorrelationBackend::Ring,
189        sequential,
190        equilibration_diagnostic: false,
191        conditional_link_overlap: true,
192    };
193
194    println!(
195        "Lattice: {l}^{n_dims}  |  Temps: {n_temps}  |  Replicas: {n_replicas}  |  Sweeps: {n_sweeps}  |  Realizations: {n_realizations}"
196    );
197    println!(
198        "Config: {couplings_kind}, mode={mode}, sequential={sequential}, generic_lattice={generic_lattice}"
199    );
200    println!("{}", "-".repeat(70));
201
202    let t0 = Instant::now();
203    let result = run_sweep_parallel(
204        &lattice,
205        &mut realizations,
206        n_replicas,
207        n_temps,
208        &config,
209        &interrupted,
210        &|| {},
211    )
212    .unwrap();
213    let elapsed = t0.elapsed().as_secs_f64();
214
215    let mut state_hash = DefaultHasher::new();
216    for realization in &realizations {
217        realization.spins.hash(&mut state_hash);
218        realization.system_ids.hash(&mut state_hash);
219        for rng in realization.rngs.iter().chain(&realization.pair_rngs) {
220            let mut rng = rng.clone();
221            rng.next_u64().hash(&mut state_hash);
222        }
223    }
224    for values in [
225        &result.mags,
226        &result.mags2,
227        &result.mags4,
228        &result.energies,
229        &result.energies2,
230    ] {
231        for value in values {
232            value.to_bits().hash(&mut state_hash);
233        }
234    }
235    for values in [
236        &result.overlap_stats.overlap,
237        &result.overlap_stats.overlap2,
238        &result.overlap_stats.overlap4,
239        &result.overlap_stats.link_overlap,
240        &result.overlap_stats.link_overlap2,
241        &result.overlap_stats.link_overlap4,
242    ] {
243        for value in values {
244            value.to_bits().hash(&mut state_hash);
245        }
246    }
247    result.overlap_stats.histogram.hash(&mut state_hash);
248    result
249        .overlap_stats
250        .per_sample_histogram
251        .hash(&mut state_hash);
252    for samples in [
253        &result.overlap_stats.ql_at_q_sum,
254        &result.overlap_stats.ql2_at_q_sum,
255    ] {
256        for bins in samples {
257            for value in bins {
258                value.to_bits().hash(&mut state_hash);
259            }
260        }
261    }
262    for samples in [
263        &result.overlap_stats.per_sample_ql_at_q_sum,
264        &result.overlap_stats.per_sample_ql2_at_q_sum,
265    ] {
266        for temperatures in samples {
267            for bins in temperatures {
268                for value in bins {
269                    value.to_bits().hash(&mut state_hash);
270                }
271            }
272        }
273    }
274
275    let per_sweep = elapsed / n_sweeps as f64 * 1000.0;
276    println!("Total: {:.3} s  |  {:.3} ms/sweep", elapsed, per_sweep);
277    println!("State checksum: {:016x}", state_hash.finish());
278}
bench/bench.rs

bench/
bench.rs
