Can GPU-initiated discrete simulated bifurcation handle low-latency requests and streaming dense couplings?

Two bottlenecks, one communication path

The paper frames GPU-based optimization as facing two distinct scaling limits that share the same communication path. At small problem sizes, the model fits in GPU memory, but repeated host coordination can cap how many concurrent decisions the service completes. At large problem sizes, the dense coupling matrix needs O(n²) storage while the O(Bn) state of B solver replicas can stay resident, so scaling becomes a coefficient-delivery problem rather than a request-rate problem.

The proposed architecture addresses both through NVIDIA DOCA GPUNetIO, which exposes receive and transmit operations to CUDA kernels. For resident models, a persistent service receives field updates, executes each solve inside one GPU thread block, and returns the result. Exact integer coupling sums, GPU work queues, and batched transmission keep the receive-solve-reply path on the device without a dedicated CPU data-path core. For models larger than device memory, a streaming solver retains dynamical state on the GPU and reuses incoming coupling tiles across replicas.

The distinction matters because the two modes stress different resources. The resident mode is judged on scheduling overhead across many concurrent short solves; the streaming mode is judged on sustained bandwidth while solver state stays put. The paper evaluates them separately rather than reporting one blended number.

One request per block keeps solves independent

Each request is solved by one thread block of 256 threads. With n padded to a multiple of 32, the B·n_p (replica, row) pairs are distributed cyclically across threads, and each thread keeps the position and momentum of its at most 16 pairs in registers for the whole trajectory, so B·n_p ≤ 4096. Because n_p and the thread count are multiples of 32, the 32 lanes of a warp always hold 32 consecutive rows of the same replica, so the new signs of a warp form exactly one 32-bit word of that replica's sign vector via __ballot_sync, with no atomics or shuffles.

The server fixes one shared-memory footprint S for all workers, which sets how many blocks are resident per streaming multiprocessor and therefore how many requests are in service at once. With S = 48 KiB, two blocks fit on each of the 188 SMs of the server GPU, giving 376 concurrent requests. Four block slots are left free for driver work.

The paper reports that in a first version every spinning block read the shutdown flag on every iteration; 372 blocks issuing PCIe reads in a tight loop starved the NIC's own DMA and doorbells and raised round-trip time to milliseconds. Each block now polls the flag no more than once every 0.5 ms, reducing PCIe traffic. This is a concrete example of how host monitoring and GPU occupancy interact with latency.

Exact integer sums make GPU and CPU trajectories match

The solver's coupling sums are computed exactly rather than approximately. For ±1 couplings stored as bit planes, the sum over j of J_ij·s_j equals n − 2·popcount(u_i ⊕ v) + s_i. For ternary couplings in {−1, 0, 1}, the sum equals 2·(popcount(P_i ∧ v) − popcount(M_i ∧ v)) − (|P_i| − |M_i|). For int8 couplings with |J_ij| ≤ 127, |h_i| < 2^15 and n ≤ 2048, the int32 DP4A sum is exact and |a_i| < 2^24, so it is also exact in FP32.

The consequence the paper draws is that if the remaining arithmetic is evaluated with the same IEEE round-to-nearest operations in the same order, the GPU and a CPU reference produce identical trajectories, spins, and energies. The GPU update uses explicitly rounded multiplication, addition, and subtraction (__fmul_rn, __fadd_rn, __fsub_rn) to prevent fused multiply-add contraction, and the host reference is compiled with -ffp-contract=off.

This matters for readers comparing implementations: the claim is not that the GPU result is close to the reference, but that the coupling sum is the same integer on both sides and every remaining operation is a correctly rounded IEEE operation on identical operands. The exactness bound is stated for n ≤ 2048, which is a scope condition rather than a general guarantee.

Latency gains appear under concurrent load

The headline comparison is against socket-based servers using the same solver. At 800 thousand requests per second offered, GPUNetIO completes 618–628 kreq/s and the tuned proxy 620–668 kreq/s without loss, with the client limiting injection at this load. GPUNetIO's median round-trip time is then 189 µs against 288 µs for the tuned proxy, and its 99th-percentile round-trip time is 218 µs against 609 µs — reductions by factors of 1.5 and 2.8. At 600 kreq/s the corresponding values are 188 against 240 µs and 216 against 327 µs; at 400 kreq/s the two are still close, at 178 against 193 µs and 207 against 229 µs.

The scaling behavior is the more decision-relevant number. Between 400 and 800 kreq/s offered, GPUNetIO's median round-trip time rises by 6% while the proxy's rises by 50%. Server-side timestamps help explain this: over the same range, GPUNetIO's median server residence stays at 156–159 µs and its 99th percentile at 180–186 µs, whereas the proxy's rise from 173 to 235 µs and from 208 to 317 µs.

The paper is explicit that the unloaded latency is close to the tuned proxy's. Without load, the main gains are the freed server CPU core and lower host sensitivity; the larger latency difference appears under concurrent load. Each CPU-mediated baseline dedicates a polling core to receiving requests and returning results, while GPUNetIO performs that work on the device after setup.

Batched replies avoid a send-side bottleneck

The paper reports a failure mode that batched transmission is designed to avoid. At the highest offered load, per-reply sends complete only 334 kreq/s with a median round-trip time of 24 ms and 41% request loss. Batched descriptor posting avoids this bottleneck at the tested rates, keeping reply submission from limiting request service.

The comparison also has a boundary the paper states directly. The launch server reaches its four-stream limit much earlier, near 24 kreq/s, so the launch path is not the relevant baseline at the highest tested loads. At the two highest rates, where the client is saturated, individual runs of both servers occasionally show millisecond-scale 99th percentiles; at 700 kreq/s this places the median 99th percentile for GPUNetIO (423.5 µs) slightly above that of the proxy (381 µs).

Across the 19 settings in the closed-loop grid, the GPU-initiated service has the lowest recorded median in 15, improving on the faster tuned baseline by 1.5–9.2 µs; the remaining four settings differ by 1.5–2.5 µs in favor of a baseline. Its recorded server residence is lower throughout the grid. The paper notes that the residuals use different timestamp boundaries, which is why it reports residence and the rest of the round trip separately.

Streaming trades memory for bandwidth

For models larger than device memory, the streaming solver receives the coupling matrix tile by tile while positions, momenta, signs, and row accumulators remain on the receiving GPU. A dense one-bit representation has logical size n²/8 bytes. Frames carry a 64-byte header including Ethernet and 4032 payload bytes, representing 32,256 couplings from one row, and the final column block is padded.

For n = 10^7, the logical matrix is 12.5 TB, while a pass carries 12.53952 TB in 3.11×10^9 frames. The receiver has one cyclic queue of 16,384 slots of 4096 bytes and publishes bursts to 372 persistent worker blocks of 512 threads each. At n = 10^7 and B = 8, the row accumulators require 320 MB. The paper reports that the 12.5-TB logical matrix is consumed through a 64-MiB packet buffer.

The measured rate is approximately 307 Gb/s while performing the coupling accumulation. A 4096-byte frame occupies 4120 bytes including FCS, preamble, and interpacket gap, giving an ideal payload ceiling of 391.46 Gb/s on the 400-GbE link, so the observed rate reaches about 78% of that ceiling. The paper attributes the gap to the current PCIe 5.0 server's GPU-to-NIC paths traversing host interconnects, and names a direct PCIe-switch path, sufficient upstream bandwidth, and balanced queue placement as candidates for closing it.

Validation rests on planted instances and reference runs

Validation combines stream accounting, planted-solution checks, and reference execution. The receiver compares parsed stream-frame totals with end markers and records dropped bursts. Mattis instances provide the exact ground-state value 2E = −n(n−1) for checking the accumulated energy and solver integration. At n = 10^5, the Mattis and random energy and checksum records match a resident GPU reference byte for byte; at n = 10^4, the random records also match a CPU reference. These records contain the best-replica energy and a 64-bit spin-bitmap checksum at each recorded step.

The streaming runs use seed 42, Δt = 1.1, and ξ = 0.5/√n. The n = 10^7, B = 1 run uses two updates; its B = 8 runs use ten. Reference comparisons at n = 10^4 and 10^5 use 20 updates. A run of K updates requires K + 1 matrix passes, with the last pass evaluating the final energy without advancing the state.

The reported receive-pass interval ends when the host observes the end markers, which it polls every 2 ms. Full iteration time also includes draining outstanding work, updating the state, and acknowledging the pass, and these control operations use host polling. The paper states that the evaluated implementation uses mapped host memory for its receive and control counters, and that the cyclic packet buffer holds approximately 1.7 ms of arriving frames at the measured rate.

What the two modes do not yet cover

The request evaluation uses one platform and two socket-based baselines. The paper states that a matched DPDK/GPUdev baseline using the same solver, packet format, and batching policy is a useful next comparison to separate GPU-initiation gains from kernel-bypass gains. This is a scope condition on the latency comparison: the measured advantage is against socket-based servers, not against a kernel-bypass baseline with matched batching.

The current request bounds remain one frame per solve, one 512-KiB model slot, and the register and shared-memory limits described in the mapping. Deployment extensions the paper names include synchronized replacement of in-use models, retransmission and per-frame sequence validation, and distributing streaming work across multiple receivers.

On the streaming side, the paper notes that the measured 306–307 Gb/s corresponds to approximately 38.3 GB/s of coefficient payload, and that a Samsung 9100 PRO PCIe 5.0 NVMe SSD specifies sequential reads up to 14.8 GB/s, so feeding the measured stream from such storage would call for parallel drives. The coefficient source must sustain this bandwidth to feed the receiver continuously. The paper also observes that the similar rates for one and eight replicas suggest room to increase delivery bandwidth while retaining replica reuse.

Where the two modes diverge in practice

The paper's own split is between request-rate scaling and model-size scaling, and the two produce different decision thresholds. If the model fits in GPU memory and the workload is many concurrent short solves, the relevant comparison is the 400-to-800 kreq/s window, where GPUNetIO's median round-trip time rises 6% against the tuned proxy's 50%, and where the 99th percentile at 800 kreq/s is 218 µs against 609 µs. The unloaded latency is not where the advantage lives.

If the model exceeds device memory, the relevant number is the sustained coefficient-delivery rate: approximately 307 Gb/s against a 391.46-Gb/s payload ceiling on the 400-GbE link, with each of K updates requiring another matrix pass and a final energy-only pass. At the ideal payload ceiling of the same link, the transfer term for the n = 10^7 configuration would be approximately 256 s instead of the present roughly 327 s per receive pass; an 800-GbE link with the same frame efficiency would reduce the ideal term to approximately 128 s.

The most useful reading is therefore conditional rather than a single verdict. The GPU-initiated path pays off most where host coordination would otherwise grow with arrival rate and queue count, and the streaming path pays off where the quadratic coefficients cannot be retained but the solver state can. The paper's own framing — that the two modes place different demands on the same communication path — is the distinction worth carrying away.

Sources used in this answer