wai-quantum 0.3.22: sparse Pauli dynamics at utility scale, and identical bits on every platform #137
25 changed files with 2528 additions and 44 deletions
18
.github/workflows/wai-energy-receipts.yml
vendored
18
.github/workflows/wai-energy-receipts.yml
vendored
|
|
@ -19,9 +19,13 @@ name: WAI energy receipts
|
|||
# the set can build (one whose required features are all on), which also
|
||||
# proves those bins compile; and the one integration test file
|
||||
# tests/video_receipt.rs.
|
||||
# 2. wai-quantum: `cargo test --locked --features full --lib`: the crate's
|
||||
# unit tests, among them the pins and label tests of every quantum
|
||||
# receipt. The tests already marked #[ignore] stay ignored.
|
||||
# 2. wai-quantum: `cargo test --locked --features full --lib --tests`: the
|
||||
# crate's unit tests, among them the pins and label tests of every quantum
|
||||
# receipt, and its integration tests — tests/bit_identity.rs, which pins
|
||||
# one result of each f64 layer so every architecture here must produce
|
||||
# the same bits, and tests/utility_reference.rs, which holds the 127-qubit
|
||||
# results to published values. The tests already marked #[ignore] stay
|
||||
# ignored.
|
||||
# 3. wai-quantum-serve: `cargo test --locked --lib`, native: the verdict
|
||||
# logic behind POST /verify, for a receipt in either form, and that an
|
||||
# unmetered receipt's 0 is rendered as the unmetered marker (null), never
|
||||
|
|
@ -129,11 +133,13 @@ jobs:
|
|||
--lib --bins --test video_receipt
|
||||
|
||||
# wai-quantum is published separately; its receipts' legacy pins and the
|
||||
# labelled forms are tested in its own workspace, against its own lock.
|
||||
- name: wai-quantum — receipt pins and labels
|
||||
# labelled forms are tested in its own workspace, against its own lock,
|
||||
# with the integration tests that pin its f64 layers' bits on every
|
||||
# architecture this matrix runs.
|
||||
- name: wai-quantum — receipt pins, labels and bit identity
|
||||
working-directory: wai/wai-quantum
|
||||
shell: bash
|
||||
run: cargo test --locked --features full --lib
|
||||
run: cargo test --locked --features full --lib --tests
|
||||
|
||||
# The HTTP handler's POST /verify must accept a receipt in either form —
|
||||
# as sealed, or labelled as wai_quantum_meter emits it — and report the
|
||||
|
|
|
|||
2
wai/wai-quantum-serve/Cargo.lock
generated
2
wai/wai-quantum-serve/Cargo.lock
generated
|
|
@ -589,7 +589,7 @@ dependencies = [
|
|||
|
||||
[[package]]
|
||||
name = "wai-quantum"
|
||||
version = "0.3.21"
|
||||
version = "0.3.22"
|
||||
dependencies = [
|
||||
"blake3",
|
||||
"ed25519-dalek",
|
||||
|
|
|
|||
2
wai/wai-quantum/Cargo.lock
generated
2
wai/wai-quantum/Cargo.lock
generated
|
|
@ -415,7 +415,7 @@ checksum = "0b928f33d975fc6ad9f86c8f283853ad26bdd5b10b7f1542aa2fa15e2289105a"
|
|||
|
||||
[[package]]
|
||||
name = "wai-quantum"
|
||||
version = "0.3.21"
|
||||
version = "0.3.22"
|
||||
dependencies = [
|
||||
"blake3",
|
||||
"ed25519-dalek",
|
||||
|
|
|
|||
|
|
@ -1,10 +1,10 @@
|
|||
[package]
|
||||
name = "wai-quantum"
|
||||
version = "0.3.21"
|
||||
version = "0.3.22"
|
||||
edition = "2024"
|
||||
rust-version = "1.97"
|
||||
license = "Apache-2.0"
|
||||
description = "A deterministic quantum stack in pure Rust: byte-exact circuit simulation (statevector / stabilizer / tensor-network MPS / sparse-Pauli backends), error mitigation, qLDPC decoding, noise learning, circuit-equivalence proofs, a phasor interference-ML layer, information-theoretic limits, noisy channels and state tomography, and signed energy-accounted receipts. No QPU, no cloud, no system libraries — identical results native, in the browser, and as a WASI component at the edge."
|
||||
description = "A deterministic quantum stack in pure Rust: byte-exact circuit simulation (statevector / stabilizer / tensor-network MPS / sparse-Pauli backends), sparse Pauli dynamics at utility scale (arbitrary angles, 1024 qubits), error mitigation, qLDPC decoding, noise learning, circuit-equivalence proofs, a phasor interference-ML layer, information-theoretic limits, noisy channels and state tomography, and signed energy-accounted receipts. No QPU, no cloud, no system libraries — identical results native, in the browser, and as a WASI component at the edge."
|
||||
repository = "https://git.transaction.science/Transaction-Science/open-standards"
|
||||
homepage = "https://wai.transaction.science/quantum-ide"
|
||||
documentation = "https://docs.rs/wai-quantum"
|
||||
|
|
@ -39,6 +39,12 @@ quantum_toolchain = ["quantum"]
|
|||
quantum_mps = ["quantum"]
|
||||
quantum_stabilizer = ["quantum"]
|
||||
quantum_pauli = ["quantum"]
|
||||
# Sparse Pauli dynamics at utility scale: arbitrary-angle Pauli rotations on up
|
||||
# to 1024 qubits, reproducible bit for bit, with a rigorous truncation bound.
|
||||
# Independent of the dyadic core; `from_dyadic` appears when `quantum` is on.
|
||||
quantum_spd = []
|
||||
# Signed, energy-accounted receipts for it, re-checkable bit for bit anywhere.
|
||||
quantum_spd_receipt = ["quantum_spd", "_signing"]
|
||||
# Mixed states and the observables that distinguish quantum theory from a classical
|
||||
# account: density matrices, the partial trace, exact Pauli expectations, CHSH.
|
||||
quantum_info = ["quantum"]
|
||||
|
|
@ -78,13 +84,21 @@ quantum_phasor_meter = ["_signing"]
|
|||
|
||||
# Everything.
|
||||
full = [
|
||||
"quantum", "quantum_toolchain", "quantum_mps", "quantum_stabilizer", "quantum_pauli", "quantum_qir", "quantum_info", "quantum_qudit", "quantum_comm", "quantum_source", "quantum_channel",
|
||||
"quantum", "quantum_toolchain", "quantum_mps", "quantum_stabilizer", "quantum_pauli", "quantum_spd", "quantum_spd_receipt", "quantum_qir", "quantum_info", "quantum_qudit", "quantum_comm", "quantum_source", "quantum_channel",
|
||||
"quantum_receipt", "quantum_ops", "quantum_cal", "quantum_mitigate", "quantum_qec",
|
||||
"quantum_control", "quantum_noise", "quantum_compile", "quantum_atom", "quantum_qbom",
|
||||
"quantum_vml", "quantum_phasor", "quantum_qfhrr", "quantum_kernel", "quantum_qdata",
|
||||
"quantum_phasor_meter",
|
||||
]
|
||||
|
||||
[[example]]
|
||||
name = "utility"
|
||||
required-features = ["quantum_spd"]
|
||||
|
||||
[[example]]
|
||||
name = "utility_mz"
|
||||
required-features = ["quantum_spd"]
|
||||
|
||||
[[example]]
|
||||
name = "bell"
|
||||
required-features = ["quantum", "quantum_toolchain"]
|
||||
|
|
|
|||
|
|
@ -58,7 +58,9 @@ never creates information.
|
|||
| `quantum_toolchain` | algorithm library, backend recommendation, **OpenQASM 2/3 ingest of what the mainstream toolchains actually export** — both register and both measurement spellings, `rx`/`ry`/`rz`/`sx`/`u`/`u2`/`u3`/`crz`, and exact rational-multiple-of-pi angles (non-dyadic angles are refused, never rounded) — Bloch / entanglement / purity analysis |
|
||||
| `quantum_stabilizer` | stabilizer (CHP) tableau — `O(n²)`, scales far past statevector |
|
||||
| `quantum_mps` | tensor-network (matrix-product-state) backend |
|
||||
| `quantum_pauli` | sparse-Pauli / Heisenberg observable propagation |
|
||||
| `quantum_pauli` | sparse-Pauli / Heisenberg observable propagation over the dyadic gate set (up to 64 qubits) |
|
||||
| `quantum_spd` | **sparse Pauli dynamics at utility scale**: arbitrary-angle Pauli rotations on up to 1024 qubits, exact Clifford and quarter-turn rotations, a rigorous bound on what truncation can move a result, and the 127-qubit heavy-hex kicked-Ising circuits — whose exactly computed 5-step magnetization it reproduces on all 127 qubits |
|
||||
| `quantum_spd_receipt` | a signed, energy-accounted receipt over a sparse-Pauli-dynamics run, re-checkable bit for bit on any machine |
|
||||
| `quantum_receipt` | signed, energy-accounted receipt over a reconstruction |
|
||||
| `quantum_ops` | the operations/attestation layer |
|
||||
| `quantum_cal`, `quantum_control`, `quantum_noise` | calibration, filter-function robust control, DD noise spectroscopy, Cycle-Benchmarking noise learning |
|
||||
|
|
@ -93,13 +95,33 @@ the imaginary cross terms are multiplications by exactly zero.
|
|||
At 20 qubits that is 8.9x on `Z`, 7.1x on `X`, 4.6x on `CZ`, 4.3x on `CX`, 3.4x on
|
||||
`T`, 1.9x on `Y` and 1.8x on `H`.
|
||||
|
||||
## Utility scale
|
||||
|
||||
`quantum_spd` evolves an observable instead of a state, so a circuit costs what
|
||||
its non-Clifford rotations make it cost rather than `2ⁿ`:
|
||||
|
||||
```sh
|
||||
cargo run --release --example utility --features quantum_spd -- 0.2 20 1e-4
|
||||
```
|
||||
|
||||
computes ⟨Z₆₂⟩ after 20 Trotter steps of the 127-qubit heavy-hex kicked-Ising
|
||||
circuit in a fraction of a second, with the bound on how far truncation can have moved
|
||||
it. `--example utility_mz` checks the 5-step magnetization against the exactly
|
||||
computed values for every angle, and `tests/utility_reference.rs` holds the
|
||||
published references.
|
||||
|
||||
## Determinism
|
||||
|
||||
Simulation is **byte-exact**: amplitudes are dyadic fixed-point, so a circuit
|
||||
reconstructs to the same `statevector_hash` on every machine and every target.
|
||||
The approximate methods (MPS truncation, the phasor/ML layer) are
|
||||
*reproducible f64* — IEEE-754-strict, identical given the same inputs and seed —
|
||||
and are documented as such rather than claimed byte-exact.
|
||||
The approximate methods (MPS truncation, sparse Pauli dynamics, the phasor/ML
|
||||
layer) are *reproducible f64*: IEEE-754-strict, and identical given the same
|
||||
inputs and seed on every platform. They are documented as such rather than
|
||||
claimed byte-exact. Their `sin`, `cos`, `ln` and `atan2` are computed in-crate
|
||||
from `+ − × ÷`, because platform math libraries differ in the last place — enough
|
||||
that four of these layers once gave different bits natively and as a WASI
|
||||
component. `tests/bit_identity.rs` pins one result from each, so any platform
|
||||
whose arithmetic differs fails there.
|
||||
|
||||
A receipt separates the two halves of a cost honestly: **work** is
|
||||
portable-exact (`n_ops · 2^n` amplitude updates, recomputed by the verifier, so a
|
||||
|
|
|
|||
31
wai/wai-quantum/examples/utility.rs
Normal file
31
wai/wai-quantum/examples/utility.rs
Normal file
|
|
@ -0,0 +1,31 @@
|
|||
//! The 127-qubit heavy-hex kicked-Ising experiment, on a laptop.
|
||||
//!
|
||||
//! `cargo run --release --example utility --features quantum_spd -- [θ_h] [steps] [threshold]`
|
||||
//!
|
||||
//! Computes ⟨Z₆₂⟩ after `steps` Trotter steps of `RX(θ_h)` on every qubit and
|
||||
//! `RZZ(−π/2)` on every edge of the heavy-hex map. Prints the value, the bound
|
||||
//! on what truncation can have moved it, and the work done.
|
||||
|
||||
use std::time::Instant;
|
||||
use wai_quantum::quantum_spd::{SpdConfig, expect_z, heavy_hex_127, kicked_ising};
|
||||
|
||||
fn main() {
|
||||
let args: Vec<f64> = std::env::args().skip(1).map(|a| a.parse().expect("a number")).collect();
|
||||
let theta_h = args.first().copied().unwrap_or(std::f64::consts::FRAC_PI_4);
|
||||
let steps = args.get(1).copied().unwrap_or(20.0) as u32;
|
||||
let threshold = args.get(2).copied().unwrap_or(1e-3);
|
||||
|
||||
let circuit = kicked_ising(127, &heavy_hex_127(), theta_h, steps);
|
||||
let started = Instant::now();
|
||||
let r = expect_z(&circuit, 62, &SpdConfig { threshold, max_weight: None }).expect("a valid circuit");
|
||||
let secs = started.elapsed().as_secs_f64();
|
||||
|
||||
println!("theta_h {theta_h}");
|
||||
println!("steps {steps}");
|
||||
println!("threshold {threshold:e}");
|
||||
println!("<Z_62> {:.6}", r.value);
|
||||
println!("bound {:.3e} (|exact - value| is at most this)", r.truncation_bound);
|
||||
println!("peak terms {}", r.peak_terms);
|
||||
println!("updates {}", r.term_updates);
|
||||
println!("seconds {secs:.2}");
|
||||
}
|
||||
37
wai/wai-quantum/examples/utility_mz.rs
Normal file
37
wai/wai-quantum/examples/utility_mz.rs
Normal file
|
|
@ -0,0 +1,37 @@
|
|||
//! The 5-step magnetization of the heavy-hex kicked-Ising experiment, against
|
||||
//! the exactly computed values published for it (θ_h = kπ/32).
|
||||
//!
|
||||
//! `cargo run --release --example utility_mz --features quantum_spd -- [threshold]`
|
||||
|
||||
use std::time::Instant;
|
||||
use wai_quantum::quantum_spd::{Pauli, SpdConfig, expectation, heavy_hex_127, kicked_ising};
|
||||
|
||||
/// Exact Mz after 5 steps at θ_h = kπ/32, k = 1…16.
|
||||
const EXACT: [f64; 16] = [
|
||||
0.995131, 0.979836, 0.951490, 0.904843, 0.833134, 0.731464, 0.601104, 0.452804, 0.305680, 0.180317, 0.090301,
|
||||
0.037263, 0.012260, 0.002936, 0.000324, 0.0,
|
||||
];
|
||||
|
||||
fn main() {
|
||||
let threshold: f64 = std::env::args().nth(1).map_or(0.0, |a| a.parse().expect("a number"));
|
||||
let kmax: u32 = std::env::args().nth(2).map_or(16, |a| a.parse().expect("a number"));
|
||||
let edges = heavy_hex_127();
|
||||
let mz: Vec<(f64, Vec<(u32, Pauli)>)> = (0..127).map(|q| (1.0 / 127.0, vec![(q, Pauli::Z)])).collect();
|
||||
let mut worst: f64 = 0.0;
|
||||
for (k, exact) in (1..=kmax).zip(EXACT) {
|
||||
let theta_h = f64::from(k) * std::f64::consts::PI / 32.0;
|
||||
let started = Instant::now();
|
||||
let r = expectation(&kicked_ising(127, &edges, theta_h, 5), &mz, &[], &SpdConfig { threshold, max_weight: None })
|
||||
.expect("a valid circuit");
|
||||
let err = (r.value - exact).abs();
|
||||
worst = worst.max(err);
|
||||
println!(
|
||||
"k={k:2} Mz={:.6} exact={exact:.6} |err|={err:.1e} bound={:.1e} peak={} {:.2}s",
|
||||
r.value,
|
||||
r.truncation_bound,
|
||||
r.peak_terms,
|
||||
started.elapsed().as_secs_f64()
|
||||
);
|
||||
}
|
||||
println!("worst |err| = {worst:.2e} (the published values carry 6 decimals)");
|
||||
}
|
||||
|
|
@ -13,8 +13,20 @@
|
|||
//!
|
||||
//! Everything is behind a feature; enable only what you use.
|
||||
|
||||
// Every receipt seals over a Merkle root; nothing else needs one.
|
||||
#[cfg(feature = "_signing")]
|
||||
mod merkle;
|
||||
|
||||
#[cfg(any(
|
||||
feature = "quantum_spd",
|
||||
feature = "quantum_pauli",
|
||||
feature = "quantum_mps",
|
||||
feature = "quantum_vml",
|
||||
feature = "quantum_phasor",
|
||||
feature = "quantum_kernel"
|
||||
))]
|
||||
mod repro;
|
||||
|
||||
/// Re-exported so consumers get a **version-matched** signing type: the receipt
|
||||
/// sealing APIs take an `ed25519_dalek::SigningKey`, and a caller who depended on
|
||||
/// their own copy could easily pick an incompatible version.
|
||||
|
|
@ -60,6 +72,12 @@ pub mod quantum_ops;
|
|||
#[cfg(feature = "quantum_pauli")]
|
||||
pub mod quantum_pauli;
|
||||
|
||||
#[cfg(feature = "quantum_spd")]
|
||||
pub mod quantum_spd;
|
||||
|
||||
#[cfg(feature = "quantum_spd_receipt")]
|
||||
pub mod quantum_spd_receipt;
|
||||
|
||||
#[cfg(feature = "quantum_info")]
|
||||
pub mod quantum_info;
|
||||
|
||||
|
|
|
|||
|
|
@ -34,7 +34,7 @@ fn u01(s: &mut u64) -> f64 {
|
|||
fn gauss(s: &mut u64) -> f64 {
|
||||
let u1 = u01(s).max(1e-300);
|
||||
let u2 = u01(s);
|
||||
(-2.0 * u1.ln()).sqrt() * (2.0 * PI * u2).cos()
|
||||
(-2.0 * crate::repro::ln(u1)).sqrt() * crate::repro::sin_cos(2.0 * PI * u2).1
|
||||
}
|
||||
|
||||
/// A random-Fourier-feature (phasor) map approximating a Gaussian kernel
|
||||
|
|
@ -71,7 +71,7 @@ pub fn features(rff: &Rff, x: &[f64]) -> Vec<f64> {
|
|||
.zip(&rff.phases)
|
||||
.map(|(w, &b)| {
|
||||
let dot: f64 = w.iter().zip(x).map(|(wi, xi)| wi * xi).sum();
|
||||
rff.scale * (dot + b).cos()
|
||||
rff.scale * crate::repro::sin_cos(dot + b).1
|
||||
})
|
||||
.collect()
|
||||
}
|
||||
|
|
@ -86,7 +86,7 @@ pub fn features_quantized(rff: &Rff, x: &[f64], bits: u32) -> Vec<f64> {
|
|||
let dot: f64 = w.iter().zip(x).map(|(wi, xi)| wi * xi).sum();
|
||||
let t = dot + b;
|
||||
let tq = (t / (2.0 * PI) * levels).round() / levels * 2.0 * PI;
|
||||
rff.scale * tq.cos()
|
||||
rff.scale * crate::repro::sin_cos(tq).1
|
||||
})
|
||||
.collect()
|
||||
}
|
||||
|
|
@ -202,7 +202,8 @@ pub fn circles(n: usize, seed: u64) -> (Vec<Vec<f64>>, Vec<f64>) {
|
|||
let inner = i % 2 == 0;
|
||||
let ang = u01(&mut st) * 2.0 * PI;
|
||||
let r = if inner { 0.35 * u01(&mut st) } else { 0.9 + 0.35 * u01(&mut st) };
|
||||
x.push(vec![r * ang.cos() + 0.03 * gauss(&mut st), r * ang.sin() + 0.03 * gauss(&mut st)]);
|
||||
let (sin, cos) = crate::repro::sin_cos(ang);
|
||||
x.push(vec![r * cos + 0.03 * gauss(&mut st), r * sin + 0.03 * gauss(&mut st)]);
|
||||
y.push(if inner { -1.0 } else { 1.0 });
|
||||
}
|
||||
(x, y)
|
||||
|
|
|
|||
|
|
@ -10,9 +10,9 @@
|
|||
//! bounded `χ` is **approximate** — the whole point is to discard small Schmidt
|
||||
//! weight — so byte-exact-vs-analytic is not the promise here. What IS preserved is
|
||||
//! **determinism**: the contraction and truncation use only IEEE-754 `+ − × ÷ √`
|
||||
//! (no FMA). The phase gate's angle is computed with the platform `cos` and `sin`,
|
||||
//! so the result is bit-reproducible only where the math library returns the same
|
||||
//! bits for those calls. The engine reports the **max bond dimension** it used and the **retained weight**
|
||||
//! (no FMA), and the phase gate's `cos` and `sin` are computed in-crate rather than
|
||||
//! by the platform's math library, so the result is bit-identical on every machine
|
||||
//! and target. The engine reports the **max bond dimension** it used and the **retained weight**
|
||||
//! (the standard truncation-fidelity heuristic). Verified against the byte-exact dense
|
||||
//! simulator: exact agreement when `χ` is unbounded, graceful, measured degradation
|
||||
//! when it is squeezed.
|
||||
|
|
@ -367,7 +367,8 @@ fn m1(base: BaseGate, param: u16) -> [[C; 2]; 2] {
|
|||
let no = C::new(-1.0, 0.0);
|
||||
let phase = |k: u16| {
|
||||
let th = 2.0 * std::f64::consts::PI / (1u64 << k) as f64;
|
||||
C::new(th.cos(), th.sin())
|
||||
let (s, c) = crate::repro::sin_cos(th);
|
||||
C::new(c, s)
|
||||
};
|
||||
match base {
|
||||
BaseGate::I => [[o, z], [z, o]],
|
||||
|
|
|
|||
|
|
@ -15,10 +15,15 @@
|
|||
//! Real arithmetic throughout in the Hermitian-Pauli basis (`X, Y, Z` real; the
|
||||
//! rotation coefficients `cos φ, sin φ` real), and the final expectation is
|
||||
//! `Σ_{x-support = 0} c_P` (only all-`I`/`Z` strings survive `⟨0|·|0⟩`). The
|
||||
//! coefficients are `f64` and the truncation is deterministic, so the result is
|
||||
//! reproducible; its honest contract, like the MPS backend, is determinism + a
|
||||
//! stated truncation threshold, not the dense sim's byte-exact guarantee. Verified
|
||||
//! against the dense simulator's exact expectation values. Supports `n ≤ 64`.
|
||||
//! coefficients are `f64`. Terms are visited in sorted order wherever
|
||||
//! coefficients are added, and `cos φ, sin φ` are computed in-crate rather than by
|
||||
//! the platform's math library, so the result is bit-identical from run to run
|
||||
//! and across platforms. Its contract, like the MPS backend's, is that
|
||||
//! reproducibility plus a stated truncation threshold, not the dense sim's
|
||||
//! byte-exact guarantee. Verified against the dense simulator's exact
|
||||
//! expectation values. Supports `n ≤ 64` and the dyadic gate set; arbitrary
|
||||
//! angles and registers up to 1024 qubits are
|
||||
//! [`quantum_spd`](https://docs.rs/wai-quantum/latest/wai_quantum/quantum_spd/)'s.
|
||||
|
||||
use crate::quantum::{BaseGate, Circuit, Gate};
|
||||
use std::collections::HashMap;
|
||||
|
|
@ -31,6 +36,15 @@ pub struct Unsupported(pub String);
|
|||
/// A qubit's factor is `I` (00), `X` (10), `Z` (01), `Y` (11).
|
||||
pub type PauliSum = HashMap<(u64, u64), f64>;
|
||||
|
||||
/// The terms of `sum` in key order. A `HashMap` iterates in an order seeded
|
||||
/// afresh for every map; adding floating-point coefficients in that order
|
||||
/// would change the last bits of the result from one run to the next.
|
||||
fn sorted(sum: &PauliSum) -> Vec<((u64, u64), f64)> {
|
||||
let mut v: Vec<((u64, u64), f64)> = sum.iter().map(|(&k, &c)| (k, c)).collect();
|
||||
v.sort_unstable_by_key(|&(k, _)| k);
|
||||
v
|
||||
}
|
||||
|
||||
#[inline]
|
||||
fn bit(v: u64, q: u8) -> u64 {
|
||||
(v >> q) & 1
|
||||
|
|
@ -118,8 +132,8 @@ fn conjugate(sum: &PauliSum, g: &Gate, out: &mut PauliSum) -> Result<(), Unsuppo
|
|||
// diagonal rotation: split
|
||||
if let Some(phi) = is_diag_rotation(g) {
|
||||
let q = g.target;
|
||||
let (cs, sn) = (phi.cos(), phi.sin());
|
||||
for (&(x, z), &c) in sum {
|
||||
let (sn, cs) = crate::repro::sin_cos(phi);
|
||||
for ((x, z), c) in sorted(sum) {
|
||||
if bit(x, q) == 0 {
|
||||
add(x, z, c, out); // commutes
|
||||
} else {
|
||||
|
|
@ -135,14 +149,14 @@ fn conjugate(sum: &PauliSum, g: &Gate, out: &mut PauliSum) -> Result<(), Unsuppo
|
|||
// gates by control count
|
||||
match (g.base, g.controls.as_slice()) {
|
||||
(BaseGate::I, []) => {
|
||||
for (&k, &c) in sum {
|
||||
for (k, c) in sorted(sum) {
|
||||
add(k.0, k.1, c, out);
|
||||
}
|
||||
}
|
||||
(BaseGate::H | BaseGate::S | BaseGate::Sdg | BaseGate::X | BaseGate::Y | BaseGate::Z, [])
|
||||
| (BaseGate::P, []) => {
|
||||
let q = g.target;
|
||||
for (&(x, z), &c) in sum {
|
||||
for ((x, z), c) in sorted(sum) {
|
||||
let (nx, nz, s) = conj1(g.base, g.param, bit(x, q), bit(z, q));
|
||||
let (x2, z2) = set_factor(x, z, q, nx, nz);
|
||||
add(x2, z2, c * s, out);
|
||||
|
|
@ -151,7 +165,7 @@ fn conjugate(sum: &PauliSum, g: &Gate, out: &mut PauliSum) -> Result<(), Unsuppo
|
|||
(BaseGate::X, [ctrl]) => {
|
||||
// CX(c,t): x_t ^= x_c; z_c ^= z_t; sign (-1)^{x_c z_t (x_t ⊕ z_c ⊕ 1)}
|
||||
let (c, t) = (*ctrl, g.target);
|
||||
for (&(x, z), &co) in sum {
|
||||
for ((x, z), co) in sorted(sum) {
|
||||
let (xc, zc, xt, zt) = (bit(x, c), bit(z, c), bit(x, t), bit(z, t));
|
||||
let s = if (xc & zt & (xt ^ zc ^ 1)) == 1 { -1.0 } else { 1.0 };
|
||||
let (x2, _) = set_factor(x, z, t, xt ^ xc, zt);
|
||||
|
|
@ -199,7 +213,7 @@ pub fn propagate(
|
|||
|
||||
/// `⟨0…0| (sum) |0…0⟩` — the sum of coefficients of all-`I`/`Z` strings.
|
||||
pub fn expectation_zero(sum: &PauliSum) -> f64 {
|
||||
sum.iter().filter(|((x, _), _)| *x == 0).map(|(_, c)| c).sum()
|
||||
sorted(sum).into_iter().filter(|((x, _), _)| *x == 0).map(|(_, c)| c).sum()
|
||||
}
|
||||
|
||||
/// `⟨0…0| U† Z_q U |0…0⟩` — the single-qubit magnetization, via propagation.
|
||||
|
|
@ -341,6 +355,32 @@ mod tests {
|
|||
assert!(peak < 1 << 20, "Pauli sum stayed bounded: {peak}");
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn results_are_bit_identical_from_run_to_run() {
|
||||
// Each HashMap is seeded afresh, so a sum accumulated in a map's
|
||||
// iteration order would already differ between calls in one process.
|
||||
let n = 20u8;
|
||||
let mut c = Circuit::new(n);
|
||||
for layer in 0..4u8 {
|
||||
for q in 0..n {
|
||||
c.h(q);
|
||||
c.t(q);
|
||||
}
|
||||
for q in (layer % 2..n - 1).step_by(2) {
|
||||
c.cx(q, q + 1);
|
||||
}
|
||||
for q in 0..n {
|
||||
c.t(q);
|
||||
}
|
||||
}
|
||||
let (a, pa) = expect_z(&c, 10, 1e-4).unwrap();
|
||||
assert!(pa > 500, "enough terms to merge in many orders: {pa}");
|
||||
for _ in 0..6 {
|
||||
let (b, pb) = expect_z(&c, 10, 1e-4).unwrap();
|
||||
assert_eq!((a.to_bits(), pa), (b.to_bits(), pb));
|
||||
}
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn rejects_unsupported() {
|
||||
// a non-Clifford controlled-phase (CP(2)) isn't handled.
|
||||
|
|
|
|||
|
|
@ -46,7 +46,8 @@ impl C {
|
|||
}
|
||||
#[inline]
|
||||
fn expi(t: f64) -> C {
|
||||
C::new(t.cos(), t.sin())
|
||||
let (s, c) = crate::repro::sin_cos(t);
|
||||
C::new(c, s)
|
||||
}
|
||||
|
||||
// ---------------------------------------------------------------------------
|
||||
|
|
@ -56,7 +57,7 @@ fn expi(t: f64) -> C {
|
|||
/// Trainable block `W = RZ(φ)·RY(θ)` applied to a single-qubit amplitude pair.
|
||||
fn apply_w(a: &mut [C; 2], theta: f64, phi: f64) {
|
||||
// RY(θ)
|
||||
let (c, s) = ((theta * 0.5).cos(), (theta * 0.5).sin());
|
||||
let (s, c) = crate::repro::sin_cos(theta * 0.5);
|
||||
let (a0, a1) = (a[0], a[1]);
|
||||
let r0 = C::new(a0.re * c - a1.re * s, a0.im * c - a1.im * s);
|
||||
let r1 = C::new(a0.re * s + a1.re * c, a0.im * s + a1.im * c);
|
||||
|
|
|
|||
|
|
@ -222,6 +222,7 @@ fn is_zero(row: &[u64]) -> bool {
|
|||
/// [`Gf2Basis`] stops at row-echelon form — a pivot may still appear in earlier
|
||||
/// rows — which is enough to test membership but not to read a row off as an
|
||||
/// independent generator. Encoding needs the stronger form.
|
||||
#[cfg(feature = "quantum")]
|
||||
fn rref(rows_cols: &[Vec<usize>], n: usize) -> Vec<(usize, Vec<usize>)> {
|
||||
let mut m: Vec<Vec<bool>> = rows_cols
|
||||
.iter()
|
||||
|
|
|
|||
|
|
@ -42,14 +42,15 @@ fn qfhrr(a: C, bits: u32) -> C {
|
|||
if m <= 0.0 {
|
||||
return C { re: 0.0, im: 0.0 };
|
||||
}
|
||||
let theta = a.im.atan2(a.re);
|
||||
let theta = crate::repro::atan2(a.im, a.re);
|
||||
let levels = (1u64 << bits) as f64;
|
||||
let tq = (theta / (2.0 * PI) * levels).round() / levels * 2.0 * PI;
|
||||
C { re: m * tq.cos(), im: m * tq.sin() }
|
||||
let (s, c) = crate::repro::sin_cos(tq);
|
||||
C { re: m * c, im: m * s }
|
||||
}
|
||||
|
||||
fn apply_ry(amps: &mut [C], q: u8, theta: f64, bits: u32) {
|
||||
let (c, s) = ((theta * 0.5).cos(), (theta * 0.5).sin());
|
||||
let (s, c) = crate::repro::sin_cos(theta * 0.5);
|
||||
let bit = 1usize << q;
|
||||
let mut i = 0;
|
||||
while i < amps.len() {
|
||||
|
|
@ -63,7 +64,7 @@ fn apply_ry(amps: &mut [C], q: u8, theta: f64, bits: u32) {
|
|||
}
|
||||
}
|
||||
fn apply_rz(amps: &mut [C], q: u8, theta: f64, bits: u32) {
|
||||
let (cm, sm) = ((theta * 0.5).cos(), (theta * 0.5).sin());
|
||||
let (sm, cm) = crate::repro::sin_cos(theta * 0.5);
|
||||
let bit = 1usize << q;
|
||||
for (i, a) in amps.iter_mut().enumerate() {
|
||||
let (er, ei) = if i & bit == 0 { (cm, -sm) } else { (cm, sm) };
|
||||
|
|
|
|||
1268
wai/wai-quantum/src/quantum_spd.rs
Normal file
1268
wai/wai-quantum/src/quantum_spd.rs
Normal file
File diff suppressed because it is too large
Load diff
427
wai/wai-quantum/src/quantum_spd_receipt.rs
Normal file
427
wai/wai-quantum/src/quantum_spd_receipt.rs
Normal file
|
|
@ -0,0 +1,427 @@
|
|||
//! Signed, energy-accounted receipts for sparse Pauli dynamics
|
||||
//! (`wai.quantum.spd`).
|
||||
//!
|
||||
//! A receipt binds three things:
|
||||
//!
|
||||
//! - **the job**: the circuit, the observable, the initial basis state and the
|
||||
//! truncation settings. Angles and coefficients are hashed as their exact
|
||||
//! IEEE-754 bits.
|
||||
//! - **the result**: the expectation value, the truncation bound, the term
|
||||
//! counts and the work done.
|
||||
//! - **the energy** the sink measured while computing it (`joules_micro`, `0`
|
||||
//! if unmetered).
|
||||
//!
|
||||
//! The job is run inside [`SpdReceipt::seal`], under a meter the caller
|
||||
//! supplies, so the sealed result is never a caller's claim. Because the
|
||||
//! computation is reproducible bit for bit on every platform,
|
||||
//! [`SpdReceipt::verify_computation`] re-runs the job anywhere and demands the
|
||||
//! identical result: value, bound, term counts and work. The energy stays what
|
||||
//! it always is, attested by the signer's key. The receipt is a JWP profile like
|
||||
//! the crate's others (Merkle root over leaf hashes, Ed25519), and can carry an
|
||||
//! energy-class label (`quantum_energy::Labelled`).
|
||||
|
||||
use crate::merkle::merkle_root;
|
||||
use crate::quantum_spd::{Op, Pauli, PauliString, RotCircuit, SpdConfig, SpdError, SpdResult, expectation};
|
||||
use ed25519_dalek::{Signature, Signer, SigningKey, Verifier, VerifyingKey};
|
||||
|
||||
const DOMAIN_RECEIPT: &[u8] = b"wai:spd-receipt\x01";
|
||||
const DOMAIN_RECEIPT_ID: &[u8] = b"wai:spd-receipt-id\x01";
|
||||
const DOMAIN_CIRCUIT: &[u8] = b"wai:spd-circuit\x01";
|
||||
const DOMAIN_QUERY: &[u8] = b"wai:spd-query\x01";
|
||||
const DOMAIN_RESULT: &[u8] = b"wai:spd-result\x01";
|
||||
|
||||
/// Everything that determines a result: the circuit, the observable `Σ c·P`,
|
||||
/// the initial basis state `|ones⟩`, and the truncation settings.
|
||||
#[derive(Clone, Debug, PartialEq)]
|
||||
pub struct SpdJob {
|
||||
pub circuit: RotCircuit,
|
||||
pub observable: Vec<(f64, PauliString)>,
|
||||
pub ones: Vec<u32>,
|
||||
pub config: SpdConfig,
|
||||
}
|
||||
|
||||
impl SpdJob {
|
||||
pub fn run(&self) -> Result<SpdResult, SpdError> {
|
||||
expectation(&self.circuit, &self.observable, &self.ones, &self.config)
|
||||
}
|
||||
|
||||
/// BLAKE3 over the circuit's canonical encoding.
|
||||
pub fn circuit_hash(&self) -> [u8; 32] {
|
||||
let mut o = DOMAIN_CIRCUIT.to_vec();
|
||||
o.extend_from_slice(&self.circuit.n.to_be_bytes());
|
||||
o.extend_from_slice(&(self.circuit.ops.len() as u64).to_be_bytes());
|
||||
for op in &self.circuit.ops {
|
||||
encode_op(op, &mut o);
|
||||
}
|
||||
*blake3::hash(&o).as_bytes()
|
||||
}
|
||||
|
||||
/// BLAKE3 over the observable, the initial state and the settings.
|
||||
pub fn query_hash(&self) -> [u8; 32] {
|
||||
let mut o = DOMAIN_QUERY.to_vec();
|
||||
o.extend_from_slice(&(self.observable.len() as u64).to_be_bytes());
|
||||
for (c, s) in &self.observable {
|
||||
o.extend_from_slice(&c.to_bits().to_be_bytes());
|
||||
encode_string(s, &mut o);
|
||||
}
|
||||
o.extend_from_slice(&(self.ones.len() as u64).to_be_bytes());
|
||||
for q in &self.ones {
|
||||
o.extend_from_slice(&q.to_be_bytes());
|
||||
}
|
||||
o.extend_from_slice(&self.config.threshold.to_bits().to_be_bytes());
|
||||
match self.config.max_weight {
|
||||
Some(w) => {
|
||||
o.push(1);
|
||||
o.extend_from_slice(&w.to_be_bytes());
|
||||
}
|
||||
None => o.push(0),
|
||||
}
|
||||
*blake3::hash(&o).as_bytes()
|
||||
}
|
||||
}
|
||||
|
||||
fn encode_string(s: &[(u32, Pauli)], o: &mut Vec<u8>) {
|
||||
o.extend_from_slice(&(s.len() as u64).to_be_bytes());
|
||||
for (q, p) in s {
|
||||
o.extend_from_slice(&q.to_be_bytes());
|
||||
o.push(match p {
|
||||
Pauli::X => b'X',
|
||||
Pauli::Y => b'Y',
|
||||
Pauli::Z => b'Z',
|
||||
});
|
||||
}
|
||||
}
|
||||
|
||||
fn encode_op(op: &Op, o: &mut Vec<u8>) {
|
||||
let one = |tag: u8, q: u32, o: &mut Vec<u8>| {
|
||||
o.push(tag);
|
||||
o.extend_from_slice(&q.to_be_bytes());
|
||||
};
|
||||
let two = |tag: u8, a: u32, b: u32, o: &mut Vec<u8>| {
|
||||
o.push(tag);
|
||||
o.extend_from_slice(&a.to_be_bytes());
|
||||
o.extend_from_slice(&b.to_be_bytes());
|
||||
};
|
||||
match op {
|
||||
Op::H(q) => one(1, *q, o),
|
||||
Op::S(q) => one(2, *q, o),
|
||||
Op::Sdg(q) => one(3, *q, o),
|
||||
Op::SX(q) => one(4, *q, o),
|
||||
Op::SXdg(q) => one(5, *q, o),
|
||||
Op::X(q) => one(6, *q, o),
|
||||
Op::Y(q) => one(7, *q, o),
|
||||
Op::Z(q) => one(8, *q, o),
|
||||
Op::CX(a, b) => two(9, *a, *b, o),
|
||||
Op::CZ(a, b) => two(10, *a, *b, o),
|
||||
Op::Swap(a, b) => two(11, *a, *b, o),
|
||||
Op::Rot { axis, theta } => {
|
||||
o.push(12);
|
||||
encode_string(axis, o);
|
||||
o.extend_from_slice(&theta.to_bits().to_be_bytes());
|
||||
}
|
||||
Op::QuarterRot { axis, k } => {
|
||||
o.push(13);
|
||||
encode_string(axis, o);
|
||||
o.push(*k);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
fn result_hash(r: &SpdResult) -> [u8; 32] {
|
||||
let mut o = DOMAIN_RESULT.to_vec();
|
||||
o.extend_from_slice(&r.value.to_bits().to_be_bytes());
|
||||
o.extend_from_slice(&r.truncation_bound.to_bits().to_be_bytes());
|
||||
o.extend_from_slice(&(r.peak_terms as u64).to_be_bytes());
|
||||
o.extend_from_slice(&(r.final_terms as u64).to_be_bytes());
|
||||
o.extend_from_slice(&r.term_updates.to_be_bytes());
|
||||
*blake3::hash(&o).as_bytes()
|
||||
}
|
||||
|
||||
fn hx(b: &[u8]) -> String {
|
||||
b.iter().map(|x| format!("{x:02x}")).collect()
|
||||
}
|
||||
|
||||
fn unhex<const N: usize>(s: &str) -> Option<[u8; N]> {
|
||||
if s.len() != 2 * N || !s.bytes().all(|b| b.is_ascii_hexdigit()) {
|
||||
return None;
|
||||
}
|
||||
let mut out = [0u8; N];
|
||||
for (i, c) in s.as_bytes().chunks(2).enumerate() {
|
||||
out[i] = u8::from_str_radix(std::str::from_utf8(c).ok()?, 16).ok()?;
|
||||
}
|
||||
Some(out)
|
||||
}
|
||||
|
||||
/// A signed record of one sparse-Pauli-dynamics computation and its energy.
|
||||
#[derive(Clone, Debug, PartialEq)]
|
||||
pub struct SpdReceipt {
|
||||
pub circuit_hash: [u8; 32],
|
||||
pub query_hash: [u8; 32],
|
||||
pub result: SpdResult,
|
||||
pub result_hash: [u8; 32],
|
||||
/// Merkle root over (circuit, query, result) hashes.
|
||||
pub merkle_root: [u8; 32],
|
||||
/// Measured microjoules the sink spent on the run (`0` = unmetered).
|
||||
/// Attested by the signature, not verifiable by anyone else.
|
||||
pub joules_micro: u64,
|
||||
pub parent_receipt_hash: Option<[u8; 32]>,
|
||||
pub signer_pubkey: [u8; 32],
|
||||
pub signer_id: String,
|
||||
pub sig: [u8; 64],
|
||||
}
|
||||
|
||||
impl SpdReceipt {
|
||||
fn leaves(&self) -> Vec<[u8; 32]> {
|
||||
vec![self.circuit_hash, self.query_hash, self.result_hash]
|
||||
}
|
||||
|
||||
fn payload_with(&self, label: Option<&crate::quantum_energy::EnergyClass>) -> Vec<u8> {
|
||||
let mut o = crate::quantum_energy::labelled_domain(DOMAIN_RECEIPT, label);
|
||||
o.extend_from_slice(&self.circuit_hash);
|
||||
o.extend_from_slice(&self.query_hash);
|
||||
o.extend_from_slice(&self.result_hash);
|
||||
o.extend_from_slice(&self.merkle_root);
|
||||
o.extend_from_slice(&self.result.term_updates.to_be_bytes());
|
||||
o.extend_from_slice(&self.joules_micro.to_be_bytes());
|
||||
crate::quantum_energy::write_label(label, &mut o);
|
||||
match self.parent_receipt_hash {
|
||||
Some(h) => {
|
||||
o.push(1);
|
||||
o.extend_from_slice(&h);
|
||||
}
|
||||
None => o.push(0),
|
||||
}
|
||||
o.extend_from_slice(&self.signer_pubkey);
|
||||
o.extend_from_slice(self.signer_id.as_bytes());
|
||||
o
|
||||
}
|
||||
|
||||
/// Run `job` under `meter` and seal what it computed. `meter` is handed a
|
||||
/// closure that performs the run; it must call it once, and returns the
|
||||
/// microjoules it measured (`|run| { run(); 0 }` for an unmetered seal).
|
||||
/// The result sealed is the one that closure produced.
|
||||
pub fn seal(
|
||||
signer: &SigningKey,
|
||||
signer_id: impl Into<String>,
|
||||
job: &SpdJob,
|
||||
parent_receipt_hash: Option<[u8; 32]>,
|
||||
meter: impl FnOnce(&mut dyn FnMut()) -> u64,
|
||||
) -> Result<SpdReceipt, SpdError> {
|
||||
let mut out = None;
|
||||
let joules_micro = meter(&mut || out = Some(job.run()));
|
||||
let result = out.ok_or(SpdError::NotRun)??;
|
||||
let mut r = SpdReceipt {
|
||||
circuit_hash: job.circuit_hash(),
|
||||
query_hash: job.query_hash(),
|
||||
result,
|
||||
result_hash: result_hash(&result),
|
||||
merkle_root: [0; 32],
|
||||
joules_micro,
|
||||
parent_receipt_hash,
|
||||
signer_pubkey: signer.verifying_key().to_bytes(),
|
||||
signer_id: signer_id.into(),
|
||||
sig: [0; 64],
|
||||
};
|
||||
r.merkle_root = merkle_root(&r.leaves());
|
||||
r.sig = signer.sign(&r.payload_with(None)).to_bytes();
|
||||
Ok(r)
|
||||
}
|
||||
|
||||
fn body_verifies(&self) -> bool {
|
||||
result_hash(&self.result) == self.result_hash && merkle_root(&self.leaves()) == self.merkle_root
|
||||
}
|
||||
|
||||
/// The receipt is internally consistent and its signature is valid. Does
|
||||
/// not re-run the job: see [`Self::verify_computation`].
|
||||
pub fn verify(&self) -> bool {
|
||||
if !self.body_verifies() {
|
||||
return false;
|
||||
}
|
||||
let Ok(k) = VerifyingKey::from_bytes(&self.signer_pubkey) else {
|
||||
return false;
|
||||
};
|
||||
k.verify(&self.payload_with(None), &Signature::from_bytes(&self.sig)).is_ok()
|
||||
}
|
||||
|
||||
/// The re-run half of [`Self::verify_computation`].
|
||||
pub(crate) fn computation_matches(&self, job: &SpdJob) -> bool {
|
||||
job.circuit_hash() == self.circuit_hash
|
||||
&& job.query_hash() == self.query_hash
|
||||
&& job.run().is_ok_and(|r| result_hash(&r) == self.result_hash)
|
||||
}
|
||||
|
||||
/// The strong check: `job` is the job this receipt names, re-running it
|
||||
/// here reproduces the sealed result bit for bit, and the receipt verifies.
|
||||
pub fn verify_computation(&self, job: &SpdJob) -> bool {
|
||||
self.computation_matches(job) && self.verify()
|
||||
}
|
||||
|
||||
/// Stable id, for chaining a derived receipt via `parent_receipt_hash`.
|
||||
pub fn receipt_hash(&self) -> [u8; 32] {
|
||||
let mut h = blake3::Hasher::new();
|
||||
h.update(DOMAIN_RECEIPT_ID);
|
||||
h.update(&self.payload_with(None));
|
||||
h.update(&self.sig);
|
||||
*h.finalize().as_bytes()
|
||||
}
|
||||
|
||||
/// Canonical JSON, keys sorted. Floating-point figures appear twice: as a
|
||||
/// readable number, and as their exact bits (`*_bits`), which are what a
|
||||
/// verifier reads back.
|
||||
pub fn to_json(&self) -> String {
|
||||
let parent = match self.parent_receipt_hash {
|
||||
Some(h) => format!("\"{}\"", hx(&h)),
|
||||
None => "null".into(),
|
||||
};
|
||||
let r = &self.result;
|
||||
format!(
|
||||
"{{\"kind\":\"quantum-spd\",\"circuit_hash\":\"{}\",\"final_terms\":{},\"joules_micro\":{},\
|
||||
\"parent_receipt_hash\":{},\"peak_terms\":{},\"query_hash\":\"{}\",\"receipt_hash\":\"{}\",\
|
||||
\"result_hash\":\"{}\",\"root_hash\":\"{}\",\"sig\":\"{}\",\"signer_id\":{},\
|
||||
\"signer_pubkey\":\"{}\",\"term_updates\":{},\"truncation_bound\":{},\
|
||||
\"truncation_bound_bits\":\"{:016x}\",\"value\":{},\"value_bits\":\"{:016x}\"}}",
|
||||
hx(&self.circuit_hash),
|
||||
r.final_terms,
|
||||
self.joules_micro,
|
||||
parent,
|
||||
r.peak_terms,
|
||||
hx(&self.query_hash),
|
||||
hx(&self.receipt_hash()),
|
||||
hx(&self.result_hash),
|
||||
hx(&self.merkle_root),
|
||||
hx(&self.sig),
|
||||
serde_json::to_string(&self.signer_id).unwrap(),
|
||||
hx(&self.signer_pubkey),
|
||||
r.term_updates,
|
||||
serde_json::to_string(&r.truncation_bound).unwrap(),
|
||||
r.truncation_bound.to_bits(),
|
||||
serde_json::to_string(&r.value).unwrap(),
|
||||
r.value.to_bits(),
|
||||
)
|
||||
}
|
||||
|
||||
pub fn from_json(s: &str) -> Option<SpdReceipt> {
|
||||
let v: serde_json::Value = serde_json::from_str(s).ok()?;
|
||||
let o = v.as_object()?;
|
||||
if o.get("kind")?.as_str()? != "quantum-spd" {
|
||||
return None;
|
||||
}
|
||||
let u = |k: &str| o.get(k).and_then(|x| x.as_u64());
|
||||
let s = |k: &str| o.get(k).and_then(|x| x.as_str());
|
||||
let bits = |k: &str| u64::from_str_radix(s(k)?, 16).ok().map(f64::from_bits);
|
||||
let parent = match o.get("parent_receipt_hash") {
|
||||
Some(serde_json::Value::String(h)) => Some(unhex::<32>(h)?),
|
||||
_ => None,
|
||||
};
|
||||
Some(SpdReceipt {
|
||||
circuit_hash: unhex(s("circuit_hash")?)?,
|
||||
query_hash: unhex(s("query_hash")?)?,
|
||||
result: SpdResult {
|
||||
value: bits("value_bits")?,
|
||||
truncation_bound: bits("truncation_bound_bits")?,
|
||||
peak_terms: usize::try_from(u("peak_terms")?).ok()?,
|
||||
final_terms: usize::try_from(u("final_terms")?).ok()?,
|
||||
term_updates: u("term_updates")?,
|
||||
},
|
||||
result_hash: unhex(s("result_hash")?)?,
|
||||
merkle_root: unhex(s("root_hash")?)?,
|
||||
joules_micro: u("joules_micro")?,
|
||||
parent_receipt_hash: parent,
|
||||
signer_pubkey: unhex(s("signer_pubkey")?)?,
|
||||
signer_id: s("signer_id")?.to_string(),
|
||||
sig: unhex(s("sig")?)?,
|
||||
})
|
||||
}
|
||||
}
|
||||
|
||||
crate::quantum_energy::labellable!(SpdReceipt, DOMAIN_RECEIPT_ID);
|
||||
|
||||
#[cfg(test)]
|
||||
mod tests {
|
||||
use super::*;
|
||||
use crate::quantum_spd::{heavy_hex_127, kicked_ising};
|
||||
|
||||
fn key(b: u8) -> SigningKey {
|
||||
SigningKey::from_bytes(&[b; 32])
|
||||
}
|
||||
|
||||
fn job() -> SpdJob {
|
||||
SpdJob {
|
||||
circuit: kicked_ising(127, &heavy_hex_127(), 0.4, 4),
|
||||
observable: vec![(1.0, vec![(62, Pauli::Z)])],
|
||||
ones: vec![],
|
||||
config: SpdConfig { threshold: 1e-4, max_weight: None },
|
||||
}
|
||||
}
|
||||
|
||||
fn unmetered(run: &mut dyn FnMut()) -> u64 {
|
||||
run();
|
||||
0
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn a_sealed_receipt_verifies_and_reproduces() {
|
||||
let j = job();
|
||||
let r = SpdReceipt::seal(&key(1), "sink-a", &j, None, |run| {
|
||||
run();
|
||||
4_200
|
||||
})
|
||||
.unwrap();
|
||||
assert!(r.verify());
|
||||
assert!(r.verify_computation(&j));
|
||||
assert_eq!(r.joules_micro, 4_200);
|
||||
assert_eq!(r.result, j.run().unwrap());
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn the_json_round_trips_to_the_same_bytes() {
|
||||
let r = SpdReceipt::seal(&key(2), "sink \"b\"", &job(), Some([7; 32]), unmetered).unwrap();
|
||||
let back = SpdReceipt::from_json(&r.to_json()).unwrap();
|
||||
assert_eq!(back, r);
|
||||
assert_eq!(back.to_json(), r.to_json());
|
||||
assert!(back.verify());
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn tampering_is_caught() {
|
||||
let j = job();
|
||||
let r = SpdReceipt::seal(&key(3), "sink", &j, None, unmetered).unwrap();
|
||||
let mut a = r.clone();
|
||||
a.joules_micro = 1;
|
||||
assert!(!a.verify(), "the energy is signed");
|
||||
let mut b = r.clone();
|
||||
b.result.value += 1e-12;
|
||||
assert!(!b.verify(), "the result is hashed");
|
||||
let mut c = r.clone();
|
||||
c.result.term_updates -= 1;
|
||||
assert!(!c.verify(), "the work is hashed");
|
||||
// A different job — one angle moved by one ulp — is not the job sealed.
|
||||
let mut other = j.clone();
|
||||
if let Op::Rot { theta, .. } = &mut other.circuit.ops[5] {
|
||||
*theta = f64::from_bits(theta.to_bits() + 1);
|
||||
}
|
||||
assert!(!r.verify_computation(&other));
|
||||
let mut looser = j.clone();
|
||||
looser.config.threshold = 2e-4;
|
||||
assert!(!r.verify_computation(&looser), "the settings are part of the job");
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn a_meter_that_never_runs_the_job_cannot_seal() {
|
||||
assert_eq!(SpdReceipt::seal(&key(4), "sink", &job(), None, |_run| 99), Err(SpdError::NotRun));
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn an_energy_class_label_seals_over_it() {
|
||||
use crate::quantum_energy::{EnergyClass, Labelled};
|
||||
let r = SpdReceipt::seal(&key(5), "sink", &job(), None, |run| {
|
||||
run();
|
||||
1_000
|
||||
})
|
||||
.unwrap();
|
||||
let class = EnergyClass::ModelBased;
|
||||
let labelled = Labelled::seal(r, class, &key(5)).unwrap();
|
||||
assert!(labelled.verify());
|
||||
assert!(!labelled.receipt.verify(), "the inner signature now covers the label");
|
||||
}
|
||||
}
|
||||
|
|
@ -65,7 +65,7 @@ fn rand_angle(s: &mut u64) -> f64 {
|
|||
|
||||
#[inline]
|
||||
fn apply_ry(amps: &mut [C], n: u8, q: u8, theta: f64) {
|
||||
let (c, s) = ((theta * 0.5).cos(), (theta * 0.5).sin());
|
||||
let (s, c) = crate::repro::sin_cos(theta * 0.5);
|
||||
let bit = 1usize << q;
|
||||
let dim = 1usize << n;
|
||||
let mut i = 0;
|
||||
|
|
@ -81,7 +81,7 @@ fn apply_ry(amps: &mut [C], n: u8, q: u8, theta: f64) {
|
|||
}
|
||||
#[inline]
|
||||
fn apply_rz(amps: &mut [C], q: u8, theta: f64) {
|
||||
let (cm, sm) = ((theta * 0.5).cos(), (theta * 0.5).sin());
|
||||
let (sm, cm) = crate::repro::sin_cos(theta * 0.5);
|
||||
let e0 = C { re: cm, im: -sm }; // e^{-iθ/2}
|
||||
let e1 = C { re: cm, im: sm }; // e^{+iθ/2}
|
||||
let bit = 1usize << q;
|
||||
|
|
@ -148,7 +148,7 @@ fn kl(target: &[f64], model: &[f64]) -> f64 {
|
|||
let mut s = 0.0;
|
||||
for (&q, &p) in target.iter().zip(model) {
|
||||
if q > 0.0 {
|
||||
s += q * (q / p.max(1e-12)).ln();
|
||||
s += q * crate::repro::ln(q / p.max(1e-12));
|
||||
}
|
||||
}
|
||||
s
|
||||
|
|
|
|||
387
wai/wai-quantum/src/repro.rs
Normal file
387
wai/wai-quantum/src/repro.rs
Normal file
|
|
@ -0,0 +1,387 @@
|
|||
//! Transcendental functions that give the same bits on every platform.
|
||||
//!
|
||||
//! The platform's libm may differ from another's in the last place. That is
|
||||
//! enough to break a result promised bit-identical natively, in the browser and
|
||||
//! at the edge, and it did: four layers of this crate computed different bits
|
||||
//! natively and as a WASI component until they moved here. Everything below is
|
||||
//! fdlibm's algorithm, as musl carries it, using only `+ − × ÷`, which Rust
|
||||
//! never fuses.
|
||||
|
||||
/// `sin` and `cos` of `x`, bit-identical on every IEEE-754 platform. The
|
||||
/// algorithm is fdlibm's (as musl carries it): Cody–Waite reduction by π/2 in
|
||||
/// three parts, then minimax polynomials on [−π/4, π/4]. Only `+ − × ÷` are
|
||||
/// used, and Rust never fuses them, so the bits do not depend on the target.
|
||||
/// Accurate to within an ulp of the true value.
|
||||
///
|
||||
/// For `|x| > 2²⁰` (and for non-finite `x`) both are NaN: the reduction used
|
||||
/// here is exact only below that, and a NaN is better than a value that
|
||||
/// depends on the platform.
|
||||
pub fn sin_cos(x: f64) -> (f64, f64) {
|
||||
if !(x.abs() <= 1_048_576.0) {
|
||||
return (f64::NAN, f64::NAN);
|
||||
}
|
||||
let ix = (x.to_bits() >> 32) as u32 & 0x7fff_ffff;
|
||||
if ix <= 0x3fe9_21fb {
|
||||
// |x| ≤ π/4
|
||||
if ix < 0x3e40_0000 {
|
||||
// |x| < 2⁻²⁷: sin x = x and cos x = 1 to the last bit.
|
||||
return (x, 1.0);
|
||||
}
|
||||
return (k_sin(x, 0.0, false), k_cos(x, 0.0));
|
||||
}
|
||||
let (n, y0, y1) = rem_pio2(x, ix);
|
||||
let (s, c) = (k_sin(y0, y1, true), k_cos(y0, y1));
|
||||
match n & 3 {
|
||||
0 => (s, c),
|
||||
1 => (c, -s),
|
||||
2 => (-s, -c),
|
||||
_ => (-c, s),
|
||||
}
|
||||
}
|
||||
|
||||
const S1: f64 = -1.666_666_666_666_663_243_48e-1;
|
||||
const S2: f64 = 8.333_333_333_322_489_461_24e-3;
|
||||
const S3: f64 = -1.984_126_982_985_794_931_34e-4;
|
||||
const S4: f64 = 2.755_731_370_707_006_767_89e-6;
|
||||
const S5: f64 = -2.505_076_025_340_686_341_95e-8;
|
||||
const S6: f64 = 1.589_690_995_211_550_102_21e-10;
|
||||
|
||||
fn k_sin(x: f64, y: f64, tail: bool) -> f64 {
|
||||
let z = x * x;
|
||||
let w = z * z;
|
||||
let r = S2 + z * (S3 + z * S4) + z * w * (S5 + z * S6);
|
||||
let v = z * x;
|
||||
if !tail {
|
||||
x + v * (S1 + z * r)
|
||||
} else {
|
||||
x - ((z * (0.5 * y - v * r) - y) - v * S1)
|
||||
}
|
||||
}
|
||||
|
||||
const C1: f64 = 4.166_666_666_666_660_190_37e-2;
|
||||
const C2: f64 = -1.388_888_888_887_410_957_49e-3;
|
||||
const C3: f64 = 2.480_158_728_947_672_941_78e-5;
|
||||
const C4: f64 = -2.755_731_435_139_066_330_35e-7;
|
||||
const C5: f64 = 2.087_572_321_298_174_827_90e-9;
|
||||
const C6: f64 = -1.135_964_755_778_819_482_65e-11;
|
||||
|
||||
fn k_cos(x: f64, y: f64) -> f64 {
|
||||
let z = x * x;
|
||||
let w = z * z;
|
||||
let r = z * (C1 + z * (C2 + z * C3)) + w * w * (C4 + z * (C5 + z * C6));
|
||||
let hz = 0.5 * z;
|
||||
let w = 1.0 - hz;
|
||||
w + (((1.0 - w) - hz) + (z * r - x * y))
|
||||
}
|
||||
|
||||
const TOINT: f64 = 1.5 / f64::EPSILON;
|
||||
const PIO4: f64 = 7.853_981_633_974_482_789_99e-1;
|
||||
const INVPIO2: f64 = 6.366_197_723_675_813_824_33e-1;
|
||||
const PIO2_1: f64 = 1.570_796_326_734_125_614_17;
|
||||
const PIO2_1T: f64 = 6.077_100_506_506_192_249_32e-11;
|
||||
const PIO2_2: f64 = 6.077_100_506_303_965_976_60e-11;
|
||||
const PIO2_2T: f64 = 2.022_266_248_795_950_631_54e-21;
|
||||
const PIO2_3: f64 = 2.022_266_248_711_166_455_80e-21;
|
||||
const PIO2_3T: f64 = 8.478_427_660_368_899_569_97e-32;
|
||||
|
||||
/// `x = n·π/2 + (y0 + y1)` with `|y0 + y1| ≤ π/4`, for `|x| ≤ 2²⁰`.
|
||||
fn rem_pio2(x: f64, ix: u32) -> (i64, f64, f64) {
|
||||
let mut f = (x * INVPIO2 + TOINT) - TOINT;
|
||||
let mut r = x - f * PIO2_1;
|
||||
let mut w = f * PIO2_1T;
|
||||
// Matters only under directed rounding; kept so the port is exact.
|
||||
if r - w < -PIO4 {
|
||||
f -= 1.0;
|
||||
r = x - f * PIO2_1;
|
||||
w = f * PIO2_1T;
|
||||
} else if r - w > PIO4 {
|
||||
f += 1.0;
|
||||
r = x - f * PIO2_1;
|
||||
w = f * PIO2_1T;
|
||||
}
|
||||
let mut y0 = r - w;
|
||||
let ex = (ix >> 20) as i32;
|
||||
let ey = ((y0.to_bits() >> 52) & 0x7ff) as i32;
|
||||
if ex - ey > 16 {
|
||||
let t = r;
|
||||
w = f * PIO2_2;
|
||||
r = t - w;
|
||||
w = f * PIO2_2T - ((t - r) - w);
|
||||
y0 = r - w;
|
||||
let ey = ((y0.to_bits() >> 52) & 0x7ff) as i32;
|
||||
if ex - ey > 49 {
|
||||
let t = r;
|
||||
w = f * PIO2_3;
|
||||
r = t - w;
|
||||
w = f * PIO2_3T - ((t - r) - w);
|
||||
y0 = r - w;
|
||||
}
|
||||
}
|
||||
let y1 = (r - y0) - w;
|
||||
(f as i64, y0, y1)
|
||||
}
|
||||
|
||||
// ---------------------------------------------------------------------------
|
||||
// Natural logarithm
|
||||
// ---------------------------------------------------------------------------
|
||||
|
||||
#[cfg(any(feature = "quantum_vml", feature = "quantum_kernel"))]
|
||||
const LN2_HI: f64 = 6.931_471_803_691_238_164_90e-1;
|
||||
#[cfg(any(feature = "quantum_vml", feature = "quantum_kernel"))]
|
||||
const LN2_LO: f64 = 1.908_214_929_270_587_700_02e-10;
|
||||
#[cfg(any(feature = "quantum_vml", feature = "quantum_kernel"))]
|
||||
const LG1: f64 = 6.666_666_666_666_735_130e-1;
|
||||
#[cfg(any(feature = "quantum_vml", feature = "quantum_kernel"))]
|
||||
const LG2: f64 = 3.999_999_999_940_941_908e-1;
|
||||
#[cfg(any(feature = "quantum_vml", feature = "quantum_kernel"))]
|
||||
const LG3: f64 = 2.857_142_874_366_239_149e-1;
|
||||
#[cfg(any(feature = "quantum_vml", feature = "quantum_kernel"))]
|
||||
const LG4: f64 = 2.222_219_843_214_978_396e-1;
|
||||
#[cfg(any(feature = "quantum_vml", feature = "quantum_kernel"))]
|
||||
const LG5: f64 = 1.818_357_216_161_805_012e-1;
|
||||
#[cfg(any(feature = "quantum_vml", feature = "quantum_kernel"))]
|
||||
const LG6: f64 = 1.531_383_769_920_937_332e-1;
|
||||
#[cfg(any(feature = "quantum_vml", feature = "quantum_kernel"))]
|
||||
const LG7: f64 = 1.479_819_860_511_658_591e-1;
|
||||
|
||||
/// The natural logarithm, within an ulp of the true value: `−∞` at zero, NaN
|
||||
/// below it.
|
||||
#[cfg(any(feature = "quantum_vml", feature = "quantum_kernel"))]
|
||||
pub fn ln(x: f64) -> f64 {
|
||||
let mut bits = x.to_bits();
|
||||
let mut hx = (bits >> 32) as u32;
|
||||
let mut k: i32 = 0;
|
||||
let mut x = x;
|
||||
if hx < 0x0010_0000 || hx >> 31 != 0 {
|
||||
if bits << 1 == 0 {
|
||||
return f64::NEG_INFINITY;
|
||||
}
|
||||
if hx >> 31 != 0 {
|
||||
return f64::NAN;
|
||||
}
|
||||
// subnormal: scale up
|
||||
k -= 54;
|
||||
x *= 18_014_398_509_481_984.0; // 2⁵⁴
|
||||
bits = x.to_bits();
|
||||
hx = (bits >> 32) as u32;
|
||||
} else if hx >= 0x7ff0_0000 {
|
||||
return x;
|
||||
} else if hx == 0x3ff0_0000 && bits << 32 == 0 {
|
||||
return 0.0;
|
||||
}
|
||||
// reduce x into [√2/2, √2]
|
||||
hx = hx.wrapping_add(0x3ff0_0000 - 0x3fe6_a09e);
|
||||
k += (hx >> 20) as i32 - 0x3ff;
|
||||
hx = (hx & 0x000f_ffff) + 0x3fe6_a09e;
|
||||
let x = f64::from_bits((u64::from(hx) << 32) | (bits & 0xffff_ffff));
|
||||
let f = x - 1.0;
|
||||
let hfsq = 0.5 * f * f;
|
||||
let s = f / (2.0 + f);
|
||||
let z = s * s;
|
||||
let w = z * z;
|
||||
let t1 = w * (LG2 + w * (LG4 + w * LG6));
|
||||
let t2 = z * (LG1 + w * (LG3 + w * (LG5 + w * LG7)));
|
||||
let r = t2 + t1;
|
||||
let dk = f64::from(k);
|
||||
s * (hfsq + r) + dk * LN2_LO - hfsq + f + dk * LN2_HI
|
||||
}
|
||||
|
||||
// ---------------------------------------------------------------------------
|
||||
// Arctangent
|
||||
// ---------------------------------------------------------------------------
|
||||
|
||||
#[cfg(feature = "quantum_qfhrr")]
|
||||
const ATANHI: [f64; 4] = [
|
||||
4.636_476_090_008_060_935_15e-1,
|
||||
7.853_981_633_974_482_789_99e-1,
|
||||
9.827_937_232_473_290_540_82e-1,
|
||||
1.570_796_326_794_896_558_00,
|
||||
];
|
||||
#[cfg(feature = "quantum_qfhrr")]
|
||||
const ATANLO: [f64; 4] = [
|
||||
2.269_877_745_296_168_709_24e-17,
|
||||
3.061_616_997_868_383_017_93e-17,
|
||||
1.390_331_103_123_099_845_16e-17,
|
||||
6.123_233_995_736_766_035_87e-17,
|
||||
];
|
||||
#[cfg(feature = "quantum_qfhrr")]
|
||||
const AT: [f64; 11] = [
|
||||
3.333_333_333_333_293_180_27e-1,
|
||||
-1.999_999_999_987_648_324_76e-1,
|
||||
1.428_571_427_250_346_637_11e-1,
|
||||
-1.111_111_040_546_235_578_80e-1,
|
||||
9.090_887_133_436_506_561_96e-2,
|
||||
-7.691_876_205_044_829_994_95e-2,
|
||||
6.661_073_137_387_531_206_69e-2,
|
||||
-5.833_570_133_790_573_486_45e-2,
|
||||
4.976_877_994_615_932_360_17e-2,
|
||||
-3.653_157_274_421_691_552_70e-2,
|
||||
1.628_582_011_536_578_236_23e-2,
|
||||
];
|
||||
|
||||
/// The arctangent, within an ulp of the true value.
|
||||
#[cfg(feature = "quantum_qfhrr")]
|
||||
pub fn atan(x: f64) -> f64 {
|
||||
let ix = (x.to_bits() >> 32) as u32;
|
||||
let negative = ix >> 31 != 0;
|
||||
let ix = ix & 0x7fff_ffff;
|
||||
if ix >= 0x4410_0000 {
|
||||
// |x| ≥ 2⁶⁶
|
||||
if x.is_nan() {
|
||||
return x;
|
||||
}
|
||||
let z = ATANHI[3];
|
||||
return if negative { -z } else { z };
|
||||
}
|
||||
let mut x = x;
|
||||
let id: i32;
|
||||
if ix < 0x3fdc_0000 {
|
||||
// |x| < 0.4375
|
||||
if ix < 0x3e40_0000 {
|
||||
return x;
|
||||
}
|
||||
id = -1;
|
||||
} else {
|
||||
x = x.abs();
|
||||
if ix < 0x3ff3_0000 {
|
||||
if ix < 0x3fe6_0000 {
|
||||
id = 0;
|
||||
x = (2.0 * x - 1.0) / (2.0 + x);
|
||||
} else {
|
||||
id = 1;
|
||||
x = (x - 1.0) / (x + 1.0);
|
||||
}
|
||||
} else if ix < 0x4003_8000 {
|
||||
id = 2;
|
||||
x = (x - 1.5) / (1.0 + 1.5 * x);
|
||||
} else {
|
||||
id = 3;
|
||||
x = -1.0 / x;
|
||||
}
|
||||
}
|
||||
let z = x * x;
|
||||
let w = z * z;
|
||||
let s1 = z * (AT[0] + w * (AT[2] + w * (AT[4] + w * (AT[6] + w * (AT[8] + w * AT[10])))));
|
||||
let s2 = w * (AT[1] + w * (AT[3] + w * (AT[5] + w * (AT[7] + w * AT[9]))));
|
||||
if id < 0 {
|
||||
return x - x * (s1 + s2);
|
||||
}
|
||||
let i = id as usize;
|
||||
let z = ATANHI[i] - (x * (s1 + s2) - ATANLO[i] - x);
|
||||
if negative { -z } else { z }
|
||||
}
|
||||
|
||||
#[cfg(feature = "quantum_qfhrr")]
|
||||
const PI: f64 = 3.141_592_653_589_793_116;
|
||||
#[cfg(feature = "quantum_qfhrr")]
|
||||
const PI_LO: f64 = 1.224_646_799_147_353_177_2e-16;
|
||||
|
||||
/// The angle of the point `(x, y)`, in `(−π, π]`, within an ulp of the true
|
||||
/// value.
|
||||
#[cfg(feature = "quantum_qfhrr")]
|
||||
pub fn atan2(y: f64, x: f64) -> f64 {
|
||||
if x.is_nan() || y.is_nan() {
|
||||
return x + y;
|
||||
}
|
||||
let (xb, yb) = (x.to_bits(), y.to_bits());
|
||||
let (ix, lx) = ((xb >> 32) as u32, xb as u32);
|
||||
let (iy, ly) = ((yb >> 32) as u32, yb as u32);
|
||||
if (ix.wrapping_sub(0x3ff0_0000) | lx) == 0 {
|
||||
return atan(y); // x = 1
|
||||
}
|
||||
let m = ((iy >> 31) & 1) | ((ix >> 30) & 2);
|
||||
let (ix, iy) = (ix & 0x7fff_ffff, iy & 0x7fff_ffff);
|
||||
if (iy | ly) == 0 {
|
||||
return match m {
|
||||
0 | 1 => y,
|
||||
2 => PI,
|
||||
_ => -PI,
|
||||
};
|
||||
}
|
||||
if (ix | lx) == 0 {
|
||||
return if m & 1 != 0 { -PI / 2.0 } else { PI / 2.0 };
|
||||
}
|
||||
if ix == 0x7ff0_0000 {
|
||||
return if iy == 0x7ff0_0000 {
|
||||
match m {
|
||||
0 => PI / 4.0,
|
||||
1 => -PI / 4.0,
|
||||
2 => 3.0 * PI / 4.0,
|
||||
_ => -3.0 * PI / 4.0,
|
||||
}
|
||||
} else {
|
||||
match m {
|
||||
0 => 0.0,
|
||||
1 => -0.0,
|
||||
2 => PI,
|
||||
_ => -PI,
|
||||
}
|
||||
};
|
||||
}
|
||||
if ix.wrapping_add(64 << 20) < iy || iy == 0x7ff0_0000 {
|
||||
return if m & 1 != 0 { -PI / 2.0 } else { PI / 2.0 };
|
||||
}
|
||||
let z = if m & 2 != 0 && iy.wrapping_add(64 << 20) < ix { 0.0 } else { atan((y / x).abs()) };
|
||||
match m {
|
||||
0 => z,
|
||||
1 => -z,
|
||||
2 => PI - (z - PI_LO),
|
||||
_ => (z - PI_LO) - PI,
|
||||
}
|
||||
}
|
||||
|
||||
#[cfg(test)]
|
||||
mod tests {
|
||||
#[allow(unused_imports)]
|
||||
use super::*;
|
||||
|
||||
#[allow(dead_code)]
|
||||
fn ulps(a: f64, b: f64) -> u64 {
|
||||
if a == b {
|
||||
return 0;
|
||||
}
|
||||
(a.to_bits() as i64).abs_diff(b.to_bits() as i64)
|
||||
}
|
||||
|
||||
/// A deterministic spread of inputs over many magnitudes.
|
||||
#[allow(dead_code)]
|
||||
fn inputs() -> impl Iterator<Item = f64> {
|
||||
let mut s = 0x2545_f491_4f6c_dd1du64;
|
||||
(0..200_000).map(move |i| {
|
||||
s ^= s << 13;
|
||||
s ^= s >> 7;
|
||||
s ^= s << 17;
|
||||
let u = (s >> 11) as f64 / (1u64 << 53) as f64;
|
||||
let scale = [1e-300, 1e-8, 1e-3, 0.5, 1.0, 3.0, 1e3, 1e12, 1e300][i % 9];
|
||||
u * scale
|
||||
})
|
||||
}
|
||||
|
||||
#[cfg(any(feature = "quantum_vml", feature = "quantum_kernel"))]
|
||||
#[test]
|
||||
fn ln_is_within_an_ulp_of_the_platform() {
|
||||
for x in inputs().filter(|&x| x > 0.0) {
|
||||
assert!(ulps(ln(x), x.ln()) <= 1, "ln({x:e}) = {:e} vs {:e}", ln(x), x.ln());
|
||||
}
|
||||
assert_eq!(ln(1.0), 0.0);
|
||||
assert_eq!(ln(0.0), f64::NEG_INFINITY);
|
||||
assert!(ln(-1.0).is_nan());
|
||||
assert_eq!(ln(f64::INFINITY), f64::INFINITY);
|
||||
assert!(ulps(ln(f64::MIN_POSITIVE / 8.0), (f64::MIN_POSITIVE / 8.0).ln()) <= 1, "subnormal");
|
||||
}
|
||||
|
||||
#[cfg(feature = "quantum_qfhrr")]
|
||||
#[test]
|
||||
fn atan2_is_within_an_ulp_of_the_platform() {
|
||||
let xs: Vec<f64> = inputs().take(4000).collect();
|
||||
for (i, &a) in xs.iter().enumerate() {
|
||||
let b = xs[(i * 7 + 3) % xs.len()];
|
||||
for (y, x) in [(a, b), (-a, b), (a, -b), (-a, -b)] {
|
||||
assert!(ulps(atan2(y, x), y.atan2(x)) <= 1, "atan2({y:e}, {x:e}) = {:e} vs {:e}", atan2(y, x), y.atan2(x));
|
||||
}
|
||||
}
|
||||
assert_eq!(atan2(0.0, -1.0), PI);
|
||||
assert_eq!(atan2(1.0, 0.0), PI / 2.0);
|
||||
assert_eq!(atan(f64::INFINITY), ATANHI[3]);
|
||||
}
|
||||
}
|
||||
109
wai/wai-quantum/tests/bit_identity.rs
Normal file
109
wai/wai-quantum/tests/bit_identity.rs
Normal file
|
|
@ -0,0 +1,109 @@
|
|||
//! The crate promises identical results natively, in the browser and at the
|
||||
//! edge. Before every transcendental moved in-crate, four of these layers gave
|
||||
//! different bits natively and under `wasm32-wasip2`. These digests were taken
|
||||
//! natively and confirmed under `wasm32-wasip2`; any platform whose math
|
||||
//! differs fails here, not in a user's diff.
|
||||
//!
|
||||
//! `CARGO_TARGET_WASM32_WASIP2_RUNNER=wasmtime cargo test --release --target wasm32-wasip2 --features full`
|
||||
#![cfg(all(
|
||||
feature = "quantum",
|
||||
feature = "quantum_vml",
|
||||
feature = "quantum_phasor",
|
||||
feature = "quantum_qfhrr",
|
||||
feature = "quantum_kernel",
|
||||
feature = "quantum_mps",
|
||||
feature = "quantum_pauli",
|
||||
feature = "quantum_spd"
|
||||
))]
|
||||
|
||||
use wai_quantum::*;
|
||||
|
||||
/// FNV-1a over every value's IEEE-754 bits.
|
||||
fn digest(xs: &[f64]) -> u64 {
|
||||
let mut h: u64 = 0xcbf2_9ce4_8422_2325;
|
||||
for x in xs {
|
||||
for b in x.to_bits().to_le_bytes() {
|
||||
h = (h ^ u64::from(b)).wrapping_mul(0x100_0000_01b3);
|
||||
}
|
||||
}
|
||||
h
|
||||
}
|
||||
|
||||
fn trained() -> quantum_vml::BornFit {
|
||||
quantum_vml::train_born(&quantum_vml::parity_target(3), 3, 2, 40, 0.1, 7)
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn born_machine_training() {
|
||||
let fit = trained();
|
||||
let mut v = fit.params.clone();
|
||||
v.extend(&fit.loss_curve);
|
||||
v.extend(&fit.final_probs);
|
||||
v.push(fit.final_kl);
|
||||
assert_eq!(digest(&v), 0xd8e5_bd02_2201_2e5d);
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn phasor_bridge() {
|
||||
let p = quantum_phasor::random_params(3, 5);
|
||||
let (a, b, c, d) = quantum_phasor::compare(3, &p, 64);
|
||||
let mut v = a;
|
||||
v.extend(b);
|
||||
v.extend(c);
|
||||
v.push(d);
|
||||
assert_eq!(digest(&v), 0x08be_30dd_7ce1_9561);
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn qfhrr_born_machine() {
|
||||
let fit = trained();
|
||||
let mut v = quantum_qfhrr::born_probs(3, 2, &fit.params, 4);
|
||||
v.push(quantum_qfhrr::tv_vs_exact(3, 2, &fit.params, 4));
|
||||
assert_eq!(digest(&v), 0x5407_7c4d_1be0_d985);
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn phasor_kernel() {
|
||||
let (x, y) = quantum_kernel::circles(200, 3);
|
||||
let rff = quantum_kernel::sample_rff(2, 64, 1.0, 9);
|
||||
let w = quantum_kernel::fit_ridge(&rff, &x, &y, 1e-3);
|
||||
let mut v: Vec<f64> = x.iter().flatten().copied().collect();
|
||||
v.extend(&w);
|
||||
assert_eq!(digest(&v), 0x58ba_1666_5ef6_fd9e);
|
||||
}
|
||||
|
||||
fn brickwork() -> quantum::Circuit {
|
||||
let mut c = quantum::Circuit::new(8);
|
||||
for l in 0..4u8 {
|
||||
for q in 0..8 {
|
||||
c.h(q);
|
||||
c.t(q);
|
||||
}
|
||||
for q in (l % 2..7).step_by(2) {
|
||||
c.cx(q, q + 1);
|
||||
}
|
||||
}
|
||||
c
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn mps_truncation() {
|
||||
let c = brickwork();
|
||||
let (mps, _) = quantum_mps::run(&c, 4).unwrap();
|
||||
assert_eq!(digest(&[quantum_mps::fidelity(&mps, &c.simulate().unwrap())]), 0xa807_44ca_a97f_bff8);
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn sparse_pauli_propagation() {
|
||||
let (e, _) = quantum_pauli::expect_z(&brickwork(), 3, 1e-6).unwrap();
|
||||
assert_eq!(digest(&[e]), 0x8cc6_5829_1fb1_8909);
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn sparse_pauli_dynamics() {
|
||||
let chain: Vec<(u32, u32)> = (0..29).map(|q| (q, q + 1)).collect();
|
||||
let circuit = quantum_spd::kicked_ising(30, &chain, 0.5, 6);
|
||||
let cfg = quantum_spd::SpdConfig { threshold: 1e-6, max_weight: None };
|
||||
let r = quantum_spd::expect_z(&circuit, 7, &cfg).unwrap();
|
||||
assert_eq!(digest(&[r.value, r.truncation_bound]), 0xfc44_f4b6_0e70_68b2);
|
||||
}
|
||||
37
wai/wai-quantum/tests/utility_reference.rs
Normal file
37
wai/wai-quantum/tests/utility_reference.rs
Normal file
|
|
@ -0,0 +1,37 @@
|
|||
//! The 127-qubit heavy-hex kicked-Ising experiment against published values.
|
||||
//!
|
||||
//! - The 5-step magnetization `Mz = (1/127) Σ ⟨Z_q⟩` at `θ_h = kπ/32` was
|
||||
//! computed exactly and published as `exact.csv` with arXiv:2308.05077.
|
||||
//! With the truncation threshold at 0 this crate reproduces `k = 1` to
|
||||
//! 3.5e-7 (0.995131 against 0.995131, 357M terms); here a threshold of
|
||||
//! 1e-7 keeps the test fast, and the error is checked against the bound.
|
||||
//! - ⟨Z₆₂⟩ after 20 steps at `θ_h = 2π/32`: the converged classical values
|
||||
//! published with the same paper (sparse Pauli dynamics, a mixed method,
|
||||
//! PEPS, PEPO) span 0.9776–0.9796.
|
||||
#![cfg(feature = "quantum_spd")]
|
||||
|
||||
use wai_quantum::quantum_spd::{Pauli, SpdConfig, expect_z, expectation, heavy_hex_127, kicked_ising};
|
||||
|
||||
fn theta(k: u32) -> f64 {
|
||||
f64::from(k) * std::f64::consts::PI / 32.0
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn five_step_magnetization_matches_the_exact_values() {
|
||||
let mz: Vec<(f64, Vec<(u32, Pauli)>)> = (0..127).map(|q| (1.0 / 127.0, vec![(q, Pauli::Z)])).collect();
|
||||
for (k, exact) in [(1, 0.995131), (2, 0.979836)] {
|
||||
let c = kicked_ising(127, &heavy_hex_127(), theta(k), 5);
|
||||
let r = expectation(&c, &mz, &[], &SpdConfig { threshold: 1e-7, max_weight: None }).unwrap();
|
||||
let err = (r.value - exact).abs();
|
||||
// The published values carry six decimals: allow their rounding.
|
||||
assert!(err <= r.truncation_bound + 5e-7, "k = {k}: error {err:e} exceeds the bound {:e}", r.truncation_bound);
|
||||
assert!(err < 3e-5, "k = {k}: Mz = {} against {exact}", r.value);
|
||||
}
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn twenty_step_z62_lies_among_the_converged_classical_values() {
|
||||
let c = kicked_ising(127, &heavy_hex_127(), theta(2), 20);
|
||||
let r = expect_z(&c, 62, &SpdConfig { threshold: 1e-4, max_weight: None }).unwrap();
|
||||
assert!((0.9776..=0.9796).contains(&r.value), "⟨Z₆₂⟩ = {}", r.value);
|
||||
}
|
||||
2
wai/wai-rs/Cargo.lock
generated
2
wai/wai-rs/Cargo.lock
generated
|
|
@ -3900,7 +3900,7 @@ dependencies = [
|
|||
|
||||
[[package]]
|
||||
name = "wai-quantum"
|
||||
version = "0.3.21"
|
||||
version = "0.3.22"
|
||||
dependencies = [
|
||||
"blake3",
|
||||
"ed25519-dalek",
|
||||
|
|
|
|||
|
|
@ -335,6 +335,10 @@ quantum_stabilizer = ["quantum", "wai-quantum/quantum_stabilizer"]
|
|||
# Pauli dynamics). Computes ⟨O⟩ on large
|
||||
# low-magic circuits where no state vector exists; grows with magic, not qubits.
|
||||
quantum_pauli = ["quantum", "wai-quantum/quantum_pauli"]
|
||||
# Sparse Pauli dynamics at utility scale (arbitrary angles, up to 1024 qubits),
|
||||
# and its signed receipts.
|
||||
quantum_spd = ["quantum", "wai-quantum/quantum_spd"]
|
||||
quantum_spd_receipt = ["quantum_spd", "wai-quantum/quantum_spd_receipt"]
|
||||
quantum_qir = ["quantum_toolchain", "wai-quantum/quantum_qir"]
|
||||
# Variational quantum machine learning (extensions/quantum-sim § VML,
|
||||
# `wai.quantum.vml`): the interference-based ML the thermodynamic-computing pivot
|
||||
|
|
@ -460,6 +464,7 @@ required-features = [
|
|||
"quantum_channel",
|
||||
"quantum_comm",
|
||||
"quantum_qudit",
|
||||
"quantum_spd_receipt",
|
||||
]
|
||||
|
||||
[[bin]]
|
||||
|
|
|
|||
|
|
@ -27,6 +27,14 @@
|
|||
//! wai-quantum 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 utility --theta 0.2 --steps 20 # the 127-qubit heavy-hex experiment
|
||||
//! wai-quantum utility --receipt | wai-quantum verify
|
||||
//! ```
|
||||
//!
|
||||
//! 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
|
||||
|
|
@ -42,6 +50,7 @@ use wai::quantum_info::{chsh, DensityMatrix};
|
|||
use wai::quantum_qudit::QuditState;
|
||||
use wai::quantum_energy::{EnergyClass, MaybeLabelled};
|
||||
use wai::quantum_receipt::QuantumReceipt;
|
||||
use wai::quantum_spd_receipt::{SpdJob, SpdReceipt};
|
||||
use wai::quantum_source::{holevo_bound, schumacher_limit, von_neumann_entropy};
|
||||
use wai::quantum_toolchain as qt;
|
||||
|
||||
|
|
@ -115,6 +124,8 @@ fn usage() -> ! {
|
|||
wai-quantum bb84 [--rounds N] [--seed S] [--eavesdrop] [--json]\n \
|
||||
wai-quantum chsh [--json]\n \
|
||||
wai-quantum qutrit [-n QUTRITS] [--json]\n\n\
|
||||
UTILITY SCALE:\n \
|
||||
wai-quantum utility [--theta T] [--steps N] [--threshold D] [-q QUBIT] [--json] [--receipt]\n\n\
|
||||
Channels: identity, bitflip, dephasing, depolarizing, damping.\n\
|
||||
Run `wai-quantum list` for the algorithm catalog."
|
||||
);
|
||||
|
|
@ -257,6 +268,24 @@ fn main() {
|
|||
// carries its acquisition class inside the signature
|
||||
// (`quantum_energy::Labelled`, as `wai_quantum_meter` emits).
|
||||
let src = read_stdin();
|
||||
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 {} uJ", 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)");
|
||||
|
|
@ -285,6 +314,49 @@ fn main() {
|
|||
// 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_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);
|
||||
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") {
|
||||
// The same deterministic demo key as `receipt`: the run is
|
||||
// bit-identical everywhere, so is its receipt.
|
||||
let key = SigningKey::from_bytes(&[42u8; 32]);
|
||||
let r = SpdReceipt::seal(&key, "did:key:wai-quantum-cli", &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);
|
||||
}
|
||||
}
|
||||
|
||||
"entropy" => {
|
||||
let c = circuit_from(args, algo);
|
||||
let sv = c.simulate().unwrap_or_else(|e| die(&format!("simulate: {e:?}")));
|
||||
|
|
|
|||
|
|
@ -581,6 +581,12 @@ pub use wai_quantum::quantum_stabilizer;
|
|||
// large low-magic circuits where no 2^n state vector exists. Verified vs dense.
|
||||
#[cfg(feature = "quantum_pauli")]
|
||||
pub use wai_quantum::quantum_pauli;
|
||||
|
||||
#[cfg(feature = "quantum_spd")]
|
||||
pub use wai_quantum::quantum_spd;
|
||||
|
||||
#[cfg(feature = "quantum_spd_receipt")]
|
||||
pub use wai_quantum::quantum_spd_receipt;
|
||||
#[cfg(feature = "quantum_qir")]
|
||||
pub use wai_quantum::quantum_qir;
|
||||
|
||||
|
|
|
|||
2
wai/wai-web/Cargo.lock
generated
2
wai/wai-web/Cargo.lock
generated
|
|
@ -734,7 +734,7 @@ dependencies = [
|
|||
|
||||
[[package]]
|
||||
name = "wai-quantum"
|
||||
version = "0.3.21"
|
||||
version = "0.3.22"
|
||||
dependencies = [
|
||||
"blake3",
|
||||
"ed25519-dalek",
|
||||
|
|
|
|||
Loading…
Reference in a new issue