//! The tracker's operator-Loschmidt-echo instances by sparse Pauli dynamics, //! against the device values and the classical values published for them. //! //! `cargo run --release --example tracker_ole --features quantum_engine,quantum_receipt,quantum_spd_receipt,quantum_bptn,quantum_mpsim -- [setting ...]` //! //! `instance` is `49x648`, `49x1296` or `70x1872`. Each circuit `W` is written //! `U · V · U†` on a 156-qubit register, `U` the Trotterised Ising evolution of //! the instance and `V` the perturbation; its echo is //! `f = 2⁻ⁿ Σ_z ⟨z|O|z⟩ ⟨z|W† O W|z⟩` for `O = Z₅₂ Z₅₉ Z₇₂`, which the device //! estimates from random basis states. Two methods, each swept: //! //! - `spd [threshold ...]`: sparse Pauli dynamics ([`echo`]) carries `O` //! through the last third only and takes the overlap after the //! perturbation, an exact identity, so each value comes with a bound on what //! truncation moved. //! - `bptn [χ ...]` and `bptn-gauged [χ ...]`: the echo as one state on twice //! the qubits (each qubit's ancilla in a Bell pair with it, `W` on the //! qubits, the ancillas' `Z`s folded onto theirs), evolved as a tensor //! network on the coupling graph by simple update (`bptn-gauged` re-gauges //! into the belief-propagation gauge after every layer), and read as `Z₅₂ Z₅₉ //! Z₇₂` with the 12-qubit loop that holds the three contracted exactly. //! //! - `mpo [χ ...]`: the circuit's operator as a matrix-product state on row //! and column sites, grown outward from the perturbation's centre, so that //! it starts as `V` (a product, near the identity) and grows only by `V`'s //! dressing; the echo is `⟨⟨W|(O ⊗ O)|W⟩⟩` on it. //! //! - `fast [χ ...]`: `bptn` with the dense blocked linear algebra //! (`Tns::with_blocked_linalg`), for bonds in the hundreds; read around the //! loop (while its bonds are small enough), and along the loop's arc through //! the three sites with belief propagation's environment and with the bond //! weights'. `fast-zero` is the same circuit without its middle, whose exact //! echo is 1: the divisor of the published "rescaled" values. //! //! Settings for these two, from the environment: `OLE_THREADS` (workers), //! `OLE_BP_SWEEPS` and `OLE_BP_TOL` (the readout's belief propagation; 0 //! sweeps reads with the bond weights alone), `OLE_LOOP_MAX` (the largest //! bond the loop readout is run at). //! //! A value is converged only where its sweep stops moving it. //! //! Each run writes a signed receipt to `RECEIPTS_DIR` (default `receipts`). #[path = "tracker/common.rs"] mod common; #[path = "tracker/echo.rs"] mod echo_net; #[path = "tracker/operator.rs"] mod operator; use serde_json::json; use wai_quantum::quantum_bptn::Tns; use wai_quantum::quantum_spd::{Pauli, SpdConfig, echo, echo_split}; /// Published values: (what produced it, the value, the tracker submission it is from). type Published = &'static [(&'static str, f64, u32)]; const INSTANCES: &[(&str, &str, Published)] = &[ ( "49x648", "data/observable-estimations/circuit-models/operator_loschmidt_echo/49Q_OLE_circuit_L_3_b_0.25_delta0.15.qasm", &[ ("device, global rescaling", 0.824, 38), ("belief-propagation tensor network, bond 192", 0.8202512915, 15), ("belief-propagation tensor network, bond 512", 0.821658489, 18), ("single-path Monte Carlo", 0.808, 63), ], ), ( "49x1296", "data/observable-estimations/circuit-models/operator_loschmidt_echo/49Q_OLE_circuit_L_6_b_0.25_delta0.15.qasm", &[ ("device, global rescaling", 0.662, 37), ("device, second layout", 0.657, 77), ("second device", 0.649, 78), ("belief-propagation tensor network, bond 512", 0.88157984, 19), ("the same, rescaled", 0.94257142, 20), ("single-path Monte Carlo", 0.619, 64), ], ), ( "70x1872", "data/observable-estimations/circuit-models/operator_loschmidt_echo/70Q_OLE_circuit_L_6_b_0.25_delta0.15.qasm", &[ ("device, global rescaling", 0.673, 36), ("belief-propagation tensor network, bond 512", 0.873827508, 21), ("the same, rescaled", 0.95247989, 22), ("single-path Monte Carlo", 0.614, 65), ], ), ]; /// The observable every instance measures, on the device's qubits. const OBSERVABLE: [u32; 3] = [52, 59, 72]; fn main() { let mut args = std::env::args().skip(1); let root = common::tracker_root(args.next()); let name = args.next().unwrap_or_else(|| "49x648".into()); let method = args.next().unwrap_or_else(|| "spd".into()); let settings: Vec = args.collect(); let (_, path, published) = INSTANCES.iter().find(|(n, ..)| *n == name).unwrap_or_else(|| panic!("no instance {name}")); let src = String::from_utf8(common::read_pinned(&root, path)).expect("UTF-8"); let (circuit, map) = common::compact(&common::rotations(&src)); let (tail, middle) = echo_split(&circuit).expect("the circuit is written U · V · U†"); let circuit_hash = common::hex( &wai_quantum::quantum_spd_receipt::SpdJob { circuit: circuit.clone(), observable: vec![], ones: vec![], config: SpdConfig::default() }.circuit_hash(), ); let threads = std::thread::available_parallelism().map_or(1, |t| t.get()); println!("instance operator_loschmidt_echo_{name}"); println!("qubits {} of the device's register", circuit.n); println!("gates {} ({} in the tail, {} in the middle)", circuit.ops.len(), tail.ops.len(), middle.ops.len()); println!("circuit hash {circuit_hash}"); for (label, value, issue) in published.iter() { println!("published {value:<12} {label} (submission {issue})"); } let request = |extra: serde_json::Value| { let mut r = json!({ "instance": format!("operator_loschmidt_echo_{name}"), "source": common::source(path), "observable": OBSERVABLE.iter().map(|q| format!("Z{q}")).collect::>(), "circuit_hash": circuit_hash, "threads": threads, "published": published.iter().map(|(label, value, issue)| json!({ "label": label, "value": value, "submission": issue })).collect::>(), }); for (k, v) in extra.as_object().expect("an object") { r[k] = v.clone(); } r }; let joules = |e: &serde_json::Value| e.get("joules_micro").map_or("unmetered".to_string(), |j| format!("{:.1} J", j.as_u64().unwrap_or(0) as f64 / 1e6)); match method.as_str() { "spd" => { let obs = common::z_string(&OBSERVABLE, &map); let thresholds: Vec = settings.iter().map(|a| a.parse().expect("a threshold")).collect(); for threshold in if thresholds.is_empty() { vec![1e-3, 3e-4, 1e-4] } else { thresholds } { let cfg = SpdConfig { threshold, max_weight: None }; let (r, secs, energy) = common::metered(|| echo(&circuit, &[(1.0, obs.clone())], &cfg, threads).expect("a valid circuit")); println!( "threshold {threshold:<8e} f = {:.6} bound {:.3e} A {} terms peak {} final {} dropped {:.3e} + {:.3e} {secs:.1} s {}", r.value, r.truncation_bound, r.tail_terms, r.peak_terms, r.final_terms, r.tail_dropped, r.middle_dropped, joules(&energy), ); let req = request(json!({ "method": "sparse Pauli dynamics, echo overlap (quantum_spd::echo)", "threshold": threshold, "threshold_bits": common::bits(threshold), })); let result = json!({ "value": r.value, "value_bits": common::bits(r.value), "truncation_bound": r.truncation_bound, "tail_terms": r.tail_terms, "peak_terms": r.peak_terms, "final_terms": r.final_terms, "term_updates": r.term_updates, "tail_dropped": r.tail_dropped, "middle_dropped": r.middle_dropped, "wall_seconds": secs, "energy_figure": energy, }); let out = common::receipt_path(&format!("ole_{name}_spd_t{threshold:e}.json")); assert!(common::verifies(&common::seal(req, result, Some(&out)))); println!("receipt {}", out.display()); } } "bptn" | "bptn-gauged" => { let sites: Vec = common::z_string(&OBSERVABLE, &map).iter().map(|(q, _)| *q).collect(); let (purified, edges) = common::purified_echo(&circuit, &sites); let ring = common::loop_through(&common::coupling(&circuit), &sites).expect("the observable's qubits lie on one loop"); println!("loop {:?} (device qubits)", ring.iter().map(|&q| map[q as usize]).collect::>()); let chis: Vec = settings.iter().map(|a| a.parse().expect("a bond dimension")).collect(); for chi in if chis.is_empty() { vec![8, 16, 32] } else { chis } { let ((f, discarded, bond, sweeps, residual), secs, energy) = common::metered(|| { let mut t = Tns::new(purified.n, &edges, chi).expect("the graph"); if method == "bptn" { t.run(&purified).expect("gates on the graph"); } else { t.run_gauged(&purified, 50, 1e-10).expect("gates on the graph"); } let env = t.belief_propagation(500, 1e-12); let ops: Vec<(u32, Pauli)> = sites.iter().map(|&q| (q, Pauli::Z)).collect(); let f = t.expect_loop_bp(&env, &ring, &ops).expect("a loop"); (f, t.discarded(), t.max_bond(), env.sweeps, env.residual) }); println!( "chi {chi:<5} f = {f:.6} discarded {discarded:.3e} bond {bond} belief propagation {sweeps} sweeps, residual {residual:.1e} {secs:.1} s {}", joules(&energy), ); let req = request(json!({ "method": if method == "bptn" { "tensor network on the coupling graph by simple update, the echo purified (quantum_bptn::Tns::run), read with the observable's loop contracted exactly (expect_loop_bp)" } else { "tensor network on the coupling graph by simple update re-gauged after every layer, the echo purified (quantum_bptn::Tns::run_gauged), read with the observable's loop contracted exactly (expect_loop_bp)" }, "chi": chi, "loop": ring.iter().map(|&q| map[q as usize]).collect::>(), })); let result = json!({ "value": f, "value_bits": common::bits(f), "discarded": discarded, "max_bond": bond, "bp_sweeps": sweeps, "bp_residual": residual, "wall_seconds": secs, "energy_figure": energy, }); let out = common::receipt_path(&format!("ole_{name}_{method}_chi{chi}.json")); assert!(common::verifies(&common::seal(req, result, Some(&out)))); println!("receipt {}", out.display()); } } "mpo" => { let sites: Vec = common::z_string(&OBSERVABLE, &map).iter().map(|(q, _)| *q).collect(); let edges = common::coupling(&circuit); // V's centre: the middle of the gates the split leaves between U and U†. let mid = tail.ops.len() + middle.ops.len() / 2; let chis: Vec = settings.iter().map(|a| a.parse().expect("a bond dimension")).collect(); for chi in if chis.is_empty() { vec![16, 32, 64] } else { chis } { let ((f, report), secs, energy) = common::metered(|| operator::echo_by_operator(&circuit, &edges, chi, mid, &sites)); println!( "chi {chi:<5} f = {f:.6} discarded {:.3e} fidelity estimate {:.4} bond {} updates {} {secs:.1} s {}", report.discarded, report.fidelity_estimate, report.max_bond_used, report.updates, joules(&energy), ); let req = request(json!({ "method": "the circuit's operator, a matrix-product state on row and column sites grown outward from the perturbation's centre (quantum_mpsim::Chain), read as ⟨⟨W|(O ⊗ O)|W⟩⟩", "chi": chi, "junction": mid, })); let result = json!({ "value": f, "value_bits": common::bits(f), "discarded": report.discarded, "fidelity_estimate": report.fidelity_estimate, "max_bond": report.max_bond_used, "updates": report.updates, "truncations": report.truncations, "wall_seconds": secs, "energy_figure": energy, }); let out = common::receipt_path(&format!("ole_{name}_mpo_chi{chi}.json")); assert!(common::verifies(&common::seal(req, result, Some(&out)))); println!("receipt {}", out.display()); } } "fast" | "fast-zero" => { let sites: Vec = common::z_string(&OBSERVABLE, &map).iter().map(|(q, _)| *q).collect(); let head_len = circuit.ops.len() - tail.ops.len() - middle.ops.len(); let parts = echo_net::purified(&circuit, head_len, middle.ops.len(), &sites); let env_usize = |k: &str, d: usize| std::env::var(k).ok().and_then(|v| v.parse().ok()).unwrap_or(d); let workers = env_usize("OLE_THREADS", threads); let bp_sweeps = env_usize("OLE_BP_SWEEPS", 30); let bp_tol: f64 = std::env::var("OLE_BP_TOL").ok().and_then(|v| v.parse().ok()).unwrap_or(1e-9); let loop_max = env_usize("OLE_LOOP_MAX", 64); println!("path {:?} (device qubits)", parts.path.iter().map(|&q| map[q as usize]).collect::>()); let chis: Vec = settings.iter().map(|a| a.parse().expect("a bond dimension")).collect(); for chi in if chis.is_empty() { vec![16, 32, 64] } else { chis } { let (phases_before, _) = wai_quantum::quantum_svd::phases(); let steps_before = wai_quantum::quantum_bptn::blocked_phases(); let ((r, evolve_seconds), secs, energy) = common::metered(|| { let (net, evolve) = echo_net::whole(&parts, chi, workers, method == "fast"); (echo_net::read(&net, &parts, bp_sweeps, bp_tol, loop_max), evolve) }); let peak = common::peak_rss_bytes(); // Where the dense routines spent their time (summed over calls). let (phases_after, _) = wai_quantum::quantum_svd::phases(); let phases: Vec<(String, f64)> = phases_after.iter().zip(&phases_before).map(|((n, a), (_, b))| (n.to_string(), a - b)).collect(); // The headline: the loop closed exactly where it was run, else the arc. let value = r.bp_loop.or(r.bp_path).unwrap_or(r.vidal_path); let show = |x: Option| x.map_or("-".to_string(), |x| format!("{x:.6}")); println!( "chi {chi:<5} f = {value:.6} (BP loop {}, BP arc {}, simple-update arc {:.6}) discarded {:.3e} bond {} BP {} sweeps {:.1e} {secs:.1} s peak {} MB {}", show(r.bp_loop), show(r.bp_path), r.vidal_path, r.discarded, r.max_bond, r.bp_sweeps, r.bp_residual, peak.map_or(0, |b| b >> 20), joules(&energy), ); println!(" evolution {evolve_seconds:.1} s, readout {:.1} s; loop bonds {:?}", r.seconds, r.loop_bonds); let steps: Vec<(String, f64)> = wai_quantum::quantum_bptn::blocked_phases().iter().zip(&steps_before).map(|((n, a), (_, b))| (n.to_string(), a - b)).collect(); println!(" steps (summed over threads): {}", steps.iter().map(|(n, t)| format!("{n} {t:.1} s")).collect::>().join(", ")); println!(" inside the SVDs: {}", phases.iter().map(|(n, t)| format!("{n} {t:.1} s")).collect::>().join(", ")); let method_text = if method == "fast" { "the purified echo |W⟩⟩ on the coupling graph with an ancilla on every qubit, evolved in the circuit's order by simple update with the dense blocked linear algebra (Tns::with_blocked_linalg); read around the observable's loop with belief propagation's environment off it (expect_loop_bp), along its arc (expect_path_bp) with that environment and with the bond weights'" } else { "as fast, with the middle left out, so that the exact echo is 1: the divisor of a rescaled value" }; let req = request(json!({ "method": method_text, "chi": chi, "workers": workers, "bp_sweeps_max": bp_sweeps, "bp_tolerance": bp_tol, "loop_readout_max_bond": loop_max, "arc": parts.path.iter().map(|&q| map[q as usize]).collect::>(), "loop": parts.ring.iter().map(|&q| map[q as usize]).collect::>(), })); let opt = |x: Option| x.map_or(json!(null), |x| json!(x)); let opt_bits = |x: Option| x.map_or(json!(null), |x| json!(common::bits(x))); let result = json!({ "value": value, "value_bits": common::bits(value), "bp_loop": opt(r.bp_loop), "bp_loop_bits": opt_bits(r.bp_loop), "bp_arc": opt(r.bp_path), "bp_arc_bits": opt_bits(r.bp_path), "simple_update_arc": r.vidal_path, "simple_update_arc_bits": common::bits(r.vidal_path), "bp_sweeps": r.bp_sweeps, "bp_residual": if r.bp_residual.is_finite() { json!(r.bp_residual) } else { json!(null) }, "discarded": r.discarded, "max_bond": r.max_bond, "loop_bonds": r.loop_bonds, "evolve_seconds": evolve_seconds, "readout_seconds": r.seconds, "wall_seconds": secs, "peak_rss_bytes": peak, "energy_figure": energy, "step_seconds": steps.iter().map(|(n, t)| json!({ "step": n, "seconds": t })).collect::>(), "svd_phase_seconds": phases.iter().map(|(n, t)| json!({ "phase": n, "seconds": t })).collect::>(), }); let out = common::receipt_path(&format!("ole_{name}_{method}_chi{chi}.json")); assert!(common::verifies(&common::seal(req, result, Some(&out)))); println!("receipt {}", out.display()); } } "loops" | "loops-zero" => { let sites: Vec = common::z_string(&OBSERVABLE, &map).iter().map(|(q, _)| *q).collect(); let head_len = circuit.ops.len() - tail.ops.len() - middle.ops.len(); let parts = echo_net::purified(&circuit, head_len, middle.ops.len(), &sites); let env_usize = |k: &str, d: usize| std::env::var(k).ok().and_then(|v| v.parse().ok()).unwrap_or(d); let workers = env_usize("OLE_THREADS", threads); let bp_sweeps = env_usize("OLE_BP_SWEEPS", 30); let bp_tol: f64 = std::env::var("OLE_BP_TOL").ok().and_then(|v| v.parse().ok()).unwrap_or(1e-9); let max_cells = env_usize("OLE_CLUSTER_CELLS", 3); let max_work: f64 = std::env::var("OLE_CLUSTER_WORK").ok().and_then(|v| v.parse().ok()).unwrap_or(1e13); let max_bytes: f64 = std::env::var("OLE_CLUSTER_BYTES").ok().and_then(|v| v.parse().ok()).unwrap_or(32e9); let ring_len = 12; let dev = |v: u32| map[v as usize]; let chis: Vec = settings.iter().map(|a| a.parse().expect("a bond dimension")).collect(); for chi in if chis.is_empty() { vec![64] } else { chis } { let (((x, bp_arc, bp_loop, env_sweeps, env_residual, bonds, readout), evolve_seconds, net_discarded, net_bond), secs, energy) = common::metered(|| { let (net, evolve) = echo_net::whole(&parts, chi, workers, method == "loops"); let started = std::time::Instant::now(); let env = net.belief_propagation(bp_sweeps, bp_tol); let ops: Vec<(u32, Pauli)> = sites.iter().map(|&q| (q, Pauli::Z)).collect(); let bp_arc = net.expect_path_bp(&env, &parts.path, &ops).expect("the arc"); let bp_loop = net.expect_loop_bp(&env, &parts.ring, &ops).expect("the loop"); println!("chi {chi:<5} arc readout {bp_arc:.6} loop readout {bp_loop:.6} discarded {:.3e} bond {} BP {} sweeps {:.1e} evolution {evolve:.1} s", net.discarded(), net.max_bond(), env.sweeps, env.residual); let lam = net.bond_weights(); let bonds: Vec<(u32, u32, usize)> = parts.edges.iter().enumerate().filter(|(_, (a, b))| *a < parts.n && *b < parts.n).map(|(e, &(a, b))| (dev(a), dev(b), lam[e].len())).collect(); println!(" bonds {:?}", bonds); let x = echo_net::clusters(&net, &env, &parts, ring_len, max_cells, max_work, max_bytes); let readout = started.elapsed().as_secs_f64(); ((x, bp_arc, bp_loop, env.sweeps, env.residual, bonds, readout), evolve, net.discarded(), net.max_bond()) }); let echo_net::Expansion { cells, orders, reads, weights } = x; let peak = common::peak_rss_bytes(); let devs = |r: &[u32]| r.iter().map(|&v| dev(v)).collect::>(); let best = orders.iter().rev().find_map(|o| o.product); println!(" {secs:.1} s peak {} MB {}", peak.map_or(0, |b| b >> 20), joules(&energy)); let opt = |x: Option| x.map_or(json!(null), |x| json!(x)); let opt_bits = |x: Option| x.map_or(json!(null), |x| json!(common::bits(x))); let req = request(json!({ "method": if method == "loops" { "the purified echo |W⟩⟩ evolved as in fast; read by the loop cluster expansion around belief propagation (Gray, Park, Evenbly, Pancotti, Kjønstad and Chan 2025, Eq. 7 and Algorithm 1): the regions are the observable's ring joined with connected sets of the graph's other rings, closed under intersection with counting numbers c(r) = 1 − Σ_{a ⊋ r} c(a), each contracted exactly (Tns::expect_region_bp: its loops cut where the estimated work is least) with belief propagation's messages on its boundary; also the arc through the observable's sites (expect_path_bp, the tree below order 1) and each ring's own correction (Midha and Zhang 2025, Eq. 9)" } else { "as loops, with the middle left out, so that the exact echo is 1" }, "chi": chi, "workers": workers, "bp_sweeps_max": bp_sweeps, "bp_tolerance": bp_tol, "ring_length": ring_len, "max_cells": max_cells, "max_work": max_work, "max_table_bytes": max_bytes, "arc": devs(&parts.path), "loop": devs(&parts.ring), })); let result = json!({ "value": opt(best), "value_bits": opt_bits(best), "bp_arc": bp_arc, "bp_arc_bits": common::bits(bp_arc), "bp_loop": bp_loop, "bp_loop_bits": common::bits(bp_loop), "orders": orders.iter().map(|o| json!({ "cells": o.cells, "terms": o.terms.iter().map(|(r, c)| json!({ "region": devs(r), "count": c })).collect::>(), "product": opt(o.product), "product_bits": opt_bits(o.product), "sum": opt(o.sum), "sum_bits": opt_bits(o.sum), "partial": o.partial.as_ref().map_or(json!(null), |q| json!({ "left_out": q.left_out.iter().map(|r| devs(r)).collect::>(), "terms": q.terms.iter().map(|(r, c)| json!({ "region": devs(r), "count": c })).collect::>(), "product": opt(q.product), "product_bits": opt_bits(q.product), "sum": opt(q.sum), "sum_bits": opt_bits(q.sum), })), })).collect::>(), "regions": reads.iter().map(|r| json!({ "region": devs(&r.region), "cells": r.cells, "cuts": r.cuts.iter().map(|&(a, b, d)| json!([dev(a), dev(b), d])).collect::>(), "combinations": r.combinations, "work": r.work, "table_bytes": r.bytes, "value": opt(r.value), "value_bits": opt_bits(r.value), "norm_ratio": opt(r.norm_ratio), "seconds": r.seconds, })).collect::>(), "rings": cells.iter().zip(&weights).map(|(c, w)| json!({ "ring": devs(c), "norm_ratio": opt(*w), "norm_ratio_bits": opt_bits(*w) })).collect::>(), "bonds": bonds.iter().map(|&(a, b, d)| json!([a, b, d])).collect::>(), "bp_sweeps": env_sweeps, "bp_residual": env_residual, "discarded": net_discarded, "max_bond": net_bond, "evolve_seconds": evolve_seconds, "readout_seconds": readout, "wall_seconds": secs, "peak_rss_bytes": peak, "energy_figure": energy, }); let out = common::receipt_path(&format!("ole_{name}_{method}_chi{chi}.json")); assert!(common::verifies(&common::seal(req, result, Some(&out)))); println!("receipt {}", out.display()); } } other => panic!("spd, bptn, bptn-gauged, mpo, fast, fast-zero, loops or loops-zero, not {other}"), } }