open-standards/wai/wai-quantum/examples/gbs_loops/common.rs

//! What the reproductions of the 2022 time-domain Gaussian boson sampler
//! share (Madsen et al., Nature 606, 75 (2022), doi:10.1038/s41586-022-04725-x):
//! the pinned data set, a SHA-256 check of every file read, a reader of the
//! `.npy` arrays it is published in, the instance's loop program, the
//! photon-number patterns in a fixed order, the energy reading around a run,
//! and the signed receipt each run writes.
//!
//! The data set is published by the paper's authors under the Apache License,
//! Version 2.0 (see `tests/data/gbs/README.md` for the credit and the licence
//! copy). Its arrays are fetched by the reader and checked against the hashes
//! here; only the small fixtures the tests need are kept in this repository.

#![allow(dead_code, clippy::needless_range_loop)]

use std::path::Path;
use std::time::Instant;

use serde_json::{Value, json};
use wai_quantum::ed25519_dalek::SigningKey;
use wai_quantum::quantum_gbs::{C, Loops, Mat};

/// The repository that describes the data set, at the commit read, and the
/// bucket that serves its arrays.
pub const REPOSITORY: &str = "https://github.com/XanaduAI/xanadu-qca-data"; // allow-name: a repository URL
pub const COMMIT: &str = "869906e3f4ec2a5c5fd6a78d5d8c8a20371b9e6c";
pub const BUCKET: &str = "https://qca-data.s3.amazonaws.com/"; // allow-name: the data set's URL
/// The repository's `LICENSE` at [`COMMIT`]: the Apache License, Version 2.0.
pub const LICENSE_SHA256: &str = "c71d239df91726fc519c6eb72d318ec65820627232b2f796219e87dcf35d0ab4";

/// The five instances: the paper's figure each is shown in, and its modes.
pub const SETS: [(&str, usize); 5] = [("fig2", 16), ("fig3a", 216), ("fig3b", 72), ("fig4", 216), ("figS15", 288)];

/// Every object of the data set, by its path in the bucket: its size in bytes
/// and the SHA-256 of its bytes (the MD5 of each single-part object is the
/// bucket's ETag, and agreed when fetched).
pub const SOURCES: &[(&str, u64, &str)] = &[
    ("fig2/program_params.json", 4454, "17c2d13824e61307a367f111a6dd4faa5c8f43d249effa00d5f2654d3ce8c0bc"),
    ("fig2/U.npy", 4224, "18685d0b53da451ae4cbb91128ab38ef86e88e6bb5db26e01d19496cb1b8bb31"),
    ("fig2/T.npy", 4224, "bd2b1805bf3a07734dbb0ebd6bedfe7e60a8fef04c0712f79d72451c3d070e00"),
    ("fig2/r.npy", 256, "e5bea96c7a25b9b30dbf1180bacf7165e0c1f92c7f00ebd9d06bf39d245f089b"),
    ("fig2/samples.npy", 1345536128, "a24c3500857c89ad7172dfd5eee7d3d961d6d8ec15169b9f11fb7eaf1da9ecac"),
    ("fig3a/program_params.json", 45621, "ed584aa6539d8fdca56b5dabc069cfacb809bd1f03d491f8307259f02b224eee"),
    ("fig3a/U.npy", 746624, "58ebeb2e1769f9c768eec7672ec322877e43e2c20e2077c8658b1bcfccd4c178"),
    ("fig3a/T.npy", 746624, "0f8d001df2a1d86c31b2bea5289b1f7fab9e9df9ec9127de268d2eccc3943091"),
    ("fig3a/r.npy", 1856, "599c27dff53d52f831dd99816481c8a2d231a37e2d0f3132c8e77fdfaee1aef0"),
    ("fig3a/samples.npy", 225504128, "50b4b3bb37a0e4e2162ff8d188d6386f4dcd0c5b3395828b016f43c19f230725"),
    ("fig3b/program_params.json", 21265, "95d723496d328a7116f6b999a1e99e90c3c0d6b6174fca78e082222fbbd6e643"),
    ("fig3b/U.npy", 83072, "f7aa1d75e0e5945eadc123b74cd5b5d92ed6e1c01f61cb6354f01c6f738ee393"),
    ("fig3b/T.npy", 83072, "c6b796dc2478bc2786ca933255786665b3fd41e73f631ca8269f2b21115ffc54"),
    ("fig3b/r.npy", 1856, "b5e5b36bc0bc576f40c9c313afdcc0004a469d8d5acead1228795573d2d71327"),
    ("fig3b/samples.npy", 74592128, "37691c55f58549c353d0749aba71bca07bbcb0f322c1355e04b5acf72766ed65"),
    ("fig4/program_params.json", 45587, "9d7a9b02cf934ea0913f48a2f68a7285967c555af15b7de9246bec8aedae7070"),
    ("fig4/U.npy", 746624, "3654a12fbdb707b92eeda5bafc58f6ef9ca03b8fbb3e59b0ac316870f0c2fd00"),
    ("fig4/T.npy", 746624, "c7e5c10b63002ff6e270f9e049aebf36371f62dfc38e73c91725c10c2d0cddbf"),
    ("fig4/r.npy", 1856, "f593db9ace468d06fa556f5a074b583756b9e6b4c5264a23e6edcb45fd48c40e"),
    ("fig4/samples.npy", 221616128, "20a3a5a334f4d40b6a9de9abb46df3eb5390ce6f660eb15e1a630f406c88165a"),
    ("figS15/program_params.json", 59802, "6900e557752c9608d7582ec3920e68faef4d57edc997f3ba999bed3db0acb085"),
    ("figS15/U.npy", 1327232, "88338b73614f3af7034498d45c1612c6875c75e31e75fa6481d2f2409256b15a"),
    ("figS15/T.npy", 1327232, "fae02f0cd1d26abf6b152ab80386ea2799809ce592194ad1b48309038d1f4c28"),
    ("figS15/r.npy", 2432, "9309ea8d12dd68deb7fa677e039efbe52c14ed653633b10fac2569850a985905"),
    ("figS15/samples.npy", 320544128, "94e7f0c8ac2280a6cc26028a0ab75afc6bdec2b620cfb626def08a78adb51130"),
];

/// The pinned size and SHA-256 of `path`.
pub fn pinned(path: &str) -> (u64, &'static str) {
    SOURCES.iter().find(|(p, _, _)| *p == path).map(|(_, n, h)| (*n, *h)).unwrap_or_else(|| panic!("{path} is not a pinned object"))
}

/// The bytes of `path` under the data directory `root` (laid out as the
/// bucket is: `root/fig2/T.npy`), refused unless their size and SHA-256 are
/// the pinned ones.
pub fn read_pinned(root: &Path, path: &str) -> Vec<u8> {
    let bytes = std::fs::read(root.join(path)).unwrap_or_else(|e| panic!("reading {path} under {}: {e}", root.display()));
    let (size, sha) = pinned(path);
    assert_eq!(bytes.len() as u64, size, "{path}: {} bytes, not the pinned {size}", bytes.len());
    let got = hex(&sha256(&bytes));
    assert_eq!(got, sha, "{path}: SHA-256 {got} is not the pinned one; fetch it again from {BUCKET}");
    bytes
}

/// What a receipt names as an input: where it is served, its size and hash.
pub fn source(path: &str) -> Value {
    let (bytes, sha256) = pinned(path);
    json!({ "url": format!("{BUCKET}{path}"), "bytes": bytes, "sha256": sha256 })
}

pub fn hex(b: &[u8]) -> String {
    b.iter().map(|x| format!("{x:02x}")).collect()
}

/// SHA-256 (FIPS 180-4), over the bytes in place.
pub fn sha256(data: &[u8]) -> [u8; 32] {
    const K: [u32; 64] = [
        0x428a2f98, 0x71374491, 0xb5c0fbcf, 0xe9b5dba5, 0x3956c25b, 0x59f111f1, 0x923f82a4, 0xab1c5ed5, 0xd807aa98, 0x12835b01, 0x243185be,
        0x550c7dc3, 0x72be5d74, 0x80deb1fe, 0x9bdc06a7, 0xc19bf174, 0xe49b69c1, 0xefbe4786, 0x0fc19dc6, 0x240ca1cc, 0x2de92c6f, 0x4a7484aa,
        0x5cb0a9dc, 0x76f988da, 0x983e5152, 0xa831c66d, 0xb00327c8, 0xbf597fc7, 0xc6e00bf3, 0xd5a79147, 0x06ca6351, 0x14292967, 0x27b70a85,
        0x2e1b2138, 0x4d2c6dfc, 0x53380d13, 0x650a7354, 0x766a0abb, 0x81c2c92e, 0x92722c85, 0xa2bfe8a1, 0xa81a664b, 0xc24b8b70, 0xc76c51a3,
        0xd192e819, 0xd6990624, 0xf40e3585, 0x106aa070, 0x19a4c116, 0x1e376c08, 0x2748774c, 0x34b0bcb5, 0x391c0cb3, 0x4ed8aa4a, 0x5b9cca4f,
        0x682e6ff3, 0x748f82ee, 0x78a5636f, 0x84c87814, 0x8cc70208, 0x90befffa, 0xa4506ceb, 0xbef9a3f7, 0xc67178f2,
    ];
    let mut h: [u32; 8] = [0x6a09e667, 0xbb67ae85, 0x3c6ef372, 0xa54ff53a, 0x510e527f, 0x9b05688c, 0x1f83d9ab, 0x5be0cd19];
    let block = |h: &mut [u32; 8], block: &[u8]| {
        let mut w = [0u32; 64];
        for i in 0..16 {
            w[i] = u32::from_be_bytes([block[4 * i], block[4 * i + 1], block[4 * i + 2], block[4 * i + 3]]);
        }
        for i in 16..64 {
            let s0 = w[i - 15].rotate_right(7) ^ w[i - 15].rotate_right(18) ^ (w[i - 15] >> 3);
            let s1 = w[i - 2].rotate_right(17) ^ w[i - 2].rotate_right(19) ^ (w[i - 2] >> 10);
            w[i] = w[i - 16].wrapping_add(s0).wrapping_add(w[i - 7]).wrapping_add(s1);
        }
        let [mut a, mut b, mut c, mut d, mut e, mut f, mut g, mut hh] = *h;
        for i in 0..64 {
            let s1 = e.rotate_right(6) ^ e.rotate_right(11) ^ e.rotate_right(25);
            let ch = (e & f) ^ (!e & g);
            let t1 = hh.wrapping_add(s1).wrapping_add(ch).wrapping_add(K[i]).wrapping_add(w[i]);
            let s0 = a.rotate_right(2) ^ a.rotate_right(13) ^ a.rotate_right(22);
            let maj = (a & b) ^ (a & c) ^ (b & c);
            let t2 = s0.wrapping_add(maj);
            hh = g;
            g = f;
            f = e;
            e = d.wrapping_add(t1);
            d = c;
            c = b;
            b = a;
            a = t1.wrapping_add(t2);
        }
        for (x, y) in h.iter_mut().zip([a, b, c, d, e, f, g, hh]) {
            *x = x.wrapping_add(y);
        }
    };
    let whole = data.len() / 64 * 64;
    for chunk in data[..whole].chunks(64) {
        block(&mut h, chunk);
    }
    let mut tail = data[whole..].to_vec();
    tail.push(0x80);
    while tail.len() % 64 != 56 {
        tail.push(0);
    }
    tail.extend_from_slice(&((data.len() as u64) * 8).to_be_bytes());
    for chunk in tail.chunks(64) {
        block(&mut h, chunk);
    }
    let mut out = [0u8; 32];
    for (i, x) in h.iter().enumerate() {
        out[4 * i..4 * i + 4].copy_from_slice(&x.to_be_bytes());
    }
    out
}

// ---------------------------------------------------------------------------
// .npy arrays
// ---------------------------------------------------------------------------

/// A `.npy` array: its element type, memory order, shape and data bytes.
pub struct Npy<'a> {
    pub descr: String,
    pub fortran: bool,
    pub shape: Vec<usize>,
    pub data: &'a [u8],
}

/// The array in `bytes`; the header is a Python dict literal, read for its
/// three keys only (nothing in it is evaluated).
pub fn npy(bytes: &[u8]) -> Npy<'_> {
    assert!(bytes.len() >= 10 && &bytes[..6] == b"\x93NUMPY", "not an .npy array");
    let (len, start) = match bytes[6] {
        1 => (u16::from_le_bytes([bytes[8], bytes[9]]) as usize, 10),
        2 | 3 => (u32::from_le_bytes([bytes[8], bytes[9], bytes[10], bytes[11]]) as usize, 12),
        v => panic!(".npy version {v}"),
    };
    let header = std::str::from_utf8(&bytes[start..start + len]).expect("an ASCII header");
    let value = |key: &str| -> &str {
        let at = header.find(&format!("'{key}':")).unwrap_or_else(|| panic!("no {key} in the header"));
        header[at + key.len() + 3..].trim_start()
    };
    let descr = value("descr").trim_start_matches('\'').split('\'').next().unwrap().to_string();
    let fortran = value("fortran_order").starts_with("True");
    let shape_text = value("shape");
    let shape = shape_text[1..shape_text.find(')').unwrap()].split(',').filter_map(|s| s.trim().parse().ok()).collect::<Vec<usize>>();
    let count: usize = shape.iter().product();
    let width = match descr.as_str() {
        "<c16" => 16,
        "<f8" => 8,
        "|i1" | "|u1" => 1,
        d => panic!("element type {d}"),
    };
    let data = &bytes[start + len..];
    assert_eq!(data.len(), count * width, "the data does not fill the shape");
    Npy { descr, fortran, shape, data }
}

fn f64_at(b: &[u8], i: usize) -> f64 {
    f64::from_le_bytes(b[8 * i..8 * i + 8].try_into().unwrap())
}

/// A square complex matrix (`<c16`), row-major whatever the file's order.
pub fn complex_matrix(bytes: &[u8]) -> Mat {
    let a = npy(bytes);
    assert!(a.descr == "<c16" && a.shape.len() == 2 && a.shape[0] == a.shape[1], "a square complex matrix");
    let n = a.shape[0];
    let mut m = Mat::zeros(n);
    for i in 0..n {
        for j in 0..n {
            let k = if a.fortran { j * n + i } else { i * n + j };
            m.a[i * n + j] = C { re: f64_at(a.data, 2 * k), im: f64_at(a.data, 2 * k + 1) };
        }
    }
    m
}

/// A real vector (`<f8`).
pub fn real_vector(bytes: &[u8]) -> Vec<f64> {
    let a = npy(bytes);
    assert!(a.descr == "<f8" && a.shape.len() == 1, "a real vector");
    (0..a.shape[0]).map(|i| f64_at(a.data, i)).collect()
}

/// The samples (`|i1`, shape `shots × 1 × modes`) as one row of `modes` photon
/// numbers per shot, in shot order, whatever the file's memory order.
pub fn sample_rows(bytes: &[u8]) -> (usize, usize, std::borrow::Cow<'_, [u8]>) {
    let a = npy(bytes);
    assert!(a.descr == "|i1" && a.shape.len() == 3 && a.shape[1] == 1, "photon-number samples");
    let (shots, modes) = (a.shape[0], a.shape[2]);
    if !a.fortran {
        return (shots, modes, std::borrow::Cow::Borrowed(a.data));
    }
    // Column-major: mode m of shot s at s + shots·m. Transposed in blocks of shots.
    let mut rows = vec![0u8; shots * modes];
    const BLOCK: usize = 4096;
    for s0 in (0..shots).step_by(BLOCK) {
        let s1 = (s0 + BLOCK).min(shots);
        for m in 0..modes {
            let col = &a.data[m * shots + s0..m * shots + s1];
            for (k, &v) in col.iter().enumerate() {
                rows[(s0 + k) * modes + m] = v;
            }
        }
    }
    (shots, modes, std::borrow::Cow::Owned(rows))
}

// ---------------------------------------------------------------------------
// The instance's program
// ---------------------------------------------------------------------------

/// An instance's program, as `program_params.json` gives it.
pub struct Program {
    /// Squeezing per input bin, as listed.
    pub r: Vec<f64>,
    pub loops: Loops,
    /// The energy transmittance from the source to the detectors, common to every mode.
    pub common: f64,
    /// The 16 detector channels' relative efficiencies, bin `b` read by channel `b mod 16`.
    pub channels: Vec<f64>,
    /// Bins read before the first pulse arrives, dropped from the samples and matrices.
    pub crop: usize,
}

pub fn program(bytes: &[u8]) -> Program {
    let v: Value = serde_json::from_slice(bytes).expect("the program is JSON");
    let list = |k: &str| -> Vec<f64> { v[k].as_array().unwrap_or_else(|| panic!("{k}")).iter().map(|x| x.as_f64().expect("a number")).collect() };
    Program {
        r: list("r"),
        loops: Loops {
            delays: vec![1, 6, 36],
            phi: vec![list("phi_0"), list("phi_1"), list("phi_2")],
            alpha: vec![list("alpha_0"), list("alpha_1"), list("alpha_2")],
            round_trip_phase: list("loop_phases"),
            round_trip_transmission: list("loop_efficiencies"),
        },
        common: v["common_efficiency"].as_f64().expect("common_efficiency"),
        channels: list("relative_mode_efficiencies"),
        crop: v["crop"].as_u64().expect("crop") as usize,
    }
}

impl Program {
    /// The transfer matrix of the instance's `m` modes: the pulses of input
    /// bins `0..m` to the detected bins `crop..crop + m`, each detected bin `b`
    /// through the common transmittance and channel `b mod 16`'s
    /// (Methods, eq. 10). With `lossless`, the loops' and detectors' losses are
    /// left out: the instance's unitary.
    pub fn transfer(&self, m: usize, lossless: bool) -> Mat {
        let detection: Vec<f64> = (0..m).map(|i| if lossless { 1.0 } else { self.common * self.channels[(i + self.crop) % 16] }).collect();
        let loops = if lossless { Loops { round_trip_transmission: vec![1.0; 3], ..self.loops.clone() } } else { self.loops.clone() };
        loops.transfer(m, self.crop, &detection).expect("the program covers the instance")
    }
}

/// The largest `|a − φ b|` over the entries, for the unit phase `φ ∈ {1, i, −1, −i}`
/// that makes it smallest, and that phase as a quarter-turn count: a global
/// phase changes no photon-number statistic.
pub fn distance_up_to_quarter_phase(a: &Mat, b: &Mat) -> (f64, u32) {
    let phases = [C { re: 1.0, im: 0.0 }, C { re: 0.0, im: 1.0 }, C { re: -1.0, im: 0.0 }, C { re: 0.0, im: -1.0 }];
    let mut best = (f64::INFINITY, 0);
    for (k, p) in phases.iter().enumerate() {
        let d = a
            .a
            .iter()
            .zip(&b.a)
            .map(|(x, y)| {
                let (re, im) = (p.re * y.re - p.im * y.im, p.re * y.im + p.im * y.re);
                ((x.re - re) * (x.re - re) + (x.im - im) * (x.im - im)).sqrt()
            })
            .fold(0.0, f64::max);
        if d < best.0 {
            best = (d, k as u32);
        }
    }
    best
}

// ---------------------------------------------------------------------------
// Patterns
// ---------------------------------------------------------------------------

/// Every photon-number pattern of `modes` modes with `total` photons, in a
/// fixed order: the first mode's count from `total` down to 0, and within it
/// the rest in the same order.
pub fn patterns(modes: usize, total: usize) -> Vec<Vec<u8>> {
    fn go(modes: usize, total: usize, prefix: &mut Vec<u8>, out: &mut Vec<Vec<u8>>) {
        if modes == 1 {
            prefix.push(total as u8);
            out.push(prefix.clone());
            prefix.pop();
            return;
        }
        for k in (0..=total).rev() {
            prefix.push(k as u8);
            go(modes - 1, total - k, prefix, out);
            prefix.pop();
        }
    }
    let mut out = Vec::new();
    go(modes, total, &mut Vec::with_capacity(modes), &mut out);
    out
}

/// A pattern of at most 16 modes with at most 15 photons in any mode, packed
/// four bits a mode.
pub fn pack(p: &[u8]) -> u64 {
    p.iter().enumerate().fold(0u64, |acc, (i, &c)| acc | (u64::from(c) << (4 * i)))
}

/// The pattern counts file: one line per total photon number `N`, `N` then
/// the count of each pattern in [`patterns`]' order.
pub fn counts_text(modes: usize, nmax: usize, count: impl Fn(&[u8]) -> u64) -> String {
    let mut s = String::new();
    for n in 0..=nmax {
        s += &n.to_string();
        for p in patterns(modes, n) {
            s += &format!(" {}", count(&p));
        }
        s.push('\n');
    }
    s
}

/// The counts file read back: per `N`, the counts in [`patterns`]' order.
pub fn parse_counts(text: &str) -> Vec<Vec<u64>> {
    text.lines()
        .enumerate()
        .map(|(n, line)| {
            let mut it = line.split(' ');
            assert_eq!(it.next().and_then(|x| x.parse::<usize>().ok()), Some(n), "line {n} starts with its N");
            it.map(|x| x.parse().expect("a count")).collect()
        })
        .collect()
}

/// The fidelity `Σ √(pᵢ qᵢ)` and the total variation distance `Σ |pᵢ − qᵢ| / 2`
/// of two distributions over the same outcomes.
pub fn fidelity_tvd(p: &[f64], q: &[f64]) -> (f64, f64) {
    let (mut f, mut t) = (0.0, 0.0);
    for (a, b) in p.iter().zip(q) {
        f += (a * b).sqrt();
        t += (a - b).abs();
    }
    (f, t / 2.0)
}

/// A deterministic stream of uniform doubles in (0, 1) (SplitMix64).
pub struct Rng(pub u64);

impl Rng {
    pub fn next(&mut self) -> f64 {
        self.0 = self.0.wrapping_add(0x9e37_79b9_7f4a_7c15);
        let mut z = self.0;
        z = (z ^ (z >> 30)).wrapping_mul(0xbf58_476d_1ce4_e5b9);
        z = (z ^ (z >> 27)).wrapping_mul(0x94d0_49bb_1331_11eb);
        (((z ^ (z >> 31)) >> 11) as f64 + 0.5) / 9_007_199_254_740_992.0
    }
}

/// `shots` draws from the distribution `p` (normalised), as counts.
pub fn draw(p: &[f64], shots: u64, rng: &mut Rng) -> Vec<u64> {
    let mut cum = Vec::with_capacity(p.len());
    let mut s = 0.0;
    for &x in p {
        s += x;
        cum.push(s);
    }
    let mut counts = vec![0u64; p.len()];
    for _ in 0..shots {
        let u = rng.next() * s;
        let k = cum.partition_point(|&c| c <= u).min(p.len() - 1);
        counts[k] += 1;
    }
    counts
}

// ---------------------------------------------------------------------------
// Energy, time and receipts
// ---------------------------------------------------------------------------

/// Run `f`, timing it and reading an energy counter around it: the result, the
/// wall time in seconds, and the energy figure (or why there is none). On macOS
/// the counter is this process's own CPU energy; on Linux, the processor
/// packages' counters, which count every process on the host (the figure says
/// so, and gives the host's load when the run began).
pub fn metered<T>(f: impl FnOnce() -> T) -> (T, f64, Value) {
    if cfg!(target_os = "linux") {
        return metered_packages(f);
    }
    metered_process(f)
}

fn metered_packages<T>(f: impl FnOnce() -> T) -> (T, f64, Value) {
    use std::sync::Mutex;
    use std::sync::atomic::{AtomicBool, Ordering};
    use wai_quantum::engine::energy::{Accumulated, Packages, package_figure};
    let load = std::fs::read_to_string("/proc/loadavg").ok().and_then(|s| s.split_whitespace().next().map(str::to_string));
    let packages = match Packages::open() {
        Ok(p) => p,
        Err(why) => {
            let started = Instant::now();
            let out = f();
            return (out, started.elapsed().as_secs_f64(), json!({ "unmetered": why }));
        }
    };
    let Some(first) = packages.read() else {
        let started = Instant::now();
        let out = f();
        return (out, started.elapsed().as_secs_f64(), json!({ "unmetered": "a package counter did not read" }));
    };
    let acc = Mutex::new(Accumulated::new(&packages, first));
    let done = AtomicBool::new(false);
    let started = Instant::now();
    let (out, ok) = std::thread::scope(|s| {
        // Read every 250 ms, so that no counter wraps between two reads.
        s.spawn(|| {
            while !done.load(Ordering::Relaxed) {
                std::thread::sleep(std::time::Duration::from_millis(250));
                let _ = acc.lock().unwrap().read();
            }
        });
        let out = f();
        done.store(true, Ordering::Relaxed);
        let ok = acc.lock().unwrap().read().is_some();
        (out, ok)
    });
    let secs = started.elapsed().as_secs_f64();
    let mut figure = if ok {
        package_figure(Some(&acc.lock().unwrap().reading(0)), "")
    } else {
        json!({ "unmetered": "a package counter did not read at the end of the run" })
    };
    if figure.get("joules_micro").is_some() {
        figure["covers"] = json!("the energy of the host's processor packages over the run: every process on the host in that window, this one among them, and the packages' idle draw");
        figure["source"] = json!("powercap energy_uj of every package zone, read before the run, every 250 ms during it and after it, each wrap corrected at the zone's max_energy_range_uj");
        figure["host_load_1m_at_start"] = json!(load);
        if let Some(m) = figure.as_object_mut() {
            m.remove("concurrent_jobs");
        }
    }
    (out, secs, figure)
}

fn metered_process<T>(f: impl FnOnce() -> T) -> (T, f64, Value) {
    use wai_quantum::engine::energy::own;
    let before = own();
    let started = Instant::now();
    let out = f();
    let secs = started.elapsed().as_secs_f64();
    let after = own();
    let figure = match (before, after) {
        (Some(a), Some(b)) => {
            let nj = |r: &wai_quantum::engine::energy::Reading| r.energy_nj.saturating_add(r.billed_nj).saturating_sub(r.serviced_nj);
            let delta = nj(&b).saturating_sub(nj(&a));
            let joules_micro = delta.saturating_add(500) / 1000;
            if joules_micro == 0 {
                json!({ "unmetered": "the process's energy counter did not advance over the run" })
            } else {
                json!({
                    "joules_micro": joules_micro,
                    "energy_provenance": "OnChipCounter",
                    "energy_counter_filtering": "Unknown",
                    "attribution": "process",
                    "coverage": "partial",
                    "covers": "the CPU energy of this process's own threads over the run, with what system services billed to it",
                    "excludes": ["memory", "the host's idle draw", "energy of threads still on a core at the second read"],
                    "cpu_ms": b.cpu_ns.saturating_sub(a.cpu_ns) / 1_000_000,
                    "source": "proc_pid_rusage ri_energy_nj + ri_billed_energy − ri_serviced_energy of this process, read before and after the run",
                })
            }
        }
        _ => json!({ "unmetered": "this platform gives no per-process energy counter here" }),
    };
    (out, secs, figure)
}

/// The signing key: 64 hex digits in the file `GBS_SIGNING_KEY` names, or,
/// with none, a key every reader knows, which the receipt then says it was
/// sealed with.
pub fn signing_key() -> (SigningKey, &'static str) {
    match std::env::var("GBS_SIGNING_KEY") {
        Ok(path) => {
            let text = std::fs::read_to_string(&path).unwrap_or_else(|e| panic!("reading the key {path}: {e}"));
            let text = text.trim();
            assert_eq!(text.len(), 64, "the key file holds 64 hex digits");
            let mut b = [0u8; 32];
            for (i, x) in b.iter_mut().enumerate() {
                *x = u8::from_str_radix(&text[2 * i..2 * i + 2], 16).expect("hex digits");
            }
            (SigningKey::from_bytes(&b), "the operator's key")
        }
        Err(_) => (SigningKey::from_bytes(&[42u8; 32]), "a demonstration key every reader knows (set GBS_SIGNING_KEY to sign with your own)"),
    }
}

/// Seal `result` of `request` and write the three, as one JSON document, to
/// `out` (when given). Returns the document.
pub fn seal(request: Value, mut result: Value, out: Option<&Path>) -> Value {
    let (key, which) = signing_key();
    result["signed_with"] = json!(which);
    let receipt = wai_quantum::engine::receipt::seal(&key, &request, &result);
    let doc = json!({ "request": request, "result": result, "receipt": receipt });
    if let Some(path) = out {
        if let Some(dir) = path.parent() {
            std::fs::create_dir_all(dir).unwrap_or_else(|e| panic!("creating {}: {e}", dir.display()));
        }
        std::fs::write(path, serde_json::to_string_pretty(&doc).unwrap() + "\n").unwrap_or_else(|e| panic!("writing {}: {e}", path.display()));
    }
    doc
}

/// Whether a written receipt document's seal holds over its request and result.
pub fn verifies(doc: &Value) -> bool {
    wai_quantum::engine::receipt::verify(&doc["request"], &doc["result"], &doc["receipt"])
}

/// The number as its exact bits, for a reader that must get the same double back.
pub fn bits(x: f64) -> String {
    format!("{:016x}", x.to_bits())
}