710cb954a6
- `wai-quantum estimate prog.qasm [--budget E] [--qubit NAME|all|FILE.json] [--calibration exact_matching|union_find|assumed] [--max-slowdown N]`: the program's counts, the ε budget and its synthesis share per rotation, the calibration, the headline and the frontier, and the receipts checked and re-run here. `--json` prints the answer sealed with the user's key; `--counts FILE` and `--routine FILE` estimate a document; `--export counts|routine` prints the program's counts as one. `lab estimate` is unchanged and its help points here. - `wai-quantum verify` on an estimate checks its seal, every resource receipt and the request re-run, and says which held. - MCP: estimate_resources (sealed with the user's key) beside compile_qasm; verify_quantum_result re-runs an estimate. The tool's description and schema are the engine's, so both servers serve the same. - smoke.sh: an estimate's counts, report, seal, verify, the counts exported and read back, a routine on every reference qubit, lab estimate as before. CI compares an estimate natively and as a WASI component byte for byte. Co-Authored-By: Claude Opus 5.5 <[email protected]>
2015 lines
105 KiB
Rust
2015 lines
105 KiB
Rust
//! `wai-quantum lab`: the research tools of the command line, the deterministic
|
||
//! quantum stack's every capability.
|
||
//!
|
||
//! Every quantum capability in this crate is pure Rust with no host bindings, so
|
||
//! this binary targets **native and cloud WASM from one source**:
|
||
//! `cargo build --release --target wasm32-wasip2` produces a component that
|
||
//! `wasmtime` runs with **byte-identical output** — the same circuit reconstructs
|
||
//! to the same statevector hash on a laptop or at the edge. That equality is the
|
||
//! point: a quantum computation you can put in a receipt and re-run anywhere.
|
||
//! The one exception is `fabric`, which describes the host it runs on: its
|
||
//! certificates name the member they ran on, and agree everywhere in parity
|
||
//! and corpus digest.
|
||
//!
|
||
//! ```text
|
||
//! wai-quantum lab list
|
||
//! wai-quantum lab run ghz --qubits 5 [--top 8] [--json]
|
||
//! wai-quantum lab run --qasm # read OpenQASM on stdin
|
||
//! wai-quantum lab qasm bell # emit OpenQASM
|
||
//! wai-quantum lab receipt ghz -n 3 [--shots 4096 --seed 7]
|
||
//! wai-quantum lab verify # read a receipt JSON on stdin
|
||
//! ```
|
||
//!
|
||
//! A circuit is built in `wai.quantum.circuit2`, whose phase table is correctly
|
||
//! rounded. `--gateset v1` builds it in the frozen `wai.quantum.circuit` instead,
|
||
//! to reproduce or re-verify content made with it. Every command that builds a
|
||
//! circuit from an algorithm, OpenQASM or QIR (`run`, `receipt`, `qasm`, `qir`,
|
||
//! `entropy`) takes it, and so does `fabric`.
|
||
//!
|
||
//! Beyond running circuits, the same binary answers what a state *carries* and
|
||
//! what transport does to it:
|
||
//!
|
||
//! ```text
|
||
//! wai-quantum lab entropy bell --keep 0 # half a Bell pair: 1 bit, and no state vector
|
||
//! wai-quantum lab channel depolarizing --prob 0.5
|
||
//! wai-quantum lab bb84 --eavesdrop # interception shows up in the error rate
|
||
//! wai-quantum lab chsh # 2√2, exactly, not sampled
|
||
//! wai-quantum lab qutrit # d = 3, not just qubits
|
||
//! ```
|
||
//!
|
||
//! At utility scale, sparse Pauli dynamics computes observables of circuits far
|
||
//! past any state vector, with arbitrary rotation angles:
|
||
//!
|
||
//! ```text
|
||
//! wai-quantum lab utility --theta 0.2 --steps 20 # the 127-qubit heavy-hex experiment
|
||
//! wai-quantum lab utility --method bptn --chi 32 # the same, by tensor networks
|
||
//! wai-quantum lab utility --receipt | wai-quantum verify
|
||
//! ```
|
||
//!
|
||
//! And for fault-tolerant machines that do not exist yet, what a program would
|
||
//! need — physical qubits and time — on six reference qubits, as a re-checkable
|
||
//! signed estimate:
|
||
//!
|
||
//! ```text
|
||
//! wai-quantum lab estimate factoring # 2048-bit factoring, all six qubits
|
||
//! wai-quantum lab estimate chemistry --qubit ns-e4 --frontier
|
||
//! wai-quantum lab estimate --qubits 100 --ccz 1e6 --budget 0.01
|
||
//! wai-quantum lab estimate factoring --qubit ns-e4 --receipt | wai-quantum verify
|
||
//! ```
|
||
//!
|
||
//! Error correction itself, at circuit level: sample a noisy Clifford circuit,
|
||
//! derive its detector error model, decode it exactly, and measure the code
|
||
//! model that estimates otherwise take on trust:
|
||
//!
|
||
//! ```text
|
||
//! wai-quantum lab qec -d 5 --p 0.003 --shots 100000
|
||
//! wai-quantum lab qec --file circuit.txt
|
||
//! wai-quantum lab qec --fit # then: estimate --prefactor A --threshold P
|
||
//! ```
|
||
//!
|
||
//! The `channel` command prints the Holevo bound before and after, which can only
|
||
//! fall: noise never creates information. Two results there are easy to miss —
|
||
//! fully damping `|1⟩` yields a *pure* output that still carries nothing (both
|
||
//! inputs land on `|0⟩`), and dephasing leaves a Z-basis ensemble's Holevo bound
|
||
//! untouched, because loss is relative to the encoding rather than absolute.
|
||
|
||
use wai_quantum::ed25519_dalek::SigningKey;
|
||
use std::io::Read;
|
||
use wai_quantum::quantum::{BaseGate, Circuit, Gateset, Measure, ONE};
|
||
use wai_quantum::quantum_channel::{tomography, Channel};
|
||
use wai_quantum::quantum_comm::{bb84, clone_fidelity};
|
||
use wai_quantum::quantum_info::{chsh, DensityMatrix};
|
||
use wai_quantum::quantum_qudit::QuditState;
|
||
use wai_quantum::quantum_energy::{EnergyClass, MaybeLabelled};
|
||
use wai_quantum::quantum_receipt::QuantumReceipt;
|
||
use wai_quantum::quantum_receipt::did_key;
|
||
use wai_quantum::quantum_resource::{self as qr, EstimateConfig, LogicalCounts, PhysicalQubit, QecModel};
|
||
use wai_quantum::quantum_resource_receipt::{Program, ResourceJob, ResourceReceipt};
|
||
use wai_quantum::quantum_spd_receipt::{SpdJob, SpdReceipt};
|
||
use wai_quantum::quantum_source::{holevo_bound, schumacher_limit, von_neumann_entropy};
|
||
use wai_quantum::quantum_toolchain as qt;
|
||
|
||
/// The three single-qubit states the channel command sends: the ground state, the
|
||
/// excited one, and the superposition whose coherence noise actually destroys.
|
||
fn one_qubit_state(which: &str) -> DensityMatrix {
|
||
let mut c = Circuit::new(1);
|
||
match which {
|
||
"zero" => {}
|
||
"one" => {
|
||
c.x(0);
|
||
}
|
||
"plus" => {
|
||
c.h(0);
|
||
}
|
||
other => die(&format!("unknown state '{other}' (zero, one, plus)")),
|
||
}
|
||
DensityMatrix::from_pure(&c.simulate().unwrap_or_else(|e| die(&format!("simulate: {e:?}"))))
|
||
}
|
||
|
||
/// The joule figure as a line shows it: `0` is the unmetered marker, shown as
|
||
/// unmetered and never as `0 uJ` (energy-measurement §2).
|
||
fn joules_line(joules_micro: u64) -> String {
|
||
if joules_micro == 0 {
|
||
"unmetered (no figure; 0 is the unmetered marker, not zero joules)".into()
|
||
} else {
|
||
format!("{joules_micro} uJ")
|
||
}
|
||
}
|
||
|
||
/// One line saying what the joule figure is: its class and evidence, or that
|
||
/// it has none.
|
||
fn energy_class_line(class: Option<&EnergyClass>, joules_micro: u64) -> String {
|
||
let Some(c) = class else {
|
||
return if joules_micro == 0 {
|
||
"none (unmetered)".into()
|
||
} else {
|
||
"unlabelled (legacy form: read as no class, not comparable with a labelled figure)".into()
|
||
};
|
||
};
|
||
let mut line = c.label().to_string();
|
||
if let Some(u) = c.uncertainty() {
|
||
line += &format!(
|
||
" (declared u {:.2} %, k=1, over {:.3} s",
|
||
u.relative_ppm as f64 / 1e4,
|
||
u.window_us as f64 / 1e6
|
||
);
|
||
match c {
|
||
EnergyClass::OnChipCounter { filtering, .. } => line += &format!(", filtering {})", filtering.label()),
|
||
EnergyClass::CalibratedInstrument { calibration_ref, .. } => line += &format!(", calibration {calibration_ref})"),
|
||
_ => line += ")",
|
||
}
|
||
}
|
||
line
|
||
}
|
||
|
||
fn hx(b: &[u8]) -> String {
|
||
b.iter().map(|x| format!("{x:02x}")).collect()
|
||
}
|
||
|
||
fn die(msg: &str) -> ! {
|
||
eprintln!("wai-quantum: {msg}");
|
||
std::process::exit(2);
|
||
}
|
||
|
||
fn usage() -> ! {
|
||
eprintln!(
|
||
"wai-quantum lab — the research tools (run `wai-quantum --help` for run, plan, lint, draw, verify and the rest)\n\n\
|
||
USAGE:\n \
|
||
wai-quantum lab list\n \
|
||
wai-quantum lab run <algorithm> [-n QUBITS] [-p PARAM] [--top K] [--json] [--gateset v1|v2]\n \
|
||
wai-quantum lab run --qasm [--top K] [--json] [--gateset v1|v2] (OpenQASM on stdin)\n \
|
||
wai-quantum lab run --qir [--top K] [--json] [--gateset v1|v2] (QIR on stdin)\n \
|
||
wai-quantum lab qasm <algorithm> [-n QUBITS] [-p PARAM] [--gateset v1|v2]\n \
|
||
wai-quantum lab qir <algorithm> [-n QUBITS] [-p PARAM] [--no-measure] [--gateset v1|v2]\n \
|
||
wai-quantum lab receipt <algorithm> [-n QUBITS] [-p PARAM] [--shots N] [--seed S] [--gateset v1|v2] [--key FILE | --demo-key]\n \
|
||
wai-quantum lab verify (receipt JSON on stdin)\n \
|
||
--gateset v2 (the default) is wai.quantum.circuit2; v1 is the frozen wai.quantum.circuit.\n\n\
|
||
INFORMATION AND TRANSPORT:\n \
|
||
wai-quantum lab entropy <algorithm> [-n QUBITS] [--keep 0,1] [--json] [--gateset v1|v2]\n \
|
||
wai-quantum lab channel <name> [--prob P] [--state zero|one|plus] [--json]\n \
|
||
wai-quantum lab bb84 [--rounds N] [--seed S] [--eavesdrop] [--json]\n \
|
||
wai-quantum lab chsh [--json]\n \
|
||
wai-quantum lab qutrit [-n QUTRITS] [--json]\n\n\
|
||
UTILITY SCALE:\n \
|
||
wai-quantum lab utility [--theta T] [--steps N] [--threshold D] [-q QUBIT] [--json] [--receipt]\n \
|
||
wai-quantum lab utility --method bptn [--chi X] [--theta T] [--steps N] [-q QUBIT] [--json]\n\n\
|
||
FAULT-TOLERANT RESOURCES:\n \
|
||
wai-quantum lab estimate [factoring|chemistry] [--qubit NAME|all] [--code surface|surface-meas|floquet]\n \
|
||
[--budget E] [--slowdown N] [--frontier] [--json] [--receipt]\n \
|
||
wai-quantum lab estimate --qubits N [--ccz N] [--t N] [--rotations N] [--depth N] [--measurements N] ...\n \
|
||
Qubits: us-e3, us-e4, ns-e3, ns-e4, meas-e4, meas-e6. A measured code: --prefactor A --threshold P.\n \
|
||
From a program's own counts, with the code model fitted to the decoders: wai-quantum estimate prog.qasm\n\n\
|
||
ERROR CORRECTION AT CIRCUIT LEVEL:\n \
|
||
wai-quantum lab qec [-d D] [--rounds R] [--p P] [--shots N] [--seed S] [--decoder exact|uf] [--json]\n \
|
||
wai-quantum lab qec --file CIRCUIT [--shots N] (the ecosystem's circuit text)\n \
|
||
wai-quantum lab qec [--circuit | --dem] [-d D | --file CIRCUIT]\n \
|
||
wai-quantum lab qec --fit [--shots N] (measure a and p* of the code model)\n\n\
|
||
GAUSSIAN BOSON SAMPLING:\n \
|
||
wai-quantum lab gbs --squeezing r1,r2,… [--eta E] [--seed S | --unitary U.json] [--json]\n \
|
||
wai-quantum lab gbs … --sample K (exact click patterns)\n \
|
||
wai-quantum lab gbs … --check SAMPLES.txt (observed patterns vs exact marginals)\n\n\
|
||
SPIN CHAINS (DMRG ground states, exact referees):\n \
|
||
wai-quantum lab dmrg [--chain ising|heisenberg] [-n N] [--j J] [--jz JZ] [--h H] [--bond D] [--sweeps S]\n \
|
||
[--cutoff C] [--exact] [--json]\n \
|
||
wai-quantum lab tdvp [--chain ising|heisenberg] [-n N] [--j J] [--jz JZ] [--h H] [--start plus|neel|up]\n \
|
||
[--bond D] [--dt DT] [--t T] [--every K] [--cutoff C] [--json] (time evolution by TDVP)\n\n\
|
||
MOLECULES (DMRG from FCIDUMP integrals, full CI as referee):\n \
|
||
wai-quantum lab chem --fcidump FILE [--bond D] [--sweeps S] [--penalty P] [--no-fci] [--json]\n\n\
|
||
NEURAL QUANTUM STATES (RBM, variational Monte Carlo, stochastic reconfiguration):\n \
|
||
wai-quantum lab nqs [--model heisenberg|j1j2|ising] [--lx L] [--ly L] [--j2 J] [--h H] [--periodic]\n \
|
||
[--complex] [--project [--kpi] [--point-group a1|a2|b1|b2 | --reflections a1|a2|b1|b2]\n \
|
||
[--odd]] [--alpha A] [--steps N] [--samples S] [--lr R]\n \
|
||
[--init W] [--guide ε] [--exact] [--seed S] [--json]\n\n\
|
||
THE FABRIC (members, parity against the CPU reference, placement):\n \
|
||
wai-quantum lab fabric [--circuits N] [--seed S] [--gateset v1|v2] [--json]\n\n\
|
||
OPTIMISATION (simulated bifurcation, annealing, exact):\n \
|
||
wai-quantum lab qio (--gset FILE | --sk N [--instance S]) [--method dsb|bsb|hdsb|hbsb|sa|exact]\n \
|
||
[--steps N] [--trials T] [--dt D] [--c1 C] [--gamma G] [--beta0 B --beta1 B] [--descend]\n \
|
||
[--target V] [--seed S] [--threads T] [--spins] [--json]\n\n\
|
||
QLDPC MEMORY (circuit level, BP-OSD):\n \
|
||
wai-quantum lab ldpc [--code 72|90|108|144|288] [--p P] [--cycles N] [--shots S] [--basis z|x|both]\n \
|
||
[--max-iter I] [--osd-order K] [--threads T] [--seed S] [--json]\n \
|
||
wai-quantum lab ldpc --circuit [--code …] [--p P] [--cycles N] [--basis z|x] (the circuit text)\n\n\
|
||
ARCHITECTURE:\n \
|
||
wai-quantum lab arch factoring [--bits 1024|1536|2048|3072|4096|6144|8192] [--json]\n \
|
||
wai-quantum lab arch storage [--target T] [--fits refit|uniform|si1000] [--json]\n \
|
||
wai-quantum lab arch cultivate [--p 5e-4|1e-3|2e-3] [--target T] [--json]\n\n\
|
||
NOISY CIRCUITS (trajectories, exact density matrix, cross-entropy):\n \
|
||
wai-quantum lab noisy --file CIRCUIT.jsonl [-n N] [--depol1 P] [--depol2 P] [--damping G] [--dephasing P]\n \
|
||
[--readout P] [--trajectories T] [--seed S] [--expect \"Z0 Z1; X3\"] [--xeb] [--exact] [--json]\n\n\
|
||
STATE VECTOR (any gate, any angle):\n \
|
||
wai-quantum lab sv --file CIRCUIT.jsonl [-n N] [--fuse W] [--top K | --sample S [--seed X] | --expect \"Z0 Z1; X3\"] [--json]\n\n\
|
||
TENSOR NETWORKS:\n \
|
||
wai-quantum lab tn --file CIRCUIT.jsonl [-n N] [--bits B] [--open q,q] [--trials T] [--seed S] [--cap-log2 L] [--whole-gates] [--plan-only] [--json]\n\n\
|
||
Channels: identity, bitflip, dephasing, depolarizing, damping.\n\
|
||
Run `wai-quantum lab list` for the algorithm catalog."
|
||
);
|
||
std::process::exit(2);
|
||
}
|
||
|
||
fn read_stdin() -> String {
|
||
let mut s = String::new();
|
||
if std::io::stdin().read_to_string(&mut s).is_err() {
|
||
die("could not read stdin");
|
||
}
|
||
s
|
||
}
|
||
|
||
/// The key a receipt is signed with ([`crate::keys`]), from `--key FILE` or
|
||
/// `--demo-key`.
|
||
fn signing_key(args: &[String]) -> SigningKey {
|
||
let given = opt(args, "--key", "").map(std::path::PathBuf::from);
|
||
let k = wai_quantum_cli::keys::load(given.as_deref(), flag(args, "--demo-key")).unwrap_or_else(|e| die(&e));
|
||
wai_quantum_cli::keys::notes(&k);
|
||
k.key
|
||
}
|
||
|
||
/// Pull `--flag VALUE` out of the argument list.
|
||
fn opt(args: &[String], long: &str, short: &str) -> Option<String> {
|
||
args.iter()
|
||
.position(|a| a == long || (!short.is_empty() && a == short))
|
||
.and_then(|i| args.get(i + 1).cloned())
|
||
}
|
||
fn flag(args: &[String], long: &str) -> bool {
|
||
args.iter().any(|a| a == long)
|
||
}
|
||
fn num<T: std::str::FromStr>(args: &[String], long: &str, short: &str, dflt: T) -> T {
|
||
match opt(args, long, short) {
|
||
Some(v) => v.parse().unwrap_or_else(|_| die("bad numeric argument")),
|
||
None => dflt,
|
||
}
|
||
}
|
||
|
||
/// What `--gateset v1` prints on standard error, never on standard output or in
|
||
/// JSON, for a circuit whose phases circuit/1 does not reproduce.
|
||
const V1_NOTE: &str =
|
||
"circuit/1 table: P(k≥8) deviates by up to 12868 units of 2^-30, P(k≥19) applies no phase; use --gateset v2";
|
||
|
||
/// The gateset `--gateset` names: `v2` (`wai.quantum.circuit2`, the default) or
|
||
/// `v1` (the frozen `wai.quantum.circuit`).
|
||
fn gateset_from(args: &[String]) -> Gateset {
|
||
match opt(args, "--gateset", "").as_deref() {
|
||
None | Some("v2") => Gateset::V2,
|
||
Some("v1") => Gateset::V1,
|
||
Some(other) => die(&format!("unknown gateset '{other}' (v1, v2)")),
|
||
}
|
||
}
|
||
|
||
/// Build the circuit a command operates on: either a named algorithm or OpenQASM
|
||
/// read from stdin (`--qasm`), in the gateset `--gateset` names. A circuit/1
|
||
/// circuit with a `P(k ≥ 8)` gets the note on standard error.
|
||
fn circuit_from(args: &[String], algo: Option<&String>) -> Circuit {
|
||
let gs = gateset_from(args);
|
||
let mut c = if flag(args, "--qasm") {
|
||
wai_quantum::qasm::dyadic_circuit(&read_stdin()).unwrap_or_else(|e| die(&format!("QASM parse: {e}")))
|
||
} else if flag(args, "--qir") {
|
||
wai_quantum::quantum_qir::from_qir(&read_stdin()).unwrap_or_else(|e| die(&format!("QIR parse: {e}")))
|
||
} else {
|
||
let id = algo.unwrap_or_else(|| usage());
|
||
let n: u8 = num(args, "--qubits", "-n", 3);
|
||
let p: u32 = num(args, "--param", "-p", 0);
|
||
qt::build_algorithm(id, n, p)
|
||
.unwrap_or_else(|| die(&format!("unknown algorithm '{id}' (try `wai-quantum list`)")))
|
||
};
|
||
// Every circuit2 circuit is a valid circuit/1 circuit, so the switch is safe.
|
||
c.gateset = gs;
|
||
if gs == Gateset::V1 && c.ops.iter().any(|g| g.base == BaseGate::P && g.param >= 8) {
|
||
eprintln!("{V1_NOTE}");
|
||
}
|
||
c
|
||
}
|
||
|
||
fn describe(c: &Circuit) -> (String, String, u64, String) {
|
||
let (backend, why) = qt::recommend_backend(c);
|
||
(
|
||
format!("{backend:?}"),
|
||
why.to_string(),
|
||
qt::statevector_work(c),
|
||
hx(&qt::circuit_id(c).unwrap_or_else(|e| die(&format!("circuit id: {e}")))),
|
||
)
|
||
}
|
||
|
||
/// The six reference qubits by command-line name.
|
||
fn reference_qubit(name: &str) -> PhysicalQubit {
|
||
match name {
|
||
"us-e3" => PhysicalQubit::us_e3(),
|
||
"us-e4" => PhysicalQubit::us_e4(),
|
||
"ns-e3" => PhysicalQubit::ns_e3(),
|
||
"ns-e4" => PhysicalQubit::ns_e4(),
|
||
"meas-e4" => PhysicalQubit::meas_ns_e4(),
|
||
"meas-e6" => PhysicalQubit::meas_ns_e6(),
|
||
other => die(&format!("unknown qubit '{other}' (us-e3, us-e4, ns-e3, ns-e4, meas-e4, meas-e6, all)")),
|
||
}
|
||
}
|
||
|
||
/// A duration in the unit a person would say it in.
|
||
pub(crate) fn human_time(ns: f64) -> String {
|
||
let s = ns / 1e9;
|
||
let (m, h, d, y) = (60.0, 3600.0, 86_400.0, 365.25 * 86_400.0);
|
||
if s < 1e-3 {
|
||
format!("{:.1} us", ns / 1e3)
|
||
} else if s < 1.0 {
|
||
format!("{:.1} ms", s * 1e3)
|
||
} else if s < m {
|
||
format!("{s:.2} s")
|
||
} else if s < h {
|
||
format!("{:.1} min", s / m)
|
||
} else if s < d {
|
||
format!("{} h {} min", (s / h) as u64, ((s % h) / m) as u64)
|
||
} else if s < y {
|
||
format!("{:.1} days", s / d)
|
||
} else {
|
||
format!("{:.1} years", s / y)
|
||
}
|
||
}
|
||
|
||
/// A circuit file: JSON lines, one gate per line, `{"q":[a,b],"m":[re,im, …]}`
|
||
/// (rows are outputs, the first qubit most significant). Returns the gates
|
||
/// and the register size they need.
|
||
fn read_gate_file(path: &str) -> (Vec<wai_quantum::gates::Gate>, u32) {
|
||
use wai_quantum::gates::{complex, Gate};
|
||
let text = std::fs::read_to_string(path).unwrap_or_else(|e| die(&format!("{path}: {e}")));
|
||
let mut gates = Vec::new();
|
||
let mut n: u32 = 0;
|
||
for (i, line) in text.lines().enumerate().filter(|(_, l)| !l.trim().is_empty()) {
|
||
let v: serde_json::Value = serde_json::from_str(line).unwrap_or_else(|_| die(&format!("line {}: not JSON", i + 1)));
|
||
let qubits: Vec<u32> = v["q"]
|
||
.as_array()
|
||
.unwrap_or_else(|| die(&format!("line {}: missing \"q\"", i + 1)))
|
||
.iter()
|
||
.map(|x| x.as_u64().unwrap_or_else(|| die("qubits are integers")) as u32)
|
||
.collect();
|
||
let m: Vec<f64> = v["m"]
|
||
.as_array()
|
||
.unwrap_or_else(|| die(&format!("line {}: missing \"m\"", i + 1)))
|
||
.iter()
|
||
.map(|x| x.as_f64().unwrap_or_else(|| die("matrix entries are numbers")))
|
||
.collect();
|
||
n = n.max(qubits.iter().max().map_or(0, |q| q + 1));
|
||
gates.push(Gate { qubits, matrix: m.chunks(2).map(|c| complex(c[0], *c.get(1).unwrap_or(&0.0))).collect() });
|
||
}
|
||
(gates, n)
|
||
}
|
||
|
||
/// A program for `estimate`: a named example, or counts from flags.
|
||
/// Returns the counts and the error budget the example is stated with.
|
||
fn estimate_program(args: &[String], name: Option<&String>) -> (LogicalCounts, f64) {
|
||
let counted = ["--qubits", "--ccz", "--t", "--rotations", "--depth", "--measurements"]
|
||
.iter()
|
||
.any(|f| opt(args, f, "").is_some());
|
||
if counted {
|
||
// Counts may be written as 1e9; they must still be whole numbers.
|
||
let count = |f: &str| -> u64 {
|
||
let v: f64 = num(args, f, "", 0.0);
|
||
if !(v >= 0.0 && v.fract() == 0.0 && v < 1.8e19) {
|
||
die(&format!("{f} must be a whole number"));
|
||
}
|
||
v as u64
|
||
};
|
||
let c = LogicalCounts {
|
||
qubits: count("--qubits"),
|
||
measurements: count("--measurements"),
|
||
rotations: count("--rotations"),
|
||
rotation_depth: count("--depth"),
|
||
t_gates: count("--t"),
|
||
ccz: count("--ccz"),
|
||
};
|
||
if c.rotations > 0 && c.rotation_depth == 0 {
|
||
die("--rotations needs --depth (non-Clifford layers containing a rotation)");
|
||
}
|
||
return (c, 0.001);
|
||
}
|
||
match name.map(String::as_str).unwrap_or("factoring") {
|
||
// 2048-bit factoring (windowed, coset representation), ε = 1/3: the
|
||
// answer is checkable, so a failed run is simply repeated.
|
||
"factoring" => (
|
||
LogicalCounts { qubits: 12_581, measurements: 1_080_000_000, rotations: 12, rotation_depth: 12, t_gates: 12, ccz: 3_730_000_000 },
|
||
1.0 / 3.0,
|
||
),
|
||
// Ground-state energy of a 76-orbital catalyst model by double-
|
||
// factorized qubitization, ε = 0.01.
|
||
"chemistry" => (
|
||
LogicalCounts {
|
||
qubits: 1318,
|
||
measurements: 1_370_000_000,
|
||
rotations: 206_000_000,
|
||
rotation_depth: 205_000_000,
|
||
t_gates: 55_300_000,
|
||
ccz: 135_000_000_000,
|
||
},
|
||
0.01,
|
||
),
|
||
other => die(&format!("unknown program '{other}' (factoring, chemistry, or counts via --qubits ...)")),
|
||
}
|
||
}
|
||
|
||
fn estimate_json(e: &qr::Estimate) -> String {
|
||
let rounds: Vec<String> = e
|
||
.factory
|
||
.iter()
|
||
.flat_map(|f| f.rounds.iter())
|
||
.map(|r| {
|
||
format!(
|
||
"{{\"copies\":{},\"distance\":{},\"unit\":\"{}\"}}",
|
||
r.copies,
|
||
r.distance,
|
||
match r.kind {
|
||
qr::UnitKind::SpaceEfficient => "space-efficient",
|
||
qr::UnitKind::ReedMuller => "reed-muller",
|
||
}
|
||
)
|
||
})
|
||
.collect();
|
||
format!(
|
||
"{{\"code\":{},\"distance\":{},\"factories\":{},\"factory_fraction\":{},\"factory_rounds\":[{}],\
|
||
\"logical_steps\":{},\"physical_qubits\":{},\"qubit\":{},\"runtime_s\":{},\"t_states\":{},\"tiles\":{}}}",
|
||
serde_json::to_string(&e.code).unwrap(),
|
||
e.distance,
|
||
e.factories,
|
||
e.factory_fraction(),
|
||
rounds.join(","),
|
||
e.steps,
|
||
e.physical_qubits,
|
||
serde_json::to_string(&e.qubit).unwrap(),
|
||
e.runtime_ns / 1e9,
|
||
e.isa.t_states,
|
||
e.isa.tiles,
|
||
)
|
||
}
|
||
|
||
/// Check a receipt of any kind this tool or the library signs: `src` is its
|
||
/// JSON. Prints what it binds, and exits 1 when it does not verify.
|
||
pub fn verify_receipt(src: &str) {
|
||
// Either form: a receipt as sealed, or one whose joule figure
|
||
// carries its acquisition class inside the signature
|
||
// (`quantum_energy::Labelled`, as `wai_quantum_meter` emits).
|
||
if let Some(any) = MaybeLabelled::<ResourceReceipt>::from_json(src) {
|
||
let ok = any.verify();
|
||
let r = any.receipt();
|
||
let e = &r.estimate;
|
||
println!("signature {}", if ok { "VALID" } else { "INVALID" });
|
||
println!("kind fault-tolerant resource estimate");
|
||
println!("job-hash {}", hx(&r.job_hash));
|
||
println!("qubit {} ({}, d = {})", e.qubit, e.code, e.distance);
|
||
println!("physical {} qubits", e.physical_qubits);
|
||
println!("runtime {}", human_time(e.runtime_ns));
|
||
println!("factories {}", e.factories);
|
||
println!("joules {}", joules_line(r.joules_micro));
|
||
println!("energy class {}", energy_class_line(any.energy_class(), r.joules_micro));
|
||
println!("signer {}", r.signer_id);
|
||
if !ok {
|
||
std::process::exit(1);
|
||
}
|
||
return;
|
||
}
|
||
if let Some(any) = MaybeLabelled::<SpdReceipt>::from_json(src) {
|
||
let ok = any.verify();
|
||
let r = any.receipt();
|
||
println!("signature {}", if ok { "VALID" } else { "INVALID" });
|
||
println!("kind sparse Pauli dynamics");
|
||
println!("circuit-hash {}", hx(&r.circuit_hash));
|
||
println!("query-hash {}", hx(&r.query_hash));
|
||
println!("value {}", r.result.value);
|
||
println!("bound {}", r.result.truncation_bound);
|
||
println!("work {} term-updates (peak {} terms)", r.result.term_updates, r.result.peak_terms);
|
||
println!("joules {}", joules_line(r.joules_micro));
|
||
println!("energy class {}", energy_class_line(any.energy_class(), r.joules_micro));
|
||
println!("signer {}", r.signer_id);
|
||
if !ok {
|
||
std::process::exit(1);
|
||
}
|
||
return;
|
||
}
|
||
match MaybeLabelled::<QuantumReceipt>::from_json(src) {
|
||
None => {
|
||
println!("INVALID (could not parse receipt JSON, or its energy class is malformed)");
|
||
std::process::exit(1);
|
||
}
|
||
Some(any) => {
|
||
let ok = any.verify();
|
||
let r = any.receipt();
|
||
println!("signature {}", if ok { "VALID" } else { "INVALID" });
|
||
println!("circuit-hash {}", hx(&r.circuit_hash));
|
||
println!("state-hash {}", hx(&r.statevector_hash));
|
||
println!("qubits {}", r.n_qubits);
|
||
println!("gates {}", r.n_ops);
|
||
println!("work {} amplitude-updates", r.work_amp_updates);
|
||
println!("joules {}", joules_line(r.joules_micro));
|
||
println!("energy class {}", energy_class_line(any.energy_class(), r.joules_micro));
|
||
println!("signer {}", r.signer_id);
|
||
if !ok {
|
||
std::process::exit(1);
|
||
}
|
||
}
|
||
}
|
||
}
|
||
|
||
/// The research tools, with `args` the words after `wai-quantum lab`.
|
||
pub fn main(args: Vec<String>) {
|
||
let argv: Vec<String> = std::iter::once("wai-quantum lab".to_string()).chain(args).collect();
|
||
if argv.len() < 2 {
|
||
usage();
|
||
}
|
||
let cmd = argv[1].as_str();
|
||
let args = &argv[2..];
|
||
let algo = args.iter().find(|a| !a.starts_with('-'));
|
||
|
||
match cmd {
|
||
"list" => {
|
||
println!("{:<14} {:<26} DESCRIPTION", "ID", "NAME");
|
||
for (id, name, desc) in qt::catalog() {
|
||
println!("{id:<14} {name:<26} {desc}");
|
||
}
|
||
}
|
||
|
||
"qasm" => {
|
||
let c = circuit_from(args, algo);
|
||
print!("{}", qt::to_qasm(&c));
|
||
}
|
||
|
||
"qir" => {
|
||
let c = circuit_from(args, algo);
|
||
print!("{}", wai_quantum::quantum_qir::to_qir(&c, !flag(args, "--no-measure")));
|
||
}
|
||
|
||
"run" => {
|
||
let c = circuit_from(args, algo);
|
||
if let Err(e) = c.validate() {
|
||
die(&format!("invalid circuit: {e:?}"));
|
||
}
|
||
let sv = c.simulate().unwrap_or_else(|e| die(&format!("simulate: {e:?}")));
|
||
let top: usize = num(args, "--top", "", 8);
|
||
let (backend, why, work, cid) = describe(&c);
|
||
let hash = hx(&sv.statevector_hash());
|
||
let outs = qt::top_outcomes(&sv, top);
|
||
let n = c.n_qubits as usize;
|
||
|
||
if flag(args, "--json") {
|
||
let items: Vec<String> = outs
|
||
.iter()
|
||
.map(|(k, p)| format!("{{\"state\":\"{k:0n$b}\",\"p\":{p}}}", n = n))
|
||
.collect();
|
||
println!(
|
||
"{{\"circuit_id\":\"{cid}\",\"capability\":\"{}\",\"n_qubits\":{},\"n_gates\":{},\"backend\":\"{backend}\",\
|
||
\"work_amp_updates\":{work},\"statevector_hash\":\"{hash}\",\"top\":[{}]}}",
|
||
c.gateset.capability(),
|
||
c.n_qubits,
|
||
c.ops.len(),
|
||
items.join(",")
|
||
);
|
||
} else {
|
||
println!("circuit-id {cid}");
|
||
println!("capability {}", c.gateset.capability());
|
||
println!("qubits {}", c.n_qubits);
|
||
println!("gates {}", c.ops.len());
|
||
println!("backend {backend} ({why})");
|
||
println!("work {work} amplitude-updates");
|
||
println!("state-hash {hash}");
|
||
println!("top outcomes");
|
||
for (k, p) in outs {
|
||
if p > 0.0 {
|
||
println!(" |{k:0n$b}> {p:.6}", n = n);
|
||
}
|
||
}
|
||
}
|
||
}
|
||
|
||
"receipt" => {
|
||
let c = circuit_from(args, algo);
|
||
let shots: u64 = num(args, "--shots", "", 0);
|
||
let seed: u64 = num(args, "--seed", "", 0);
|
||
let measure = (shots > 0).then_some(Measure { seed, shots });
|
||
let key = signing_key(args);
|
||
let r = QuantumReceipt::seal(&key, did_key(&key.verifying_key()), &c, measure, 0, None)
|
||
.unwrap_or_else(|e| die(&format!("seal: {e:?}")));
|
||
println!("{}", r.to_json());
|
||
}
|
||
|
||
"verify" => verify_receipt(&read_stdin()),
|
||
|
||
// How much quantum information a state actually carries. A state vector
|
||
// cannot answer this for a subsystem — half a Bell pair has no state
|
||
// vector — so `--keep` traces the rest out and reports the entropy that
|
||
// only a density matrix can see.
|
||
// The 127-qubit heavy-hex kicked-Ising experiment, by sparse Pauli
|
||
// dynamics: ⟨Z_q⟩ after `--steps` Trotter steps of RX(θ_h) on every
|
||
// qubit and RZZ(−π/2) on every edge, and the bound on how far
|
||
// truncation can have moved it. `--receipt` seals the run instead.
|
||
"utility" => {
|
||
use wai_quantum::quantum_spd::{heavy_hex_127, kicked_ising, Pauli, SpdConfig};
|
||
let theta: f64 = num(args, "--theta", "", std::f64::consts::FRAC_PI_4);
|
||
let steps: u32 = num(args, "--steps", "", 20);
|
||
let threshold: f64 = num(args, "--threshold", "", 1e-3);
|
||
let qubit: u32 = num(args, "--qubit", "-q", 62);
|
||
// The same circuit by belief-propagation tensor networks: an
|
||
// independent method, whose errors are unrelated to truncating a
|
||
// Pauli sum.
|
||
if opt(args, "--method", "").as_deref() == Some("bptn") {
|
||
if flag(args, "--receipt") {
|
||
die("--receipt seals sparse Pauli dynamics runs; tensor-network runs have no receipt yet");
|
||
}
|
||
let chi: usize = num(args, "--chi", "", 16);
|
||
let edges = heavy_hex_127();
|
||
let mut tns = wai_quantum::quantum_bptn::Tns::new(127, &edges, chi).unwrap_or_else(|e| die(&e.to_string()));
|
||
tns.run(&kicked_ising(127, &edges, theta, steps)).unwrap_or_else(|e| die(&e.to_string()));
|
||
let value = tns.expect1(qubit, Pauli::Z).unwrap_or_else(|e| die(&e.to_string()));
|
||
if flag(args, "--json") {
|
||
println!(
|
||
"{{\"chi\":{chi},\"discarded\":{},\"max_bond\":{},\"method\":\"bptn\",\"qubit\":{qubit},\
|
||
\"steps\":{steps},\"theta_h\":{theta},\"value\":{value}}}",
|
||
tns.discarded(),
|
||
tns.max_bond()
|
||
);
|
||
} else {
|
||
println!("<Z_{qubit}> {value:.6}");
|
||
println!("discarded {:.3e} (squared singular-value weight truncated)", tns.discarded());
|
||
println!("max bond {}", tns.max_bond());
|
||
}
|
||
return;
|
||
}
|
||
if opt(args, "--method", "").as_deref().is_some_and(|m| m != "spd") {
|
||
die("--method is spd (the default) or bptn");
|
||
}
|
||
let job = SpdJob {
|
||
circuit: kicked_ising(127, &heavy_hex_127(), theta, steps),
|
||
observable: vec![(1.0, vec![(qubit, Pauli::Z)])],
|
||
ones: vec![],
|
||
config: SpdConfig { threshold, max_weight: None },
|
||
};
|
||
if flag(args, "--receipt") {
|
||
let key = signing_key(args);
|
||
let r = SpdReceipt::seal(&key, did_key(&key.verifying_key()), &job, None, |run| {
|
||
run();
|
||
0
|
||
})
|
||
.unwrap_or_else(|e| die(&format!("seal: {e}")));
|
||
println!("{}", r.to_json());
|
||
return;
|
||
}
|
||
let r = job.run().unwrap_or_else(|e| die(&e.to_string()));
|
||
if flag(args, "--json") {
|
||
println!(
|
||
"{{\"bound\":{},\"final_terms\":{},\"peak_terms\":{},\"qubit\":{qubit},\"steps\":{steps},\
|
||
\"term_updates\":{},\"theta_h\":{theta},\"threshold\":{threshold},\"value\":{}}}",
|
||
r.truncation_bound, r.final_terms, r.peak_terms, r.term_updates, r.value
|
||
);
|
||
} else {
|
||
println!("<Z_{qubit}> {:.6}", r.value);
|
||
println!("bound {:.3e} (|exact - value| is at most this)", r.truncation_bound);
|
||
println!("peak terms {}", r.peak_terms);
|
||
println!("work {} term-updates", r.term_updates);
|
||
}
|
||
}
|
||
|
||
// What a fault-tolerant machine would need to run a program: code
|
||
// distance, T factories, physical qubits and wall-clock time, on the
|
||
// six reference qubits. `--receipt` seals one estimate, re-checkable
|
||
// bit for bit by anyone who re-runs it.
|
||
"estimate" => {
|
||
let (counts, default_budget) = estimate_program(args, algo);
|
||
let budget: f64 = num(args, "--budget", "", default_budget);
|
||
let slowdown: u64 = num(args, "--slowdown", "", 1);
|
||
let cfg = EstimateConfig { error_budget: budget, slowdown, ..EstimateConfig::new(budget) };
|
||
let code = opt(args, "--code", "").map(|c| match c.as_str() {
|
||
"surface" => QecModel::surface_gate(),
|
||
"surface-meas" => QecModel::surface_measurement(),
|
||
"floquet" => QecModel::floquet(),
|
||
other => die(&format!("unknown code '{other}' (surface, surface-meas, floquet)")),
|
||
});
|
||
// A code model with given constants, e.g. as `qec --fit` measured them.
|
||
let a = opt(args, "--prefactor", "").map(|v| v.parse::<f64>().unwrap_or_else(|_| die("bad --prefactor")));
|
||
let t = opt(args, "--threshold", "").map(|v| v.parse::<f64>().unwrap_or_else(|_| die("bad --threshold")));
|
||
let code = if a.is_some() || t.is_some() {
|
||
let base = code.unwrap_or_else(QecModel::surface_gate);
|
||
Some(QecModel {
|
||
name: format!("{} (given)", base.name),
|
||
prefactor: a.unwrap_or(base.prefactor),
|
||
threshold: t.unwrap_or(base.threshold),
|
||
..base
|
||
})
|
||
} else {
|
||
code
|
||
};
|
||
let which = opt(args, "--qubit", "").unwrap_or_else(|| {
|
||
if flag(args, "--receipt") || flag(args, "--frontier") { "ns-e4".into() } else { "all".into() }
|
||
});
|
||
let qubits: Vec<PhysicalQubit> =
|
||
if which == "all" { PhysicalQubit::presets() } else { vec![reference_qubit(&which)] };
|
||
if flag(args, "--receipt") {
|
||
if qubits.len() != 1 {
|
||
die("--receipt seals one estimate: pick a --qubit");
|
||
}
|
||
let job = ResourceJob { program: Program::Logical(counts), qubit: qubits[0].clone(), code, config: cfg };
|
||
let key = signing_key(args);
|
||
let r = ResourceReceipt::seal(&key, did_key(&key.verifying_key()), &job, None, |run| {
|
||
run();
|
||
0
|
||
})
|
||
.unwrap_or_else(|e| die(&format!("estimate: {e:?}")));
|
||
println!("{}", r.to_json());
|
||
return;
|
||
}
|
||
let run = |q: &PhysicalQubit, cfg: &EstimateConfig| match &code {
|
||
Some(c) => qr::estimate(&counts, q, c, cfg),
|
||
None => qr::estimate_best(&counts, q, cfg),
|
||
};
|
||
if flag(args, "--frontier") {
|
||
// Slow the program down to need fewer factories: each point
|
||
// uses fewer qubits than every faster one.
|
||
let q = &qubits[0];
|
||
let c = code.clone().unwrap_or_else(|| {
|
||
run(q, &cfg).map(|e| {
|
||
QecModel::for_qubit(q).into_iter().find(|m| m.name == e.code).expect("code of the estimate")
|
||
})
|
||
.unwrap_or_else(|e| die(&format!("estimate: {e:?}")))
|
||
});
|
||
let isa = qr::layout(&counts, budget);
|
||
let json = flag(args, "--json");
|
||
let points = qr::frontier(&isa, q, &c, &cfg, 1024);
|
||
if json {
|
||
let rows: Vec<String> = points.iter().map(estimate_json).collect();
|
||
println!("[{}]", rows.join(","));
|
||
return;
|
||
}
|
||
println!("{:>8} {:>3} {:>9} {:>14} runtime", "slowdown", "d", "factories", "physical");
|
||
for e in &points {
|
||
println!(
|
||
"{:>8} {:>3} {:>9} {:>14} {}",
|
||
e.steps / isa.min_steps,
|
||
e.distance,
|
||
e.factories,
|
||
e.physical_qubits,
|
||
human_time(e.runtime_ns)
|
||
);
|
||
}
|
||
return;
|
||
}
|
||
let estimates: Vec<_> = qubits.iter().map(|q| (q, run(q, &cfg))).collect();
|
||
if flag(args, "--json") {
|
||
let rows: Vec<String> = estimates
|
||
.iter()
|
||
.map(|(q, e)| match e {
|
||
Ok(e) => estimate_json(e),
|
||
Err(err) => format!("{{\"error\":\"{err:?}\",\"qubit\":{}}}", serde_json::to_string(&q.name).unwrap()),
|
||
})
|
||
.collect();
|
||
println!("[{}]", rows.join(","));
|
||
return;
|
||
}
|
||
let isa = qr::layout(&counts, budget);
|
||
println!("program {} tiles, {} logical steps, {} T states (budget {budget})", isa.tiles, isa.min_steps, isa.t_states);
|
||
println!(
|
||
"{:<12} {:<22} {:>3} {:>9} {:>8} {:>14} runtime",
|
||
"qubit", "code", "d", "factories", "share", "physical"
|
||
);
|
||
for (q, e) in &estimates {
|
||
match e {
|
||
Ok(e) => println!(
|
||
"{:<12} {:<22} {:>3} {:>9} {:>7.1}% {:>14} {}",
|
||
e.qubit,
|
||
e.code,
|
||
e.distance,
|
||
e.factories,
|
||
100.0 * e.factory_fraction(),
|
||
e.physical_qubits,
|
||
human_time(e.runtime_ns)
|
||
),
|
||
Err(err) => println!("{:<12} {err:?}", q.name),
|
||
}
|
||
}
|
||
}
|
||
|
||
// Error correction at circuit level: sample a noisy Clifford circuit,
|
||
// derive its detector error model, decode it exactly.
|
||
"qec" => {
|
||
use wai_quantum::quantum_frame::{error_model, parse, sample, surface_code_memory};
|
||
use wai_quantum::quantum_match::{fit_code_model, memory_experiment, MatchingGraph};
|
||
let shots: usize = num(args, "--shots", "", 100_000);
|
||
let seed: u64 = num(args, "--seed", "", 1);
|
||
let json = flag(args, "--json");
|
||
if flag(args, "--fit") {
|
||
// The surface code's logical error per round on a grid, and
|
||
// the code model a resource estimate would otherwise assume.
|
||
let mut pts = Vec::new();
|
||
if !json {
|
||
println!("{:>3} {:>7} {:>9} {:>9} {:>11}", "d", "p", "shots", "failures", "per round");
|
||
}
|
||
for d in [3u32, 5, 7] {
|
||
for p in [0.002, 0.003, 0.004, 0.005] {
|
||
let n = if d == 7 { shots / 4 } else { shots };
|
||
let t = memory_experiment(d, d, p, n, seed).unwrap_or_else(|e| die(&format!("{e:?}")));
|
||
if !json {
|
||
println!("{:>3} {:>7} {:>9} {:>9} {:>11.3e}", d, p, t.shots, t.failures, t.per_round());
|
||
}
|
||
pts.push(t);
|
||
}
|
||
}
|
||
let f = fit_code_model(&pts).unwrap_or_else(|| die("too few failures to fit; raise --shots"));
|
||
if json {
|
||
println!("{{\"prefactor\":{},\"threshold\":{}}}", f.prefactor, f.threshold);
|
||
} else {
|
||
println!("fit P(d) = a·(p/p*)^((d+1)/2) per round, a = {:.4}, p* = {:.5}", f.prefactor, f.threshold);
|
||
println!("estimate wai-quantum estimate factoring --prefactor {:.4} --threshold {:.5}", f.prefactor, f.threshold);
|
||
}
|
||
return;
|
||
}
|
||
let d: u32 = num(args, "--distance", "-d", 3);
|
||
let rounds: u32 = num(args, "--rounds", "", d);
|
||
let p: f64 = num(args, "--p", "", 0.001);
|
||
let text = match opt(args, "--file", "") {
|
||
Some(path) => std::fs::read_to_string(&path).unwrap_or_else(|e| die(&format!("{path}: {e}"))),
|
||
None => {
|
||
if d < 2 || rounds < 1 || !(0.0..=0.5).contains(&p) {
|
||
die("need -d ≥ 2, --rounds ≥ 1 and 0 ≤ --p ≤ 0.5");
|
||
}
|
||
surface_code_memory(d, rounds, p)
|
||
}
|
||
};
|
||
if flag(args, "--circuit") {
|
||
print!("{text}");
|
||
return;
|
||
}
|
||
let c = parse(&text).unwrap_or_else(|e| die(&format!("circuit: {e}")));
|
||
let model = error_model(&c).unwrap_or_else(|e| die(&format!("model: {e}")));
|
||
if flag(args, "--dem") {
|
||
print!("{}", model.to_text());
|
||
return;
|
||
}
|
||
let g = MatchingGraph::from_model(&model).unwrap_or_else(|e| die(&format!("matching graph: {e:?}")));
|
||
let s = sample(&c, shots, seed);
|
||
let decoder = opt(args, "--decoder", "").unwrap_or_else(|| "exact".into());
|
||
let failures = match decoder.as_str() {
|
||
"exact" => g.failures(&s),
|
||
"uf" => wai_quantum::quantum_uf::UnionFind::new(&g).failures(&s),
|
||
_ => die("--decoder: exact or uf"),
|
||
}
|
||
.unwrap_or_else(|e| die(&format!("decode: {e:?}")));
|
||
let rate = failures as f64 / shots as f64;
|
||
if json {
|
||
println!(
|
||
"{{\"decoder\":\"{decoder}\",\"detectors\":{},\"failures\":{failures},\"mechanisms\":{},\"observables\":{},\"qubits\":{},\"seed\":{seed},\"shots\":{shots}}}",
|
||
model.detectors,
|
||
model.mechanisms.len(),
|
||
model.observables,
|
||
c.qubits
|
||
);
|
||
} else {
|
||
println!("circuit {} qubits, {} detectors, {} observables", c.qubits, model.detectors, model.observables);
|
||
println!("model {} error mechanisms", model.mechanisms.len());
|
||
println!("failures {failures} of {shots} shots ({rate:.3e} per shot, {decoder} decoder)");
|
||
if opt(args, "--file", "").is_none() {
|
||
let per_round = wai_quantum::quantum_match::MemoryPoint { distance: d, rounds, p, shots: shots as u64, failures }.per_round();
|
||
println!("per round {per_round:.3e} (d = {d}, {rounds} rounds, p = {p})");
|
||
}
|
||
}
|
||
}
|
||
|
||
// A circuit's amplitudes by tensor-network contraction: a seeded
|
||
// hyper-optimised contraction tree, sliced to a memory cap, contracted
|
||
// in a fixed order. Circuits are JSON lines, one gate per line:
|
||
// {"q":[0,1],"m":[re,im, …]} (rows are outputs, first qubit most
|
||
// significant).
|
||
"tn" => {
|
||
use wai_quantum::quantum_tn::{contract_plan, find_path, Families, Network, PathOptions};
|
||
let path = opt(args, "--file", "").unwrap_or_else(|| die("--file CIRCUIT.jsonl is required"));
|
||
let (gates, n) = read_gate_file(&path);
|
||
let n: u32 = num(args, "--qubits", "-n", n);
|
||
let bits: Vec<u8> = match opt(args, "--bits", "") {
|
||
Some(b) => b.chars().map(|c| if c == '1' { 1 } else if c == '0' { 0 } else { die("--bits is a 0/1 string") }).collect(),
|
||
None => vec![0; n as usize],
|
||
};
|
||
let open: Vec<u32> = opt(args, "--open", "")
|
||
.map(|l| l.split(',').map(|q| q.trim().parse().unwrap_or_else(|_| die("bad --open qubit"))).collect())
|
||
.unwrap_or_default();
|
||
if open.len() > 20 {
|
||
die("at most 20 open qubits (2^20 amplitudes)");
|
||
}
|
||
// Gates split where they are block-diagonal (a control, or diagonal
|
||
// on a qubit), unless --whole-gates keeps each gate one tensor.
|
||
let whole = args.iter().any(|a| a == "--whole-gates");
|
||
let built = if whole { Network::amplitude(n, &gates, &bits, &open) } else { Network::amplitude_split(n, &gates, &bits, &open) };
|
||
let mut net = built.unwrap_or_else(|e| die(&format!("{e:?}")));
|
||
let raw = net.tensors.len();
|
||
// Projections fixed through diagonal gates, then absorption.
|
||
if whole {
|
||
net.simplify();
|
||
} else {
|
||
net.reduce();
|
||
}
|
||
let shape = net.shape();
|
||
let cap: Option<f64> = opt(args, "--cap-log2", "").map(|v| 2f64.powf(v.parse().unwrap_or_else(|_| die("bad --cap-log2"))));
|
||
let options = PathOptions {
|
||
trials: num(args, "--trials", "", 64),
|
||
seed: num(args, "--seed", "", 0),
|
||
max_size: cap,
|
||
reconfigure: 8,
|
||
families: Families::Both,
|
||
};
|
||
let plan = find_path(&shape, &net.dims, &net.open, &options);
|
||
let flops = plan.total_flops();
|
||
let plan_line = if flops > 0.0 {
|
||
format!("10^{:.3} operations, largest intermediate 2^{:.1}, {} slices", flops.log10(), plan.per_slice.max_size.log2(), plan.slices)
|
||
} else {
|
||
"no contractions needed".to_string()
|
||
};
|
||
if flag(args, "--plan-only") || flops > 1e13 {
|
||
if !flag(args, "--plan-only") {
|
||
eprintln!("wai-quantum: the plan costs 10^{:.1} operations; printing the plan only (pass a smaller circuit, or --plan-only)", flops.log10());
|
||
}
|
||
if flag(args, "--json") {
|
||
println!(
|
||
"{{\"flops\":{flops},\"log2_max_size\":{},\"slices\":{},\"tensors\":{}}}",
|
||
plan.per_slice.max_size.log2(),
|
||
plan.slices,
|
||
shape.len()
|
||
);
|
||
} else {
|
||
println!("network {raw} tensors, {} after simplification", shape.len());
|
||
println!("plan {plan_line}");
|
||
}
|
||
return;
|
||
}
|
||
let t = contract_plan(&net, &plan).unwrap_or_else(|e| die(&format!("{e:?}")));
|
||
if flag(args, "--json") {
|
||
let amps: Vec<String> = t.data.iter().map(|c| format!("[{},{}]", c.re, c.im)).collect();
|
||
println!(
|
||
"{{\"amplitudes\":[{}],\"flops\":{flops},\"log2_max_size\":{},\"slices\":{}}}",
|
||
amps.join(","),
|
||
plan.per_slice.max_size.log2(),
|
||
plan.slices
|
||
);
|
||
} else {
|
||
println!("network {raw} tensors, {} after simplification", shape.len());
|
||
println!("plan {plan_line}");
|
||
for (i, c) in t.data.iter().enumerate() {
|
||
println!("amp[{i}] {:+.12e} {:+.12e}i", c.re, c.im);
|
||
}
|
||
}
|
||
}
|
||
|
||
// Gaussian boson sampling, exactly: squeezed light through an
|
||
// interferometer with loss, and the click statistics it must show.
|
||
"gbs" => {
|
||
use wai_quantum::quantum_gbs::{random_unitary, GaussianState, Mat, C};
|
||
let r: Vec<f64> = opt(args, "--squeezing", "")
|
||
.unwrap_or_else(|| die("--squeezing r1,r2,… is required (one value per mode)"))
|
||
.split(',')
|
||
.map(|x| x.trim().parse().unwrap_or_else(|_| die("bad --squeezing value")))
|
||
.collect();
|
||
let n = r.len();
|
||
let eta: f64 = num(args, "--eta", "", 1.0);
|
||
let u = match opt(args, "--unitary", "") {
|
||
Some(path) => {
|
||
// JSON: an array of N rows, each N [re, im] pairs.
|
||
let text = std::fs::read_to_string(&path).unwrap_or_else(|e| die(&format!("{path}: {e}")));
|
||
let v: serde_json::Value = serde_json::from_str(&text).unwrap_or_else(|_| die("--unitary is not JSON"));
|
||
let rows = v.as_array().unwrap_or_else(|| die("--unitary: an array of rows"));
|
||
if rows.len() != n {
|
||
die("--unitary needs one row per mode");
|
||
}
|
||
let mut m = Mat::zeros(n);
|
||
for (i, row) in rows.iter().enumerate() {
|
||
let row = row.as_array().filter(|r| r.len() == n).unwrap_or_else(|| die("--unitary rows need one entry per mode"));
|
||
for (j, z) in row.iter().enumerate() {
|
||
let re = z[0].as_f64().unwrap_or_else(|| die("--unitary entries are [re, im]"));
|
||
let im = z[1].as_f64().unwrap_or(0.0);
|
||
m.a[i * n + j] = C { re, im };
|
||
}
|
||
}
|
||
m
|
||
}
|
||
None => random_unitary(n, num(args, "--seed", "", 0)),
|
||
};
|
||
let st = GaussianState::squeezed(&r, &u, eta).unwrap_or_else(|e| die(&format!("{e:?}")));
|
||
let json = flag(args, "--json");
|
||
if let Some(path) = opt(args, "--check", "") {
|
||
// Lines of 0/1, one click pattern per line.
|
||
let text = std::fs::read_to_string(&path).unwrap_or_else(|e| die(&format!("{path}: {e}")));
|
||
let samples: Vec<Vec<bool>> = text
|
||
.lines()
|
||
.filter(|l| !l.trim().is_empty())
|
||
.map(|l| l.trim().chars().map(|c| c == '1').collect())
|
||
.collect();
|
||
let c = st.check_samples(&samples).unwrap_or_else(|e| die(&format!("{e:?} (every line needs {n} characters)")));
|
||
if json {
|
||
println!(
|
||
"{{\"mean_clicks_exact\":{},\"mean_clicks_observed\":{},\"shots\":{},\"worst_pair_z\":{},\"worst_single_z\":{}}}",
|
||
c.mean_clicks.1, c.mean_clicks.0, c.shots, c.worst_pair, c.worst_single
|
||
);
|
||
} else {
|
||
println!("shots {}", c.shots);
|
||
println!("mean clicks {:.4} observed, {:.4} exact", c.mean_clicks.0, c.mean_clicks.1);
|
||
println!("worst |z| {:.2} over single modes, {:.2} over pairs", c.worst_single, c.worst_pair);
|
||
let verdict = if c.worst_single < 5.0 && c.worst_pair < 5.0 { "consistent" } else { "INCONSISTENT" };
|
||
println!("verdict {verdict} with the state's exact one- and two-mode click statistics");
|
||
}
|
||
return;
|
||
}
|
||
if let Some(k) = opt(args, "--sample", "") {
|
||
let k: usize = k.parse().unwrap_or_else(|_| die("bad --sample count"));
|
||
for s in st.sample_clicks(k, num(args, "--seed", "", 0) ^ 0x9e37) {
|
||
println!("{}", s.iter().map(|&b| if b { '1' } else { '0' }).collect::<String>());
|
||
}
|
||
return;
|
||
}
|
||
let (p1, _) = st.click_marginals();
|
||
let none = st.click_probability(&vec![false; n]).unwrap_or(f64::NAN);
|
||
let all = if n <= 22 { st.click_probability(&vec![true; n]).ok() } else { None };
|
||
if json {
|
||
let p1s: Vec<String> = p1.iter().map(|p| p.to_string()).collect();
|
||
println!(
|
||
"{{\"click_marginals\":[{}],\"modes\":{n},\"p_all_click\":{},\"p_no_click\":{none}}}",
|
||
p1s.join(","),
|
||
all.map_or("null".to_string(), |p| p.to_string())
|
||
);
|
||
} else {
|
||
println!("modes {n}, transmission {eta}");
|
||
println!("mean clicks {:.4}", p1.iter().sum::<f64>());
|
||
println!("P(no click) {none:.6e}");
|
||
if let Some(p) = all {
|
||
println!("P(all click) {p:.6e}");
|
||
}
|
||
for (i, p) in p1.iter().enumerate() {
|
||
println!("mode {i:>3} P(click) = {p:.6}");
|
||
}
|
||
}
|
||
}
|
||
|
||
// Ground states of spin chains by DMRG, against the exact energy
|
||
// where one is known (the Ising chain at any length, any chain by
|
||
// exact diagonalisation up to 22 sites).
|
||
"dmrg" => {
|
||
use wai_quantum::quantum_dmrg::*;
|
||
let n: usize = num(args, "--sites", "-n", 40);
|
||
if n < 2 {
|
||
die("-n: at least 2 sites");
|
||
}
|
||
let j: f64 = num(args, "--j", "", 1.0);
|
||
let h: f64 = num(args, "--h", "", if opt(args, "--chain", "").as_deref() == Some("heisenberg") { 0.0 } else { 1.0 });
|
||
let kind = opt(args, "--chain", "").unwrap_or_else(|| "ising".into());
|
||
let chain = match kind.as_str() {
|
||
"ising" => Chain::ising(n, j, h),
|
||
"heisenberg" => Chain::heisenberg(n, j, num(args, "--jz", "", j), h),
|
||
_ => die("--chain: ising or heisenberg"),
|
||
};
|
||
let cfg = DmrgConfig { max_bond: num(args, "--bond", "", 64), sweeps: num(args, "--sweeps", "", 5), cutoff: num(args, "--cutoff", "", 1e-12), ..DmrgConfig::default() };
|
||
let r = dmrg(&chain, &cfg);
|
||
let exact = if kind == "ising" {
|
||
Some(ising_exact_energy(n, j, h))
|
||
} else if flag(args, "--exact") || n <= 16 {
|
||
exact_ground_energy(&chain)
|
||
} else {
|
||
None
|
||
};
|
||
let max_bond = r.mps.dims.iter().copied().max().unwrap_or(1);
|
||
if flag(args, "--json") {
|
||
let ent: Vec<String> = r.entropies.iter().map(|e| e.to_string()).collect();
|
||
let sw: Vec<String> = r.sweep_energies.iter().map(|e| e.to_string()).collect();
|
||
println!(
|
||
"{{\"bond\":{max_bond},\"chain\":\"{kind}\",\"discarded\":{},\"energy\":{},\"entropies\":[{}],\"exact\":{},\"sites\":{n},\"sweeps\":[{}]}}",
|
||
r.discarded,
|
||
r.energy,
|
||
ent.join(","),
|
||
exact.map_or("null".into(), |e| e.to_string()),
|
||
sw.join(",")
|
||
);
|
||
} else {
|
||
println!("{kind} chain, {n} sites, bond ≤ {}, {} sweeps", cfg.max_bond, cfg.sweeps);
|
||
println!("energy {:.12} ({:.12} per site)", r.energy, r.energy / n as f64);
|
||
if let Some(e) = exact {
|
||
println!("exact {e:.12} (error per site {:.1e})", (r.energy - e) / n as f64);
|
||
}
|
||
println!("discarded {:.1e} largest bond {max_bond}", r.discarded);
|
||
println!("entropy, mid {:.6}", r.entropies[n / 2 - 1]);
|
||
}
|
||
}
|
||
|
||
// Time evolution of a spin chain by TDVP from a product state. The
|
||
// Ising chain started with every σˣ = +1 is held at each report to
|
||
// its free-fermion solution.
|
||
"tdvp" => {
|
||
use wai_quantum::quantum_dmrg::{Chain, Mps, Op, Scalar, C};
|
||
use wai_quantum::quantum_tdvp::*;
|
||
let n: usize = num(args, "--sites", "-n", 40);
|
||
if n < 2 {
|
||
die("-n: at least 2 sites");
|
||
}
|
||
let kind = opt(args, "--chain", "").unwrap_or_else(|| "ising".into());
|
||
let j: f64 = num(args, "--j", "", 1.0);
|
||
let h: f64 = num(args, "--h", "", if kind == "heisenberg" { 0.0 } else { 1.0 });
|
||
let chain = match kind.as_str() {
|
||
"ising" => Chain::ising(n, j, h),
|
||
"heisenberg" => Chain::heisenberg(n, j, num(args, "--jz", "", j), h),
|
||
_ => die("--chain: ising or heisenberg"),
|
||
};
|
||
let start_kind = opt(args, "--start", "").unwrap_or_else(|| if kind == "ising" { "plus".into() } else { "neel".into() });
|
||
let start = match start_kind.as_str() {
|
||
"plus" => Mps::product(&vec![[C::ONE, C::ONE]; n]),
|
||
"neel" => Mps::basis(&(0..n).map(|i| i % 2 == 1).collect::<Vec<_>>()),
|
||
"up" => Mps::basis(&vec![false; n]),
|
||
_ => die("--start: plus, neel or up"),
|
||
};
|
||
let dt: f64 = num(args, "--dt", "", 0.05);
|
||
let t_end: f64 = num(args, "--t", "", 2.0);
|
||
if dt.is_nan() || dt <= 0.0 || t_end.is_nan() || t_end < 0.0 {
|
||
die("--dt must be positive and --t not negative");
|
||
}
|
||
let cfg = TdvpConfig { dt, max_bond: num(args, "--bond", "", 32), cutoff: num(args, "--cutoff", "", 1e-20), ..TdvpConfig::default() };
|
||
let steps = (t_end / dt).round() as usize;
|
||
let every: usize = num(args, "--every", "", ((0.5 / dt).round() as usize).max(1));
|
||
// The observable at the middle site: σˣ for the Ising chain, σᶻ
|
||
// for the Heisenberg chain.
|
||
let (obs, obs_name) = if kind == "ising" { (Op::x(), "x") } else { (Op::z(), "z") };
|
||
let referee = kind == "ising" && start_kind == "plus";
|
||
let mut ev = Evolution::new(&chain, &start, cfg);
|
||
let e0 = ev.energy();
|
||
let mut discarded = 0.0;
|
||
let mut rows: Vec<String> = Vec::new();
|
||
let json = flag(args, "--json");
|
||
if !json {
|
||
println!("{kind} chain, {n} sites, start {start_kind}, bond ≤ {}, dt {dt}", cfg.max_bond);
|
||
print!(" t {obs_name}(mid) ");
|
||
if referee {
|
||
print!("exact max error rate exact rate ");
|
||
}
|
||
println!("S(mid) bond |ΔE| discarded");
|
||
}
|
||
for s in 1..=steps {
|
||
discarded += ev.step().discarded;
|
||
if s % every != 0 && s != steps {
|
||
continue;
|
||
}
|
||
let t = ev.time();
|
||
let o = ev.expect(obs);
|
||
let s_mid = ev.entropies()[n / 2 - 1];
|
||
let de = (ev.energy() - e0).abs();
|
||
let rate = -ev.overlap(&start).norm2().ln() / n as f64;
|
||
let ff = referee.then(|| ising_quench(n, j, h, t));
|
||
let max_err = ff.as_ref().map(|f| o.iter().zip(&f.x).map(|(p, q)| (p - q).abs()).fold(0.0, f64::max));
|
||
if json {
|
||
rows.push(format!(
|
||
"{{\"bond\":{},\"discarded\":{discarded},\"energy_drift\":{de},\"entropy\":{s_mid},\"exact\":{},\"exact_rate\":{},\"max_error\":{},\"mid\":{},\"rate\":{rate},\"t\":{t}}}",
|
||
ev.max_bond(),
|
||
ff.as_ref().map_or("null".into(), |f| f.x[n / 2].to_string()),
|
||
ff.as_ref().map_or("null".into(), |f| (-f.echo.ln() / n as f64).to_string()),
|
||
max_err.map_or("null".into(), |e| e.to_string()),
|
||
o[n / 2]
|
||
));
|
||
} else {
|
||
print!("{t:6.2} {:+.10} ", o[n / 2]);
|
||
if let (Some(f), Some(e)) = (&ff, max_err) {
|
||
print!("{:+.10} {e:.2e} {rate:.8} {:.8} ", f.x[n / 2], -f.echo.ln() / n as f64);
|
||
}
|
||
println!("{s_mid:.4} {:4} {de:.1e} {discarded:.1e}", ev.max_bond());
|
||
}
|
||
}
|
||
if json {
|
||
println!("{{\"chain\":\"{kind}\",\"dt\":{dt},\"rows\":[{}],\"sites\":{n},\"start\":\"{start_kind}\"}}", rows.join(","));
|
||
}
|
||
}
|
||
|
||
// A molecule's ground state by DMRG from its integrals, against
|
||
// Hartree–Fock and full configuration interaction where it fits.
|
||
"chem" => {
|
||
use wai_quantum::quantum_chem::*;
|
||
let path = opt(args, "--fcidump", "").unwrap_or_else(|| die("--fcidump FILE is required"));
|
||
let text = std::fs::read_to_string(&path).unwrap_or_else(|e| die(&format!("{path}: {e}")));
|
||
let fd = read_fcidump(&text).unwrap_or_else(|e| die(&format!("{path}: {e}")));
|
||
let cfg = ChemConfig { max_bond: num(args, "--bond", "", 64), sweeps: num(args, "--sweeps", "", 4), penalty: num(args, "--penalty", "", 1.0), ..ChemConfig::default() };
|
||
let mpo_bond = molecular_mpo(&fd, &cfg).max_bond();
|
||
let hf = hartree_fock_energy(&fd);
|
||
let r = ground_state(&fd, &cfg);
|
||
let full = if flag(args, "--no-fci") { None } else { fci(&fd) };
|
||
if flag(args, "--json") {
|
||
let sw: Vec<String> = r.sweep_energies.iter().map(|e| e.to_string()).collect();
|
||
println!(
|
||
"{{\"bond\":{},\"discarded\":{},\"electrons\":{},\"energy\":{},\"fci\":{},\"fci_determinants\":{},\"hartree_fock\":{hf},\"mpo_bond\":{mpo_bond},\"orbitals\":{},\"sweeps\":[{}]}}",
|
||
cfg.max_bond,
|
||
r.discarded,
|
||
fd.nelec,
|
||
r.energy,
|
||
full.as_ref().map_or("null".into(), |f| f.energy.to_string()),
|
||
full.as_ref().map_or("null".into(), |f| f.dim.to_string()),
|
||
fd.norb,
|
||
sw.join(",")
|
||
);
|
||
} else {
|
||
println!("{path}: {} orbitals, {} electrons, MPO bond {mpo_bond}, DMRG bond ≤ {}", fd.norb, fd.nelec, cfg.max_bond);
|
||
println!("Hartree–Fock {hf:.10}");
|
||
for (i, e) in r.sweep_energies.iter().enumerate() {
|
||
println!("sweep {:>2} {e:.10}", i + 1);
|
||
}
|
||
println!("discarded {:.1e}", r.discarded);
|
||
match &full {
|
||
Some(f) => println!("full CI {:.10} ({} determinants; DMRG − full CI {:.2e})", f.energy, f.dim, r.energy - f.energy),
|
||
None => println!("full CI not run"),
|
||
}
|
||
}
|
||
}
|
||
|
||
// A neural quantum state for a spin model on a chain or an open
|
||
// square lattice, against exact diagonalisation where it fits.
|
||
"nqs" => {
|
||
use wai_quantum::quantum_nqs::*;
|
||
let lx: usize = num(args, "--lx", "", 10);
|
||
let ly: usize = num(args, "--ly", "", 1);
|
||
let n = lx * ly;
|
||
if !(2..=400).contains(&n) {
|
||
die("--lx × --ly: 2 to 400 spins");
|
||
}
|
||
let periodic = flag(args, "--periodic");
|
||
let kind = opt(args, "--model", "").unwrap_or_else(|| "heisenberg".into());
|
||
let nn = if ly == 1 { chain_edges(n, periodic) } else { square_edges(lx, ly, periodic) };
|
||
let nnn = if ly == 1 { chain_next_nearest(n, periodic) } else { square_diagonals(lx, ly, periodic) };
|
||
let model = match kind.as_str() {
|
||
"heisenberg" => SpinModel::heisenberg(n, &nn, 1.0),
|
||
"j1j2" => SpinModel::j1j2(n, &nn, &nnn, 1.0, num(args, "--j2", "", 0.5)),
|
||
"ising" => SpinModel::ising(n, &nn, 1.0, num(args, "--h", "", 1.0)),
|
||
_ => die("--model: heisenberg, j1j2 or ising"),
|
||
};
|
||
let sector = if model.conserves_magnetisation() { Some(0) } else { None };
|
||
let exact = flag(args, "--exact");
|
||
if exact && n > 20 {
|
||
die("--exact enumerates every configuration: at most 20 spins");
|
||
}
|
||
let complex = flag(args, "--complex");
|
||
// Marshall's sign for the real network on the unfrustrated model;
|
||
// a complex network finds its own.
|
||
let marshall = kind == "heisenberg" && !complex;
|
||
let mut cfg = VmcConfig { samples: num(args, "--samples", "", 1000), exact, magnetisation: sector, lr: num(args, "--lr", "", 0.02), seed: num(args, "--seed", "", 7), ..VmcConfig::default() };
|
||
let steps: usize = num(args, "--steps", "", 300);
|
||
let alpha: usize = num(args, "--alpha", "", 2);
|
||
let seed: u64 = num(args, "--seed", "", 7);
|
||
let reference = if n <= 20 { model.exact_ground_energy(sector) } else { None };
|
||
// Projection onto translations at momentum 0 or π, on a square
|
||
// lattice optionally its rotations and reflections too, and a
|
||
// parity under the global flip; the lattice must be periodic.
|
||
let point = opt(args, "--point-group", "").map(|g| match g.as_str() {
|
||
"a1" => PointIrrep::A1,
|
||
"a2" => PointIrrep::A2,
|
||
"b1" => PointIrrep::B1,
|
||
"b2" => PointIrrep::B2,
|
||
_ => die("--point-group: a1, a2, b1 or b2"),
|
||
});
|
||
let reflections = opt(args, "--reflections", "").map(|g| match g.as_str() {
|
||
"a1" => RectIrrep::A1,
|
||
"a2" => RectIrrep::A2,
|
||
"b1" => RectIrrep::B1,
|
||
"b2" => RectIrrep::B2,
|
||
_ => die("--reflections: a1, a2, b1 or b2"),
|
||
});
|
||
let symmetry = flag(args, "--project").then(|| {
|
||
if !periodic {
|
||
die("--project needs --periodic");
|
||
}
|
||
let kpi = flag(args, "--kpi");
|
||
if kpi && (lx % 2 == 1 || (ly > 1 && ly % 2 == 1)) {
|
||
die("--kpi needs even lengths");
|
||
}
|
||
let base = match (point, reflections) {
|
||
(Some(_), Some(_)) => die("--point-group or --reflections, not both"),
|
||
(Some(_), None) if lx != ly || lx < 3 => die("--point-group needs a square lattice, --lx = --ly ≥ 3"),
|
||
(Some(irrep), None) => Symmetry::square(lx, kpi, irrep),
|
||
(None, Some(_)) if lx < 3 || ly < 3 => die("--reflections needs both lengths ≥ 3"),
|
||
(None, Some(irrep)) => Symmetry::rectangle(lx, ly, kpi, kpi, irrep),
|
||
(None, None) if ly == 1 => Symmetry::ring(n, kpi),
|
||
(None, None) => Symmetry::torus(lx, ly, kpi, kpi),
|
||
};
|
||
base.with_spin_flip(flag(args, "--odd"))
|
||
});
|
||
if (point.is_some() || reflections.is_some()) && symmetry.is_none() {
|
||
die("--point-group and --reflections need --project");
|
||
}
|
||
// Parameters start uniform in ±init. A projected network that
|
||
// starts near zero is nearly the projection of a constant, the
|
||
// ferromagnet: an eigenstate, which reconfiguration barely moves.
|
||
let init: f64 = num(args, "--init", "", if symmetry.is_some() { 0.1 } else { 0.01 });
|
||
// A projected network's amplitude is a sum that cancels near its
|
||
// nodes, where the local energy is large and |ψ|² rarely
|
||
// samples; a guide floor visits there.
|
||
cfg.guide = num(args, "--guide", "", if symmetry.is_some() { 0.1 } else { 0.0 });
|
||
// The exact ground energy in the projection's own sector.
|
||
let in_sector = match &symmetry {
|
||
Some(sym) if n <= 24 => model.exact_ground_energy_in(sector, sym),
|
||
_ => None,
|
||
};
|
||
fn run<T: Amplitude>(model: SpinModel, rbm: Rbm<T>, cfg: VmcConfig, steps: usize) -> (usize, Estimate, Estimate) {
|
||
let params = rbm.params();
|
||
let mut vmc = Vmc::new(model, rbm, cfg);
|
||
let mut last = vmc.estimate();
|
||
for _ in 0..steps {
|
||
last = vmc.step();
|
||
}
|
||
(params, last, vmc.estimate())
|
||
}
|
||
let (params, last, fin) = if complex {
|
||
let mut rbm = Rbm::complex(n, alpha, init, seed);
|
||
if marshall {
|
||
rbm = rbm.with_marshall_sign(checkerboard(lx, ly));
|
||
}
|
||
if let Some(sym) = symmetry {
|
||
rbm = rbm.with_symmetry(sym);
|
||
}
|
||
run(model, rbm, cfg, steps)
|
||
} else {
|
||
let mut rbm = Rbm::new(n, alpha, init, seed);
|
||
if marshall {
|
||
rbm = rbm.with_marshall_sign(checkerboard(lx, ly));
|
||
}
|
||
if let Some(sym) = symmetry {
|
||
rbm = rbm.with_symmetry(sym);
|
||
}
|
||
run(model, rbm, cfg, steps)
|
||
};
|
||
if flag(args, "--json") {
|
||
println!(
|
||
"{{\"complex\":{complex},\"energy\":{},\"error\":{},\"exact\":{},\"exact_in_sector\":{},\"last_step\":{},\"model\":\"{kind}\",\"parameters\":{params},\"sites\":{n},\"steps\":{steps},\"variance\":{}}}",
|
||
fin.energy,
|
||
fin.error,
|
||
reference.map_or("null".into(), |e| e.to_string()),
|
||
in_sector.map_or("null".into(), |e| e.to_string()),
|
||
last.energy,
|
||
fin.variance
|
||
);
|
||
} else {
|
||
println!("{kind}, {lx}×{ly}{}, {} RBM with {params} parameters, {steps} steps{}", if periodic { " periodic" } else { "" }, if complex { "complex" } else { "real" }, if exact { ", exact expectations" } else { "" });
|
||
println!("energy {:.8} ± {:.8} variance {:.2e}", fin.energy, fin.error, fin.variance);
|
||
if let Some(e) = reference {
|
||
println!("exact {e:.8} relative error {:.2e}", (fin.energy - e) / e.abs());
|
||
}
|
||
if let Some(e) = in_sector {
|
||
println!("in its sector {e:.8} relative error {:.2e}", (fin.energy - e) / e.abs());
|
||
}
|
||
}
|
||
}
|
||
|
||
// The fabric this process stands on: its detected members, and the
|
||
// fixed-point core's 32-bit-lane port certified against the native
|
||
// simulator on random circuits of the gateset `--gateset` names, then
|
||
// placed. The member lines and the certificate's member describe the
|
||
// host, so they differ between platforms; the parity and the corpus
|
||
// digest do not.
|
||
"fabric" => {
|
||
use wai_quantum::quantum_fabric::*;
|
||
let circuits: usize = num(args, "--circuits", "", 200);
|
||
let seed: u64 = num(args, "--seed", "", 7);
|
||
let gs = gateset_from(args);
|
||
let kernel = apply_kernel(gs);
|
||
let mut member = detect_local();
|
||
let cert = certify_lanes(gs, &member.id, circuits, seed);
|
||
member.kernels.push(KernelEntry { kernel: kernel.into(), implementation: "lanes32".into(), certificate: cert, energy: None });
|
||
let placed = place(std::slice::from_ref(&member), kernel, &Parity::BitExact);
|
||
let parity = |p: &Parity| match p {
|
||
Parity::BitExact => "bit-exact".to_string(),
|
||
Parity::Within { max_ulp, max_abs } => format!("within {max_ulp} ulp ({max_abs:.3e})"),
|
||
Parity::Statistical { fidelity, shots } => format!("statistical, fidelity {fidelity:.6} over {shots} shots"),
|
||
Parity::Uncertified => "uncertified".to_string(),
|
||
};
|
||
let hex = |b: &[u8]| b.iter().map(|x| format!("{x:02x}")).collect::<String>();
|
||
if flag(args, "--json") {
|
||
let kernels: Vec<String> = member
|
||
.kernels
|
||
.iter()
|
||
.map(|k| format!("{{\"implementation\":{:?},\"kernel\":{:?},\"parity\":{:?},\"generator\":{:?},\"cases\":{},\"corpus\":\"{}\"}}", k.implementation, k.kernel, parity(&k.certificate.parity), k.certificate.generator, k.certificate.cases, hex(&k.certificate.corpus)))
|
||
.collect();
|
||
println!(
|
||
"{{\"capability\":\"{}\",\"member\":{:?},\"class\":\"{}\",\"presence\":\"{}\",\"detail\":{:?},\"kernels\":[{}],\"placed\":{}}}",
|
||
gs.capability(),
|
||
member.id,
|
||
member.class.label(),
|
||
member.presence.label(),
|
||
member.detail,
|
||
kernels.join(","),
|
||
placed.map_or("null".to_string(), |(m, k)| format!("{:?}", format!("{} {}", m.id, k.implementation)))
|
||
);
|
||
} else {
|
||
println!("member {} ({}, {})", member.id, member.class.label(), member.presence.label());
|
||
println!("detail {}", member.detail);
|
||
println!("kernel {kernel}");
|
||
for k in member.kernels.iter().filter(|k| k.kernel == kernel) {
|
||
if k.implementation == REFERENCE {
|
||
println!(" {:<11}{} (the reference)", k.implementation, parity(&k.certificate.parity));
|
||
} else {
|
||
println!(" {:<11}{} over {} cases", k.implementation, parity(&k.certificate.parity), k.certificate.cases);
|
||
println!(" {:<11}corpus {} from {}", "", &hex(&k.certificate.corpus)[..16], k.certificate.generator);
|
||
}
|
||
}
|
||
match placed {
|
||
Some((m, k)) => println!("placed {} · {} (bit-exact required)", m.id, k.implementation),
|
||
None => println!("placed nothing meets bit-exact"),
|
||
}
|
||
}
|
||
}
|
||
|
||
// Quantum-inspired optimisation: Max-Cut (G-set text) or a spin
|
||
// glass, by simulated bifurcation, annealing, or exact enumeration.
|
||
// Max-Cut reports cuts; a spin glass reports energies.
|
||
"qio" => {
|
||
use wai_quantum::quantum_qio::*;
|
||
let (problem, graph) = if let Some(path) = opt(args, "--gset", "") {
|
||
let text = std::fs::read_to_string(&path).unwrap_or_else(|e| die(&format!("{path}: {e}")));
|
||
let g = read_gset(&text).unwrap_or_else(|e| die(&e.to_string()));
|
||
(g.ising(), Some(g))
|
||
} else if let Some(n) = opt(args, "--sk", "") {
|
||
let n: usize = n.parse().unwrap_or_else(|_| die("bad --sk size"));
|
||
(Ising::sk(n, num(args, "--instance", "", 0)), None)
|
||
} else {
|
||
die("qio needs --gset FILE or --sk N")
|
||
};
|
||
let score = |s: &[i8]| graph.as_ref().map_or_else(|| problem.energy(s), |g| g.cut(s));
|
||
let maximise = graph.is_some();
|
||
let method = opt(args, "--method", "").unwrap_or_else(|| "dsb".into());
|
||
let threads: usize = num(args, "--threads", "", std::thread::available_parallelism().map_or(1, |n| n.get()));
|
||
let seed: u64 = num(args, "--seed", "", 1);
|
||
let trials: usize = num(args, "--trials", "", 64);
|
||
let run = match method.as_str() {
|
||
"exact" => Run { trials: vec![ground_state(&problem, threads).unwrap_or_else(|e| die(&e.to_string()))] },
|
||
"sa" => {
|
||
let cfg = SaConfig { sweeps: num(args, "--steps", "", 10_000), beta0: num(args, "--beta0", "", 0.1), beta1: num(args, "--beta1", "", 3.0), trials, seed };
|
||
simulated_annealing(&problem, &cfg, threads)
|
||
}
|
||
m => {
|
||
let v = match m {
|
||
"bsb" => Variant::Ballistic,
|
||
"dsb" => Variant::Discrete,
|
||
"hbsb" => Variant::HeatedBallistic(0.0),
|
||
"hdsb" => Variant::HeatedDiscrete(0.0),
|
||
_ => die("--method: dsb, bsb, hdsb, hbsb, sa or exact"),
|
||
};
|
||
let mut cfg = SbConfig::published(v, num(args, "--steps", "", 10_000), trials, seed);
|
||
cfg.dt = num(args, "--dt", "", cfg.dt);
|
||
cfg.c1 = num(args, "--c1", "", cfg.c1);
|
||
if let Some(g) = opt(args, "--gamma", "") {
|
||
let g: f64 = g.parse().unwrap_or_else(|_| die("bad --gamma"));
|
||
cfg.variant = match cfg.variant {
|
||
Variant::HeatedBallistic(_) => Variant::HeatedBallistic(g),
|
||
Variant::HeatedDiscrete(_) => Variant::HeatedDiscrete(g),
|
||
other => other,
|
||
};
|
||
}
|
||
simulated_bifurcation(&problem, &cfg, threads)
|
||
}
|
||
};
|
||
let finals: Vec<Vec<i8>> = run.trials.iter().map(|t| if flag(args, "--descend") { descend(&problem, &t.spins).spins } else { t.spins.clone() }).collect();
|
||
let scores: Vec<f64> = finals.iter().map(|s| score(s)).collect();
|
||
let best_k = (0..scores.len()).fold(0, |b, k| if (maximise && scores[k] > scores[b]) || (!maximise && scores[k] < scores[b]) { k } else { b });
|
||
let best = scores[best_k];
|
||
let mean = scores.iter().sum::<f64>() / scores.len() as f64;
|
||
let target: Option<f64> = opt(args, "--target", "").map(|t| t.parse().unwrap_or_else(|_| die("bad --target")));
|
||
let hit = target.map(|t| scores.iter().filter(|&&v| if maximise { v >= t - 1e-9 } else { v <= t + 1e-9 }).count() as f64 / scores.len() as f64);
|
||
let bits: String = finals[best_k].iter().map(|&s| if s > 0 { '0' } else { '1' }).collect();
|
||
let what = if maximise { "cut" } else { "energy" };
|
||
if flag(args, "--json") {
|
||
let hit_s = hit.map_or("null".into(), |h| h.to_string());
|
||
let spins_s = if flag(args, "--spins") { format!(",\"spins\":\"{bits}\"") } else { String::new() };
|
||
println!("{{\"best\":{best},\"mean\":{mean},\"method\":\"{method}\",\"n\":{},\"objective\":\"{what}\",\"success\":{hit_s},\"trials\":{}{spins_s}}}", problem.n(), scores.len());
|
||
} else {
|
||
println!("{} spins, {method}, {} trial(s)", problem.n(), scores.len());
|
||
println!("best {what} {best}");
|
||
println!("mean {what} {mean:.3}");
|
||
if let (Some(t), Some(h)) = (target, hit) {
|
||
println!("reached {t} {:.2}% of trials", 100.0 * h);
|
||
}
|
||
if flag(args, "--spins") {
|
||
println!("spins {bits} (1 = −1; re-check: one pass over the edges)");
|
||
}
|
||
}
|
||
}
|
||
|
||
// qLDPC memory: a bivariate-bicycle code's depth-8 syndrome cycle,
|
||
// sampled with circuit-level noise and decoded by BP-OSD.
|
||
"ldpc" => {
|
||
use wai_quantum::quantum_frame::Basis;
|
||
use wai_quantum::quantum_ldpc::{memory_circuit, memory_experiment, BbCode, BpOsdConfig};
|
||
let code = match opt(args, "--code", "").as_deref().unwrap_or("72") {
|
||
"72" => BbCode::bb72(),
|
||
"90" => BbCode::bb90(),
|
||
"108" => BbCode::bb108(),
|
||
"144" => BbCode::bb144(),
|
||
"288" => BbCode::bb288(),
|
||
_ => die("--code: 72, 90, 108, 144 or 288"),
|
||
};
|
||
let p: f64 = num(args, "--p", "", 0.003);
|
||
let cycles: u32 = num(args, "--cycles", "", 6);
|
||
let bases: Vec<(&str, Basis)> = match opt(args, "--basis", "").as_deref().unwrap_or("both") {
|
||
"z" | "Z" => vec![("z", Basis::Z)],
|
||
"x" | "X" => vec![("x", Basis::X)],
|
||
"both" => vec![("z", Basis::Z), ("x", Basis::X)],
|
||
_ => die("--basis: z, x or both"),
|
||
};
|
||
if flag(args, "--circuit") {
|
||
print!("{}", memory_circuit(&code, cycles, p, bases[0].1));
|
||
return;
|
||
}
|
||
let shots: usize = num(args, "--shots", "", 1000);
|
||
let seed: u64 = num(args, "--seed", "", 0);
|
||
let threads: usize = num(args, "--threads", "", std::thread::available_parallelism().map_or(1, |n| n.get()));
|
||
let cfg = BpOsdConfig { max_iter: num(args, "--max-iter", "", 10_000), osd_order: num(args, "--osd-order", "", 7), ..BpOsdConfig::default() };
|
||
let mut rows = Vec::new();
|
||
for (i, (name, basis)) in bases.iter().enumerate() {
|
||
let r = memory_experiment(&code, cycles, p, *basis, shots, seed ^ ((i as u64 + 1) << 32), cfg, threads).unwrap_or_else(|e| die(&e.to_string()));
|
||
rows.push((*name, r));
|
||
}
|
||
let either = 1.0 - rows.iter().map(|(_, r)| 1.0 - r.per_shot()).product::<f64>();
|
||
let either_cycle = 1.0 - (1.0 - either).powf(1.0 / f64::from(cycles));
|
||
if flag(args, "--json") {
|
||
let items: Vec<String> = rows
|
||
.iter()
|
||
.map(|(b, r)| format!("{{\"basis\":\"{b}\",\"bp_converged\":{},\"failures\":{},\"per_cycle\":{},\"per_shot\":{},\"shots\":{}}}", r.converged, r.failures, r.per_cycle(), r.per_shot(), r.shots))
|
||
.collect();
|
||
println!("{{\"code\":[{},{}],\"cycles\":{cycles},\"either_per_cycle\":{either_cycle},\"either_per_shot\":{either},\"p\":{p},\"runs\":[{}]}}", code.n(), code.k(), items.join(","));
|
||
} else {
|
||
println!("[[{},{}]] bivariate-bicycle memory, p = {p}, {cycles} noisy cycles", code.n(), code.k());
|
||
for (b, r) in &rows {
|
||
println!(
|
||
"{b}-basis {} failures in {} shots {:.3e} per shot {:.3e} per cycle BP alone {:.1}%",
|
||
r.failures,
|
||
r.shots,
|
||
r.per_shot(),
|
||
r.per_cycle(),
|
||
100.0 * r.converged as f64 / r.shots as f64
|
||
);
|
||
}
|
||
if rows.len() == 2 {
|
||
println!("either {either:.3e} per shot {either_cycle:.3e} per cycle (the two halves as independent)");
|
||
}
|
||
}
|
||
}
|
||
|
||
// The architecture around the code: yoked storage, cultivated magic
|
||
// states, CCZ factories and the clock, priced; and the published
|
||
// factoring machine re-run region by region.
|
||
"arch" => {
|
||
use wai_quantum::quantum_arch::*;
|
||
let json = flag(args, "--json");
|
||
match algo.map(String::as_str) {
|
||
Some("factoring") | None => {
|
||
let bits: u32 = num(args, "--bits", "", 2048);
|
||
let params = FactoringParams::published(bits).unwrap_or_else(|| die("--bits: one of 1024, 1536, 2048, 3072, 4096, 6144, 8192"));
|
||
let e = factoring_estimate(¶ms, &Machine::published()).unwrap_or_else(|e| die(&e.to_string()));
|
||
let blocks: Vec<(u32, u32, u32, u64)> = e.cold.blocks.iter().map(|(b, k)| (b.rows - 1, b.distance, b.logical, *k)).collect();
|
||
if json {
|
||
let bs: Vec<String> = blocks.iter().map(|(w, d, l, k)| format!("{{\"count\":{k},\"distance\":{d},\"logical\":{l},\"width\":{w}}}")).collect();
|
||
println!(
|
||
"{{\"bits\":{bits},\"ccz_error\":{},\"ccz_rounds\":{},\"cold_blocks\":[{}],\"cold_logical\":{},\"cold_qubits\":{},\"compute_qubits\":{},\"days\":{},\"expected_shots\":{},\"hot_distance\":{},\"hot_logical\":{},\"hot_qubits\":{},\"hours_per_shot\":{},\"hours_per_shot_rounded\":{},\"physical_qubits\":{},\"shot_survival\":{},\"step_us\":{},\"t_state_d1\":{},\"t_state_volume\":{},\"toffolis_per_shot\":{}}}",
|
||
e.ccz_error, e.ccz_rounds, bs.join(","), e.cold_logical, e.cold_qubits, e.compute_qubits, e.days, e.expected_shots, e.hot_distance, e.hot_logical,
|
||
e.hot_qubits, e.hours_per_shot, e.hours_per_shot_rounded, e.physical_qubits(), e.shot_survival, e.timing.step_us(), e.t_state.d1, e.t_state.volume, e.toffolis_per_shot
|
||
);
|
||
} else {
|
||
println!("{bits}-bit factoring (arXiv:2505.15917 parameters), re-run");
|
||
println!("logical qubits {} cold + {} hot", e.cold_logical, e.hot_logical);
|
||
for (w, d, l, k) in &blocks {
|
||
println!("cold storage {k} × {w}×{w} grid-yoked block, d = {d}, {l} logical each");
|
||
}
|
||
println!("hot storage d = {}", e.hot_distance);
|
||
println!("physical qubits {} = {} cold + {} hot + {} compute", e.physical_qubits(), e.cold_qubits, e.hot_qubits, e.compute_qubits);
|
||
println!("T states cultivated at d1 = {}, {:.0} qubit·rounds each", e.t_state.d1, e.t_state.volume);
|
||
println!("CCZ states {:.1} rounds per factory, error {:.1e}", e.ccz_rounds, e.ccz_error);
|
||
println!("clock {:.1} µs per step ({:?}-bound)", e.timing.step_us(), e.timing.bound());
|
||
println!("Toffolis {:.4e} per shot", e.toffolis_per_shot);
|
||
println!("time {:.2} h per shot ({:.2} h at the clock), {:.2} shots, {:.2} days", e.hours_per_shot_rounded, e.hours_per_shot, e.expected_shots, e.days);
|
||
}
|
||
}
|
||
Some("storage") => {
|
||
let target: f64 = num(args, "--target", "", 1e-15);
|
||
let model = match opt(args, "--fits", "").as_deref() {
|
||
None | Some("refit") => StorageModel::uniform_refit(),
|
||
Some("uniform") => StorageModel::uniform(),
|
||
Some("si1000") => StorageModel::si1000(),
|
||
Some(_) => die("--fits: refit, uniform or si1000"),
|
||
};
|
||
let rows = [
|
||
("cold plain", best_block(&model, Storage::Cold, Yoke::None, target, 250)),
|
||
("cold rows", best_block(&model, Storage::Cold, Yoke::Rows, target, 250)),
|
||
("cold grid", best_block(&model, Storage::Cold, Yoke::Grid, target, 250)),
|
||
("hot plain", best_block(&model, Storage::Hot, Yoke::None, target, 250)),
|
||
("hot rows", best_block(&model, Storage::Hot, Yoke::Rows, target, 250)),
|
||
];
|
||
if json {
|
||
let items: Vec<String> = rows
|
||
.iter()
|
||
.filter_map(|(name, b)| b.map(|b| format!("{{\"distance\":{},\"logical\":{},\"name\":\"{name}\",\"qubits_per_logical\":{}}}", b.distance, b.logical, b.qubits_per_logical())))
|
||
.collect();
|
||
println!("{{\"layouts\":[{}],\"target\":{target}}}", items.join(","));
|
||
} else {
|
||
println!("storage at {target:e} per logical qubit per round");
|
||
for (name, b) in rows {
|
||
match b {
|
||
Some(b) => println!("{name:<11} d = {:>3} {:>4} × {:<4} patches {:>4} logical {:>8.2} qubits each", b.distance, b.rows, b.cols, b.logical, b.qubits_per_logical()),
|
||
None => println!("{name:<11} unreachable"),
|
||
}
|
||
}
|
||
}
|
||
}
|
||
Some("cultivate") => {
|
||
let p: f64 = num(args, "--p", "", 1e-3);
|
||
let target: f64 = num(args, "--target", "", 1e-7);
|
||
match cultivate(p, target) {
|
||
Some(t) if json => println!("{{\"d1\":{},\"error\":{},\"volume\":{}}}", t.d1, t.error, t.volume),
|
||
Some(t) => println!("T state at p = {p:e}: error {:.2e}, {:.0} qubit·rounds (injection distance {})", t.error, t.volume, t.d1),
|
||
None => die("no cultivation curve reaches that target at that noise (simulated: 5e-4, 1e-3, 2e-3)"),
|
||
}
|
||
}
|
||
Some(other) => die(&format!("arch: unknown '{other}' (factoring, storage, cultivate)")),
|
||
}
|
||
}
|
||
|
||
// Noisy circuits: trajectories on the state vector, optionally the
|
||
// exact density matrix (up to 12 qubits) and cross-entropy fidelity
|
||
// against the ideal circuit and the digital error model.
|
||
"noisy" => {
|
||
use wai_quantum::quantum_noisy::*;
|
||
let path = opt(args, "--file", "").unwrap_or_else(|| die("--file CIRCUIT.jsonl is required"));
|
||
let (gates, n) = read_gate_file(&path);
|
||
let n: u32 = num(args, "--qubits", "-n", n);
|
||
let noise = NoiseModel {
|
||
depol1: num(args, "--depol1", "", 0.0),
|
||
depol2: num(args, "--depol2", "", 0.0),
|
||
damping: num(args, "--damping", "", 0.0),
|
||
dephasing: num(args, "--dephasing", "", 0.0),
|
||
readout: num(args, "--readout", "", 0.0),
|
||
};
|
||
let count: usize = num(args, "--trajectories", "", 1000);
|
||
let seed: u64 = num(args, "--seed", "", 1);
|
||
let threads: usize = num(args, "--threads", "", std::thread::available_parallelism().map_or(1, |n| n.get()));
|
||
let mut names = Vec::new();
|
||
let mut observables = Vec::new();
|
||
if let Some(spec) = opt(args, "--expect", "") {
|
||
for term in spec.split(';').map(str::trim).filter(|t| !t.is_empty()) {
|
||
let string: Vec<(u32, char)> = term
|
||
.split_whitespace()
|
||
.map(|f| {
|
||
let (p, q) = f.split_at(1);
|
||
let p = p.chars().next().unwrap().to_ascii_uppercase();
|
||
if !"XYZ".contains(p) {
|
||
die("Pauli factors are X, Y or Z followed by a qubit, as in Z0");
|
||
}
|
||
(q.parse().unwrap_or_else(|_| die("bad qubit in --expect")), p)
|
||
})
|
||
.collect();
|
||
if string.iter().any(|&(q, _)| q >= n) {
|
||
die("--expect names a qubit outside the register");
|
||
}
|
||
names.push(term.to_string());
|
||
observables.push(vec![(1.0, string)]);
|
||
}
|
||
}
|
||
let run = trajectories(&gates, n, &noise, &observables, count, seed, threads).unwrap_or_else(|e| die(&format!("{e:?}")));
|
||
let exact = if flag(args, "--exact") || (flag(args, "--xeb") && n <= MAX_DENSITY_QUBITS) {
|
||
Some(DensityMatrix::run(&gates, n, &noise).unwrap_or_else(|e| die(&format!("{e:?} (the exact density matrix stops at {MAX_DENSITY_QUBITS} qubits)"))))
|
||
} else {
|
||
None
|
||
};
|
||
let xeb = if flag(args, "--xeb") {
|
||
let mut ideal = wai_quantum::quantum_sv::StateVector::zero(n).unwrap_or_else(|e| die(&format!("{e:?}")));
|
||
ideal.run(&gates, 2).unwrap_or_else(|e| die(&format!("{e:?}")));
|
||
let sampled = linear_xeb(&ideal, &run.samples);
|
||
let exact_x = exact.as_ref().map(|rho| exact_xeb(&ideal, rho, noise.readout));
|
||
Some((sampled, exact_x, digital_fidelity(&gates, n, &noise), ideal_xeb(&ideal)))
|
||
} else {
|
||
None
|
||
};
|
||
if flag(args, "--json") {
|
||
let mut parts = vec![format!("\"qubits\":{n}"), format!("\"trajectories\":{count}")];
|
||
let ex: Vec<String> = names
|
||
.iter()
|
||
.enumerate()
|
||
.map(|(k, t)| {
|
||
let e = run.expectations[k];
|
||
let exact_v = exact.as_ref().map_or("null".into(), |rho| rho.expectation(&observables[k]).to_string());
|
||
format!("{{\"exact\":{exact_v},\"mean\":{},\"stderr\":{},\"term\":{}}}", e.mean, e.stderr, serde_json::to_string(t).unwrap())
|
||
})
|
||
.collect();
|
||
parts.push(format!("\"expectations\":[{}]", ex.join(",")));
|
||
if let Some((s, e, d, i)) = xeb {
|
||
parts.push(format!("\"xeb\":{{\"digital_model\":{d},\"exact\":{},\"ideal\":{i},\"sampled\":{},\"stderr\":{}}}", e.map_or("null".into(), |v| v.to_string()), s.mean, s.stderr));
|
||
}
|
||
parts.sort();
|
||
println!("{{{}}}", parts.join(","));
|
||
} else {
|
||
println!("{n} qubits, {} gates, {count} trajectories", gates.len());
|
||
for (k, t) in names.iter().enumerate() {
|
||
let e = run.expectations[k];
|
||
match &exact {
|
||
Some(rho) => println!("<{t}> = {:+.6} ± {:.6} exact {:+.6}", e.mean, e.stderr, rho.expectation(&observables[k])),
|
||
None => println!("<{t}> = {:+.6} ± {:.6}", e.mean, e.stderr),
|
||
}
|
||
}
|
||
if let Some((s, e, d, i)) = xeb {
|
||
println!("XEB noiseless {i:.4} (the circuit's own; 1 when fully scrambling)");
|
||
println!("XEB sampled {:.4} ± {:.4} normalised {:.4}", s.mean, s.stderr, s.mean / i);
|
||
if let Some(v) = e {
|
||
println!("XEB exact {v:.4} normalised {:.4}", v / i);
|
||
}
|
||
println!("digital model {d:.4}");
|
||
}
|
||
}
|
||
}
|
||
|
||
// The dense state vector, any gate at any angle: run a circuit file
|
||
// (fused to --fuse qubits), then show the likeliest outcomes, draw
|
||
// samples, or evaluate Pauli observables.
|
||
"sv" => {
|
||
use wai_quantum::quantum_sv::StateVector;
|
||
let path = opt(args, "--file", "").unwrap_or_else(|| die("--file CIRCUIT.jsonl is required"));
|
||
let (gates, n) = read_gate_file(&path);
|
||
let n: u32 = num(args, "--qubits", "-n", n);
|
||
let width: usize = num(args, "--fuse", "", 2);
|
||
let mut sv = StateVector::zero(n).unwrap_or_else(|e| die(&format!("{e:?} (at most 32 qubits)")));
|
||
sv.run(&gates, width).unwrap_or_else(|e| die(&format!("{e:?}")));
|
||
let json = flag(args, "--json");
|
||
if let Some(k) = opt(args, "--sample", "") {
|
||
let k: usize = k.parse().unwrap_or_else(|_| die("bad --sample count"));
|
||
for x in sv.sample(k, num(args, "--seed", "", 0)) {
|
||
println!("{:0width$b}", x, width = n as usize);
|
||
}
|
||
return;
|
||
}
|
||
if let Some(spec) = opt(args, "--expect", "") {
|
||
// "Z0 X3; Y1 Y2": terms separated by ';', each a product of Pauli factors.
|
||
let mut values = Vec::new();
|
||
for term in spec.split(';').map(str::trim).filter(|t| !t.is_empty()) {
|
||
let string: Vec<(u32, char)> = term
|
||
.split_whitespace()
|
||
.map(|f| {
|
||
let (p, q) = f.split_at(1);
|
||
let p = p.chars().next().unwrap().to_ascii_uppercase();
|
||
if !"XYZ".contains(p) {
|
||
die("Pauli factors are X, Y or Z followed by a qubit, as in Z0");
|
||
}
|
||
(q.parse().unwrap_or_else(|_| die("bad qubit in --expect")), p)
|
||
})
|
||
.collect();
|
||
if string.iter().any(|&(q, _)| q >= n) {
|
||
die("--expect names a qubit outside the register");
|
||
}
|
||
values.push((term.to_string(), sv.expectation(&[(1.0, string)])));
|
||
}
|
||
if json {
|
||
let items: Vec<String> = values.iter().map(|(t, v)| format!("{}:{v}", serde_json::to_string(t).unwrap())).collect();
|
||
println!("{{{}}}", items.join(","));
|
||
} else {
|
||
for (t, v) in values {
|
||
println!("<{t}> = {v:+.12}");
|
||
}
|
||
}
|
||
return;
|
||
}
|
||
let top: usize = num(args, "--top", "", 8);
|
||
let mut idx: Vec<usize> = (0..sv.re.len()).collect();
|
||
let p = |x: usize| sv.re[x] * sv.re[x] + sv.im[x] * sv.im[x];
|
||
idx.sort_by(|&a, &b| p(b).partial_cmp(&p(a)).unwrap().then(a.cmp(&b)));
|
||
if json {
|
||
let rows: Vec<String> = idx
|
||
.iter()
|
||
.take(top)
|
||
.map(|&x| format!("{{\"amplitude\":[{},{}],\"bits\":\"{:0w$b}\",\"p\":{}}}", sv.re[x], sv.im[x], x, p(x), w = n as usize))
|
||
.collect();
|
||
println!("{{\"gates\":{},\"norm\":{},\"qubits\":{n},\"top\":[{}]}}", gates.len(), sv.norm2(), rows.join(","));
|
||
} else {
|
||
println!("{n} qubits, {} gates (fused to {width}), norm {:.15}", gates.len(), sv.norm2());
|
||
for &x in idx.iter().take(top) {
|
||
println!("|{:0w$b}> p = {:.10} amp = {:+.10} {:+.10}i", x, p(x), sv.re[x], sv.im[x], w = n as usize);
|
||
}
|
||
}
|
||
}
|
||
|
||
"entropy" => {
|
||
let c = circuit_from(args, algo);
|
||
let sv = c.simulate().unwrap_or_else(|e| die(&format!("simulate: {e:?}")));
|
||
let full = DensityMatrix::from_pure(&sv);
|
||
let rho = match opt(args, "--keep", "") {
|
||
None => full,
|
||
Some(list) => {
|
||
let keep: Vec<u8> = list
|
||
.split(',')
|
||
.map(|q| q.trim().parse().unwrap_or_else(|_| die("bad --keep qubit")))
|
||
.collect();
|
||
full.partial_trace(&keep)
|
||
}
|
||
};
|
||
let s = von_neumann_entropy(&rho);
|
||
let purity = rho.purity() as f64 / ONE as f64;
|
||
if flag(args, "--json") {
|
||
println!(
|
||
"{{\"qubits\":{},\"entropy_bits\":{s},\"purity\":{purity},\
|
||
\"schumacher_qubits_per_symbol\":{}}}",
|
||
rho.n_qubits,
|
||
schumacher_limit(&rho)
|
||
);
|
||
} else {
|
||
println!("qubits {}", rho.n_qubits);
|
||
println!("purity {purity:.6} (1 = pure, 1/d = maximally mixed)");
|
||
println!("entropy {s:.6} bits");
|
||
println!("schumacher {:.6} qubits per symbol", schumacher_limit(&rho));
|
||
}
|
||
}
|
||
|
||
// Transport. Push a state through a noisy channel and report what arrives,
|
||
// including the Holevo bound before and after — which can only fall.
|
||
"channel" => {
|
||
let name = algo.map(|s| s.as_str()).unwrap_or("depolarizing");
|
||
let p: f64 = num(args, "--prob", "", 0.5);
|
||
let p_fx = (p * ONE as f64) as i64;
|
||
let ch = match name {
|
||
"identity" => Channel::identity(),
|
||
"bitflip" => Channel::bit_flip(p_fx),
|
||
"dephasing" => Channel::dephasing(p_fx),
|
||
"depolarizing" => Channel::depolarizing(p_fx),
|
||
"damping" => Channel::amplitude_damping(p_fx),
|
||
other => die(&format!(
|
||
"unknown channel '{other}' \
|
||
(identity, bitflip, dephasing, depolarizing, damping)"
|
||
)),
|
||
};
|
||
let which = opt(args, "--state", "").unwrap_or_else(|| "plus".into());
|
||
let input = one_qubit_state(&which);
|
||
let out = ch.apply(&input).unwrap_or_else(|| die("channel/state dimension mismatch"));
|
||
|
||
// Tomography is how you learn what arrived; report how far the
|
||
// reconstruction strays, which is fixed-point rounding and nothing more.
|
||
let seen = tomography(&out);
|
||
let mut worst = 0i64;
|
||
for r in 0..out.d {
|
||
for c2 in 0..out.d {
|
||
let (a, b) = (out.get(r, c2), seen.get(r, c2));
|
||
worst = worst.max((a.re - b.re).abs()).max((a.im - b.im).abs());
|
||
}
|
||
}
|
||
// The law, on one line: noise cannot increase what a receiver can read.
|
||
let zero = one_qubit_state("zero");
|
||
let one = one_qubit_state("one");
|
||
let before = holevo_bound(&[(ONE / 2, zero.clone()), (ONE / 2, one.clone())]).unwrap();
|
||
let after = holevo_bound(&[
|
||
(ONE / 2, ch.apply(&zero).unwrap()),
|
||
(ONE / 2, ch.apply(&one).unwrap()),
|
||
])
|
||
.unwrap();
|
||
|
||
let f = |x: i64| x as f64 / ONE as f64;
|
||
let entropy = von_neumann_entropy(&out);
|
||
if flag(args, "--json") {
|
||
println!(
|
||
"{{\"channel\":\"{name}\",\"p\":{p},\"input\":\"{which}\",\
|
||
\"p0\":{},\"p1\":{},\"coherence_re\":{},\"coherence_im\":{},\
|
||
\"purity\":{},\"entropy_bits\":{entropy},\
|
||
\"holevo_before_bits\":{before},\"holevo_after_bits\":{after},\
|
||
\"tomography_max_deviation\":{}}}",
|
||
f(out.get(0, 0).re),
|
||
f(out.get(1, 1).re),
|
||
f(out.get(0, 1).re),
|
||
f(out.get(0, 1).im),
|
||
out.purity() as f64 / ONE as f64,
|
||
f(worst)
|
||
);
|
||
} else {
|
||
// The identity channel takes no parameter; printing one would imply
|
||
// a knob that does nothing.
|
||
if name == "identity" {
|
||
println!("channel {name}");
|
||
} else {
|
||
println!("channel {name} (p = {p:.3})");
|
||
}
|
||
println!("input |{which}>");
|
||
println!("defect {:.9} (0 = probability exactly conserved)", f(ch.trace_preservation_defect()));
|
||
println!("populations |0> {:.6} |1> {:.6}", f(out.get(0, 0).re), f(out.get(1, 1).re));
|
||
println!(
|
||
"coherence {:.6} {:+.6}i (0 = decohered)",
|
||
f(out.get(0, 1).re),
|
||
f(out.get(0, 1).im)
|
||
);
|
||
println!("purity {:.6}", out.purity() as f64 / ONE as f64);
|
||
println!("entropy {entropy:.6} bits");
|
||
println!("holevo {before:.6} -> {after:.6} bits (noise cannot raise it)");
|
||
println!("tomography reconstructs to {:.9}", f(worst));
|
||
}
|
||
}
|
||
|
||
// Key distribution. A quiet channel yields a perfect key; an eavesdropper
|
||
// announces themselves in the error rate, because measuring disturbs.
|
||
"bb84" => {
|
||
let rounds: u32 = num(args, "--rounds", "", 2048);
|
||
let seed: u64 = num(args, "--seed", "", 7);
|
||
let r = bb84(rounds, seed, flag(args, "--eavesdrop"));
|
||
if flag(args, "--json") {
|
||
println!(
|
||
"{{\"rounds\":{},\"sifted\":{},\"errors\":{},\"qber_ppm\":{},\
|
||
\"eavesdropper\":{}}}",
|
||
r.rounds, r.sifted, r.errors, r.qber_ppm, r.eavesdropper
|
||
);
|
||
} else {
|
||
println!("rounds {}", r.rounds);
|
||
println!("sifted {} (bases agreed)", r.sifted);
|
||
println!("errors {}", r.errors);
|
||
println!("qber {:.6}", r.qber_ppm as f64 / 1e6);
|
||
println!("eavesdropper {}", if r.eavesdropper { "yes" } else { "no" });
|
||
println!(
|
||
"no-cloning fidelity {:.6} on a conjugate-basis state \
|
||
(1.0 would break the protocol)",
|
||
clone_fidelity(true)
|
||
);
|
||
}
|
||
}
|
||
|
||
// The experiment that rules out local hidden variables, computed exactly
|
||
// rather than sampled — so the violation is arithmetic, not statistics.
|
||
"chsh" => {
|
||
let mut c = Circuit::new(2);
|
||
c.h(0).cx(0, 1);
|
||
let sv = c.simulate().unwrap_or_else(|e| die(&format!("simulate: {e:?}")));
|
||
let r = chsh(&sv).unwrap_or_else(|| die("CHSH needs a two-qubit state"));
|
||
if flag(args, "--json") {
|
||
println!(
|
||
"{{\"s\":{},\"classical_bound\":{},\"tsirelson_bound\":{},\
|
||
\"violates\":{}}}",
|
||
r.s,
|
||
r.classical_bound,
|
||
r.tsirelson_bound,
|
||
r.s > r.classical_bound
|
||
);
|
||
} else {
|
||
println!("state (|00> + |11>)/sqrt(2)");
|
||
println!("S {:.6}", r.s);
|
||
println!("classical {:.6} (any local hidden-variable theory)", r.classical_bound);
|
||
println!("tsirelson {:.6} (the quantum maximum)", r.tsirelson_bound);
|
||
println!(
|
||
"verdict {}",
|
||
if r.s > r.classical_bound {
|
||
"violates the classical bound"
|
||
} else {
|
||
"no violation"
|
||
}
|
||
);
|
||
}
|
||
}
|
||
|
||
// Qutrits. The Fourier gate over Z_3 spreads one basis state evenly across
|
||
// three, which is the d = 3 analogue of what H does for a qubit.
|
||
"qutrit" => {
|
||
let n: u8 = num(args, "--qubits", "-n", 2);
|
||
let mut q = QuditState::new(3, n).unwrap_or_else(|| die("bad qutrit register"));
|
||
q.fourier(0);
|
||
if n > 1 {
|
||
q.csum(0, 1);
|
||
}
|
||
let probs = q.probabilities();
|
||
if flag(args, "--json") {
|
||
let items: Vec<String> = probs.iter().map(|p| format!("{p}")).collect();
|
||
println!(
|
||
"{{\"d\":3,\"qutrits\":{n},\"hash\":\"{}\",\"probabilities\":[{}]}}",
|
||
hx(&q.hash()),
|
||
items.join(",")
|
||
);
|
||
} else {
|
||
println!("d 3 (qutrits)");
|
||
println!("register {n}");
|
||
println!("circuit F on qutrit 0{}", if n > 1 { ", then CSUM 0->1" } else { "" });
|
||
println!("state-hash {}", hx(&q.hash()));
|
||
println!("outcomes");
|
||
for (i, p) in probs.iter().enumerate() {
|
||
if *p > 1e-9 {
|
||
let digits: String =
|
||
(0..n).rev().map(|k| char::from(b'0' + q.digit(i, k))).collect();
|
||
println!(" |{digits}> {p:.6}");
|
||
}
|
||
}
|
||
}
|
||
}
|
||
|
||
"-h" | "--help" | "help" => usage(),
|
||
other => die(&format!("unknown command '{other}' (try --help)")),
|
||
}
|
||
}
|
||
|
||
#[cfg(test)]
|
||
mod tests {
|
||
use super::joules_line;
|
||
|
||
/// The unmetered marker prints as unmetered, never as `0 uJ`; a figure prints as itself.
|
||
#[test]
|
||
fn the_unmetered_marker_is_never_printed_as_zero_joules() {
|
||
assert!(joules_line(0).starts_with("unmetered"));
|
||
assert!(!joules_line(0).contains("0 uJ"));
|
||
assert_eq!(joules_line(1_250), "1250 uJ");
|
||
}
|
||
}
|