wai-quantum 0.6.0

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), belief-propagation tensor networks on the hardware graph, 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 in the default features — identical results native, in the browser, and as a WASI component at the edge; the opt-in WebGPU member reaches the platform's GPU interface.
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95

// wai-quantum fabric: the gate kernel for OpenCL devices.
//
// A transcription of `lanes::exact_row` and `lanes::apply_words` in 64-bit
// integers, which OpenCL C's `long` is on every full-profile device. Amplitudes
// are pairs of 32-bit fixed-point values (30 fractional bits). Each product
// `fx(a, b) = (a·b + 2^29) >> 30` is exact in 64 bits (a 32×32 product is at
// most 2^62 in size), and so is each row's sum of four products. The lane
// stores the low 32 bits of each sum; a sum outside 32 bits records the
// gate's op-log index in `report[0]`, the lowest such index winning, so the
// host refuses the run with the gate the CPU reference names.
//
// Nothing here is floating point, depends on the order work-items run in, or
// reads an amplitude another work-item of the same launch writes: each
// work-item owns one pair (i, j), and a launch is one gate.
//
// OpenCL C fills the vacated bits of `>>` on a negative signed value with ones,
// as Rust's `>>` on i64 does; the host's probe checks this, among the rest of
// the arithmetic, before a device is a member.

typedef long i64;
typedef ulong u64;

// The exact fixed-point product.
inline i64 fx(int a, int b) {
    return ((i64)a * (i64)b + ((i64)1 << 29)) >> 30;
}

// The low 32 bits of `v`, as the lane keeps them.
inline int low32(i64 v) {
    return (int)(uint)(u64)v;
}

// Whether `v` fits a 32-bit lane.
inline int fits32(i64 v) {
    return v >= -(i64)2147483648L && v <= (i64)2147483647L;
}

// One gate on a state of `2 * pairs` lanes, one work-item per pair.
//
// `m0..m7` are the gate's words: row 0 is (m0 + i·m1)·a0 + (m2 + i·m3)·a1,
// row 1 is (m4 + i·m5)·a0 + (m6 + i·m7)·a1. `tbit` is 1 << target; `ctrl`
// the OR of the controls' bits. `index` is the gate's position in the op
// log. The host checks the words before a launch (`apply_words` refuses the
// same malformed words), so the kernel trusts them.
__kernel void apply_gate(
    __global int2* restrict amps,
    ulong pairs,
    uint tbit,
    uint ctrl,
    int m0, int m1, int m2, int m3, int m4, int m5, int m6, int m7,
    uint index,
    __global uint* restrict report)
{
    u64 p = (u64)get_global_id(0);
    u64 stride = (u64)get_global_size(0);
    for (; p < pairs; p += stride) {
        // The p-th index whose target bit is clear: p with a zero inserted at
        // the target's position.
        u64 low = p & ((u64)tbit - 1);
        u64 i = ((p - low) << 1) | low;
        if ((i & (u64)ctrl) != (u64)ctrl) {
            continue;
        }
        u64 j = i | (u64)tbit;
        int2 a0 = amps[i];
        int2 a1 = amps[j];
        i64 o0r = fx(m0, a0.x) - fx(m1, a0.y) + fx(m2, a1.x) - fx(m3, a1.y);
        i64 o0i = fx(m0, a0.y) + fx(m1, a0.x) + fx(m2, a1.y) + fx(m3, a1.x);
        i64 o1r = fx(m4, a0.x) - fx(m5, a0.y) + fx(m6, a1.x) - fx(m7, a1.y);
        i64 o1i = fx(m4, a0.y) + fx(m5, a0.x) + fx(m6, a1.y) + fx(m7, a1.x);
        if (!(fits32(o0r) && fits32(o0i) && fits32(o1r) && fits32(o1i))) {
            atomic_min(&report[0], index);
        }
        amps[i] = (int2)(low32(o0r), low32(o0i));
        amps[j] = (int2)(low32(o1r), low32(o1i));
    }
}

// The arithmetic probe: for each case (a, b), the 64-bit two's-complement
// product a·b and the exact fx(a, b), each as (high, low) words, the outputs
// of `lanes::mul_wide` and `lanes::fxmul_wide`.
__kernel void probe(
    __global const int2* restrict cases,
    __global uint4* restrict out,
    ulong count)
{
    u64 k = (u64)get_global_id(0);
    u64 stride = (u64)get_global_size(0);
    for (; k < count; k += stride) {
        int2 c = cases[k];
        u64 prod = (u64)((i64)c.x * (i64)c.y);
        u64 f = (u64)fx(c.x, c.y);
        out[k] = (uint4)((uint)(prod >> 32), (uint)prod, (uint)(f >> 32), (uint)f);
    }
}