Source-linked AI summary
A single-precision floating-point systolic Givens-QRD Triangular Solver for MVDR Beamforming
Athi Ram R S, Alwin A, S. G. Sreejeesh, J. U. Kidav
TL;DR
MVDR weight computation is a latency-critical hardware challenge, and direct covariance inversion is resource-intensive and difficult to pipeline. The paper implements a float32 systolic Givens-QRD FPGA pipeline; at matched dispatch, the Xeon is faster in raw throughput, while the FPGA leads in power-normalised throughput and shows near-theoretical finite-precision accuracy.
Problem
Direct covariance-matrix inversion is a computational bottleneck for MVDR, with high resource use, poor scalability, and difficulty in efficient pipelining.
Source: §I · p. 1Method
The paper implements a systolic Givens-QRD pipeline in IEEE float32, mapping the arithmetic to FPGA DSP48E2 fabric.
Source: §I · p. 1, §VI · p. 17Results
At matched three-way dispatch, FPGA power-normalised throughput is 2.3× (programmable logic) to 1.05× (total on-chip) the Xeon’s, despite lower raw throughput; validation reports 4.27 εf global RMS error and a 100 % pass rate across 77,544 comparisons.
Source: §VI · p. 17Takeaways & Limitations
The float32 implementation avoids fixed-point condition-number sensitivity and range-design overhead in the reported FPGA design.
Source: §VI · p. 17
Abstract
from arXiv · showhide
Computation of adaptive beamforming weights in Minimum Variance Distortionless Response (MVDR) processing is a latency-critical operation that poses significant challenges for real-time hardware implementation. This paper presents an FPGA implementation of a systolic Givens-rotation QR decomposition pipeline for MVDR beamforming on a simulated 32-element ultrasound transducer array, using single-precision floating-point arithmetic. The design is deployed on a Zynq UltraScale+ FPGA at 100 MHz with three parallel kernel instances operating concurrently, achieving a measured throughput of 31,123 weight vectors/s at 90.9% parallel efficiency relative to the measured single-instance rate. At an estimated 2.451 W of programmable-logic power, and 5.286 W including the processing system, this corresponds to 12,698 and 5,888 weight vectors/s/W, respectively. Under a matched three-way dispatch, a 24-core Intel Xeon Gold 5220R at 2.20 GHz achieves 465,699 weight vectors/s at a measured 83.08 W package power, corresponding to 5,606 weight vectors/s/W. The FPGA therefore attains 2.3x the power-normalised throughput of the processor on a programmable-logic basis and 1.05x on a total on-chip basis. In contrast, the processor retains a raw throughput advantage of approximately 15x at this operating point. Numerical precision is validated against MATLAB float32 reference outputs from a Field II cyst phantom simulation, achieving a 100% pass rate with a root mean square error of 5.10x10^-7, confirming near-theoretical finite-precision behaviour without systematic bias.
I. INTRODUCTION
The paper targets the computational cost and numerical challenges of MVDR weight computation by replacing direct covariance inversion with Givens-rotation QR decomposition. It presents a float32 FPGA implementation and validates its weights against a software reference using Field II phantom data; imaging quality is assessed only qualitatively.
I. INTRODUCTION
- Direct covariance-matrix inversion is the fundamental MVDR bottleneck, with numerical instability, computational overhead, high resource use, poor scalability, and difficult pipelining.Constructing and inverting the matrix requires processing correlations across the P × L complex matrix.Source: §I · p. 1
- The float32 systolic Givens QRD implementation avoids fixed-point range design and CORDIC approximation, retaining explicit square-root and reciprocal operators.The paper identifies accumulated rounding in iterative fixed-point CORDIC rotations as a source of increasing error and poorer factorisation conditioning.Source: §I · p. 1
- The closest prior work,, reports 68.3 fps for a full 61-line imaging chain; this work isolates its computationally dominant adaptive-weight stage.The present implementation uses IEEE float32 arithmetic.Source: §I · pp. 1–2
A. Contributions of This Paper
- 31,123 weight vectors/s at 90.9 % parallel efficiency; three FPGA instances achieve 2.3× (programmable logic) to 1.05× (total on-chip) power-normalised throughput versus the Xeon.At matched three-way dispatch, the Xeon has higher raw throughput: 465,699 weight vectors/s.Source: A. Contributions of This Paper · p. 2
A. Role of Field II in This Work
- Field II generates phantom covariance matrices for hardware and software weight-vector comparison; the study reports no quantitative imaging-quality metrics.The B-mode panels are included only for qualitative reference.Source: A. Role of Field II in This Work · p. 2
B. Spatial Smoothing and Matrix Construction
The paper constructs float32 spatially smoothed inputs from overlapping subarrays, recovers MVDR weights without forming R^-1, and evaluates a Givens QRD reference against a dataset of 9,693 snapshots.
B. Spatial Smoothing and Matrix Construction
- M = 32, L = 8, and P = 25 define the overlapping-subarray construction of the float32 spatially smoothed input matrix U[n].Each received snapshot is partitioned into P overlapping subarrays of length L.Source: B. Spatial Smoothing and Matrix… · pp. 2–3
C. MVDR Weight Computation
- MVDR weights are recovered without forming R^-1, using the transform T[n] = S^-H[n] on the data and to form a_z[n].The weights are applied to transform-domain data rather than raw element data.Source: C. MVDR Weight Computation · p. 3
D. Software Reference and Dataset
- The software reference uses the identical Givens rotation sequence in explicit float32 arithmetic, not MATLAB's Householder-based QR function.The different diagonal sign convention makes MATLAB's built-in QR unsuitable for bit-comparison with the Givens path at float32 precision.Source: D. Software Reference and Dataset · p. 3
- 9,693 snapshots from RF scan line 1 yield 77,544 complex float32 weight-element comparisons against the software reference.The FPGA processes 400 float32 input words and outputs 16 float32 words per snapshot.Source: D. Software Reference and Dataset · p. 3
III. HARDWARE ARCHITECTURE
The diagonal processing element computes complex Givens rotation parameters and handles near-zero inputs, while off-diagonal elements apply the rotations to update matrix rows.
1) Diagonal PE:
- The diagonal PE derives rotation parameters (c, s) from the complex top and bottom elements, using magnitude accumulation and reciprocal-based normalization.The normalization factor is computed by squaring and accumulating the real and imaginary components before the reciprocal stage.Source: 1) Diagonal PE · p. 3, Figure 1 · p. 4
- For p < 10^-10, guarded logic bypasses the reciprocal path and defaults the rotation to an identity-like configuration.The diagonal PE also annihilates the bottom complex value during the local matrix update.Source: 1) Diagonal PE · p. 3
- Off-diagonal PEs receive (c, s) horizontally and stream rows u and v vertically to produce the updated row pair.Four complex MAC units implement the row-pair update.Source: Figure 1 · p. 4, 2) Off-Diagonal PE · p. 4
2) Off-Diagonal PE:
The off-diagonal PE expands the complex row update into real and imaginary components and implements the products with multiplier groups and controlled add/subtract logic.
2) Off-Diagonal PE:
- Three multiplier groups handle the c · u, s · v, and s · u products, with sign inversions absorbed into controlled add/subtract logic.This implementation preserves pipeline regularity and costs 28 FLOPs per call.Source: 2) Off-Diagonal PE · p. 4
3) Systolic Array Grid:
The systolic Givens QRD maps diagonal-to-off-diagonal dependencies onto a triangular PE array, reducing sequential latency relative to the operation count. Across synthesized configurations, cycle counts track sequential depth, while the QRD stage dominates the measured kernel latency.
3) Systolic Array Grid:
- The per-snapshot workload is 20,753 FLOPs: 19,380 for QRD and 1,373 for forward substitution, steering accumulation, and normalization.The solver term depends on L, while the QRD term depends on the full subarray geometry.Source: 3) Systolic Array Grid · p. 4, 4) Asymptotic Complexity · p. 6
- For the actual P = 25, L = 8 design, the triangular array spans eight columns, from 24 row pairs in column 0 to 17 in column 7.The number of active PE rows decreases across columns.Source: 3) Systolic Array Grid · p. 4, Figure 3 · p. 5
- Diagonal PEs chain vertically and broadcast rotation parameters (c, s) to off-diagonal PEs, while R-matrix boundary elements accumulate as columns complete.The scheduling staircase follows the dependency: an off-diagonal PE waits for its same-row diagonal PE to produce (c, s).Source: 3) Systolic Array Grid · p. 4, Figure 2 · p. 5
4) Asymptotic Complexity:
- A 16× increase in FLOPs from Case 5 to Case 3 yields only a 5.3× increase in cycles, which track sequential depth D rather than N_QRD.The solver fraction rises from about 7% when P ≫ L to about 28% as L approaches P in Case 3.Source: Table V · p. 6, 4) Asymptotic Complexity · p. 6
- 9,181 of 9,430 cycles (97.4%) are spent in load-and-QRD, leaving 249 cycles for the downstream solver chain.The reported stage counts describe where synthesized latency is spent.Source: 4) Asymptotic Complexity · p. 6
B. Arithmetic Precision: Float32 Without CORDIC
The design uses float32 Givens rotations rather than fixed-point CORDIC to address range and quantization concerns in ultrasound covariance processing. Its five-stage pipeline is dominated by load-and-QRD, and three parallel FPGA instances achieve 31,123 weight vectors/s at 90.9% of the theoretical three-instance limit.
B. Arithmetic Precision: Float32 Without CORDIC
- At high-condition snapshots, fixed-point forward substitution risks pivot overflow beyond the WL=16, FL=14 format’s approximately 2.0 maximum, whereas float32 has a 3.4 × 10^38 range.The effective fixed-point precision at the weight output degrades with covariance condition number.Source: B. Arithmetic Precision: Float32… · p. 6
- CORDIC is avoided because fixed-point rotations require float-to-integer conversion, 15–20 iterations offer no timing advantage over native square root, and its scale factor requires compensation.The direct float32 formulation produces |c|^2 + |s|^2 = 1 without correction.Source: B. Arithmetic Precision: Float32… · p. 6
C. HLS Dataflow Pipeline
- Forward substitution is fully parallelized by loop unrolling and complete array partitioning, while column-wise Givens operations remain sequential because of their dependency chain.The systolic inner loop over row pairs is pipelined.Source: C. HLS Dataflow Pipeline · p. 6
- 91,810 ns (9,181 cycles), or 97.4% of the 94,300 ns kernel latency, is spent in load-and-QRD; downstream stages overlap this interval.The five dataflow stages are load-and-QRD, forward substitution, steering accumulation, weight normalization, and output store.Source: C. HLS Dataflow Pipeline · p. 6, Figure 4 · p. 7
D. System Integration
- 31,123 weight vectors/s reaches 90.9% of the three-instance theoretical limit of 34,230 weight vectors/s.The reported gap is attributed to the sequential DMA trigger window of 2.23 µs per IP between dispatch cycles.Source: Table IX · p. 7, D. System Integration · pp. 7–8
- Three kernel instances use dedicated DMA channels and independent high-performance DDR4 ports, allowing concurrent transfers without shared-bus contention.The MMIO dispatch path removes per-engine driver overhead relative to the driver-managed DMA descriptor interface.Source: D. System Integration · pp. 6–7, Figure 5 · p. 8
IV. NUMERICAL ACCURACY ANALYSIS
The accuracy evaluation uses Field II cyst-phantom data and compares float32 with a 16-bit fixed-point path derived from the same double-precision source. The reported fixed-point weight error peaks at a high-condition specular snapshot, while float32 and fixed-point traces both track the signal envelope.
A. Float32 vs Fixed-Point Dynamic Range
- Both evaluation inputs derive from the same double-precision Field II cyst-phantom RF data, isolating numerical differences to the subsequent precision conversion.The simulation uses a 32-element linear array.Source: A. Float32 vs Fixed-Point Dynamic Range · p. 8
- Fixed-point RF lines are pre-scaled per scan line before quantization and normalized to unit peak; a 10^-6 pivot guard is used only in the fixed-point reference.The FPGA float32 pipeline factorizes U[n] directly and requires no such guard.Source: A. Float32 vs Fixed-Point Dynamic Range · p. 8
- At a high-condition specular snapshot, the fixed-point weight error peaks at approximately 2.43, while its magnitude trace is attenuated relative to float32.Both traces track the signal envelope and peak near sample 1.9 × 10^4.Source: Figure 6 · p. 9
- At specular tissue boundaries, κ(U[n]) reaches O(10^4)–O(10^6), beyond the fixed-point range’s ability to represent both the pivot and off-diagonal update.The passage reports float32’s range as eliminating overflow for physically encountered ultrasound covariance.Source: A. Float32 vs Fixed-Point Dynamic Range · p. 9
1) Empirical Weight Divergence Under Fixed-Point Arithmetic:
On identical Field II covariance inputs, fixed-point 16 produces severe MVDR weight errors, including loss of imaginary content and phase instability. Float32 errors instead show the signature of rounding, with indistinguishable traces across the three IP instances.
1) Empirical Weight Divergence Under Fixed-Point Arithmetic:
- The fixed-point 16 weight error norm is 1.43 (143%), with absolute error reaching ≈2.43 on the unit-normalised weight scale.The fixed-point output is systematically attenuated across the imaging zone, while float32 tracks the expected envelope.Source: 1) Empirical Weight Divergence Under… · p. 9, Figure 7 · p. 10
- Fixed-point relative error is centred near 100 across the imaging zone, with a maximum of ∼50× and a heavy-tailed distribution.Most absolute errors are below 0.1, with the bulk of the tail ending near 0.35 apart from a specular outlier at ≈2.43.Source: 1) Empirical Weight Divergence Under… · p. 9
- The fixed-point imaginary weights collapse to nearly zero despite float32 signal content varying across ±0.15, and phase error reaches ≈5.85 rad.The real fixed-point component also exhibits large, uncorrelated excursions across the dataset.Source: 1) Empirical Weight Divergence Under… · p. 9
- Float32 supports L = 16 without overflow, while the fixed-point 16 reconstruction is limited to L = 8 in the qualitative phantom comparison.No quantitative imaging metric is computed from these panels.Source: 1) Empirical Weight Divergence Under… · p. 9
1) Per-element error and error trace:
- Float32 error grows monotonically from w0 to w7 along the forward-substitution dependency chain, while all three IP traces are visually indistinguishable.The hardware-to-software difference is 10−13–10−14; the cumulative error distribution reaches a 100% pass rate, with median imaging-zone error of 4.51 ULP.Source: 1) Per-element error and error trace · pp. 9–10, 2) Cumulative distribution and ULP error · p. 10
3) Error attribution and validation:
The error analyses attribute float32 deviations to rounding proportional to output magnitude, with no systematic offset, and report a 100% pass rate. The paper also contrasts this behaviour with fixed-point failure and reports power-normalised comparisons under matched dispatch.
3) Error attribution and validation:
- The error-versus-output-magnitude fit has a log-log slope of 0.96 with R2 = 0.91, consistent with error scaling with magnitude and no systematic offset or bias.The table summary attributes global RMS error to approximately 24 complex Givens rotations plus the forward-substitution dependency chain.Source: 3) Error attribution and validation · p. 10, 4) Error breakdown by signal zone · p. 11
4) Error breakdown by signal zone:
- The zone breakdown resolves into three non-overlapping clusters mapping to the dataset’s three physical zones; a single diagonal band in magnitude scatter confirms error depends entirely on input amplitude, not hardware condition or zone-specific behavior.Source: 4) Error breakdown by signal zone · p. 10
A. Novelty Relative to Prior QRD-MVDR FPGA Architectures
- The implementation had no overflow for any condition number encountered in the Field II dataset (κ up to O(106)); fixed-point 16 (WL=16, FL=14) produced 143% overall weight error and 5.85 rad maximum phase error.Source: A. Novelty Relative to Prior QRD-MVDR… · p. 12
2) Cumulative distribution and ULP error:
- 4.27 εf global RMS error and a 100% pass rate across 77,544 comparisons support near-theoretical finite-precision behaviour, with no NaN or Inf outputs.The reported log |e|-vs-log |w| slope is 0.96 (R2 = 0.91).Source: 4) Error breakdown by signal zone · p. 11, A. Novelty Relative to Prior QRD-MVDR… · p. 13
A. Novelty Relative to Prior QRD-MVDR FPGA Architectures
- At matched three-way dispatch, FPGA power-normalised throughput is 2.3× the processor on a programmable-logic basis and 1.05× on a total on-chip basis.The reported comparison uses estimated FPGA power and measured processor package power; at its 48-thread throughput ceiling, the processor has higher power-normalised throughput than either FPGA basis.Source: A. Novelty Relative to Prior QRD-MVDR… · p. 13
B. Comparison with Prior Work
The paper cautions that prior FPGA throughput figures are not directly comparable across differing matrix sizes, arithmetic, pipeline scope, and application domains. Its closest comparison is prior QRD-MVDR work, which shares the algorithm, domain, and matrix dimensions but reports qualitative rather than weight-level accuracy.
B. Comparison with Prior Work
- Prior FPGA QRD implementations target smaller 4×4 real matrices and omit forward substitution and weight normalisation, so their throughput is not comparable to MVDR weight vectors per second.The table comparison spans differing dimensions, arithmetic precision, pipeline scope, and application domains.Source: B. Comparison with Prior Work · pp. 13, 15
- Zhang et al. report 10,557 fps and Kou et al. report 4.83 Gsamples/s, but both target full-chain conventional beamforming rather than isolated float32 adaptive weight computation.Their fixed-arithmetic systems therefore differ in scope from the present 25 × 8 complex covariance-input stage.Source: B. Comparison with Prior Work · p. 15
- The closest work uses the same algorithm, domain, and P = 25, L = 8 matrix dimensions, but reports qualitative accuracy rather than weight-level error characterisation.Its 0.25 ms per beam covers a complete scan line with 512 frequency-domain MVDR computations and an FFT/IFFT front-end; this work measures an isolated time-domain weight-computation stage.Source: B. Comparison with Prior Work · p. 15
C. Throughput and Power-Normalised Efficiency
The benchmark compares matched kernel boundaries and group-synchronous three-way dispatch on FPGA and CPU. At 100 MHz, the FPGA reaches 31,123 weight vectors/s, while matched-dispatch CPU throughput is higher; power-normalised ordering depends on the FPGA power basis.
Benchmarking Methodology:
- Both platforms time one complete MVDR weight computation, using a 400-word float32 input and 16-word output, with no I/O outside the kernel boundary.The CPU uses three synchronised worker threads to mirror the FPGA's group-synchronous dispatch of three instances.Source: Benchmarking Methodology · p. 15
Raw Throughput:
- 31,123 weight vectors/s is achieved by the FPGA three-IP MMIO configuration at 100 MHz, reaching 90.9% of its measured three-instance ceiling.Matched CPU-4 three-thread execution reaches 465,699 weight vectors/s; CPU-5's three independent MKL workers reach 860,504 weight vectors/s but do not use group-synchronous dispatch.Source: Raw Throughput · pp. 15–16, Benchmarking Methodology · p. 15
Power-Normalised Efficiency:
- At matched dispatch, FPGA power-normalised throughput is 2.3× the processor using 2.451 W programmable-logic power and 1.05× using 5.286 W total on-chip power.The FPGA power is estimated from switching-activity analysis, whereas the processor's 83.08 W package power is measured; at 48 threads the processor reaches 74,153 weight vectors/s/W.Source: Power-Normalised Efficiency · p. 16, Table XIV · p. 16, Table XIII · p. 16
- The matched-dispatch efficiency difference is attributed to the power floor: processor idle package power is 50.07 W, versus 5.286 W total FPGA on-chip power.The per-operation times are comparable in order of magnitude: 2.15 µs for the processor and 2.23 µs for FPGA dispatch.Source: Power-Normalised Efficiency · p. 16
D. Resource Utilisation and Scalability
Three FPGA instances deliver 31,123 weight vectors/s, while the design closes timing above its 100 MHz deployment clock. The conclusion also reports numerical agreement with a software reference and identifies broader validation and integration as future work.
Timing Closure:
- WNS = +1.305 ns with zero failing endpoints yields approximately 115 MHz Fmax, above the 100 MHz deployment clock.The timing summary is reported in Table XV.Source: Timing Closure · p. 16, Table XV · p. 17
Scalability:
- 84.0% DSP and 76.9% LUT utilisation leave no headroom for a fourth IP instance on the ZCU104.Each additional IP adds 484 DSPs; 16 BRAM tiles per IP leave ample block-memory headroom.Source: Table XVI · p. 17, Scalability · p. 17
VI. CONCLUSION
- 31,123 weight vectors/s at 90.9% parallel efficiency; the FPGA has 2.3× the processor’s power-normalised throughput on a programmable-logic basis and 1.05× on a total on-chip basis.Under matched three-way dispatch, the Xeon reaches 465,699 weight vectors/s at 83.08 W package power; the FPGA has lower raw throughput.Source: §VI · p. 17
- The hardware achieves 100% pass rate across 77,544 float32 comparisons, with global RMS error 4.27 εf and no systematic bias.The log-log error-to-magnitude slope is 0.96 (R2 = 0.91).Source: §VI · p. 17
- Future work targets full 50-line phantom validation, FFT-front-end integration, alternative arithmetic for other domains, and wavefront scheduling to reduce square-root-induced stalls.Source: §VI · p. 17