1-bit Massive MU-MIMO Precoding in VLSI

TL;DR

This paper introduces novel nonlinear 1-bit precoding algorithms and VLSI architectures for massive MU-MIMO systems with low-resolution DACs, significantly improving error rates and enabling practical high-rate downlink transmission.

Contribution

It presents the first efficient nonlinear 1-bit precoding algorithms with VLSI designs for massive MU-MIMO, outperforming linear methods in error performance.

Findings

Superior error-rate performance of proposed algorithms

VLSI architectures enable practical implementation

FPGA results demonstrate high-rate reliable transmission

Abstract

Massive multiuser (MU) multiple-input multiple-output (MIMO) will be a core technology in fifth-generation (5G) wireless systems as it offers significant improvements in spectral efficiency compared to existing multi-antenna technologies. The presence of hundreds of antenna elements at the base station (BS), however, results in excessively high hardware costs and power consumption, and requires high interconnect throughput between the baseband-processing unit and the radio unit. Massive MU-MIMO that uses low-resolution analog-to-digital and digital-to-analog converters (DACs) has the potential to address all these issues. In this paper, we focus on downlink precoding for massive MU-MIMO systems with 1-bit DACs at the BS. The objective is to design precoders that simultaneously mitigate multi-user interference (MUI) and quantization artifacts. We propose two nonlinear 1-bit precoding algorithms and corresponding very-large scale integration (VLSI) designs. Our algorithms rely on biconvex relaxation, which enables the design of efficient 1-bit precoding algorithms that achieve superior error-rate performance compared to that of linear precoding algorithms followed by quantization. To showcase the efficacy of our algorithms, we design VLSI architectures that enable efficient 1-bit precoding for massive MU-MIMO systems in which hundreds of antennas serve tens of user equipments. We present corresponding field-programmable gate array (FPGA) implementations to demonstrate that 1-bit precoding enables reliable and high-rate downlink data transmission in practical systems.

Tables

Table 1. TABLE I: Implementation results for C1PO and C2PO for MU-MIMO systems with U = 16 𝑈 16 U=16 UEs on a Xilinx Virtex-7 XC7VX690T FPGA

Algorithm	C1PO					C2PO
BS antennas $B$	$32$	$64$	$128$	$256$		$32$	$64$	$128$	$256$
Slices	2 700	5 187	10 324	21 951		3 375	6 519	12 690	24 748
LUTs	6 671	13 305	30 979	71 817		10 817	21 920	43 710	85 323
– LUTs as logic	6 031	12 025	25 939	51 897		10 069	20 424	40 718	79 339
– LUTs as memory	640	1 280	5 040	19 920		748	1 496	2 992	5 984
Flipflops	6 830	13 624	26 683	52 175		5 677	12 461	26 083	53 409
DSP48 units	128	256	512	1 024		136	272	544	1 088
Max. clock frequency [MHz]	285	264	244	205		222	206	208	193
Min. latency666The minimum latency is measured for one algorithm iteration. [clock cycles]	35	67	131	259		39	40	41	42
Max. throughput777The throughput corresponds to the total number of symbols precoded per unit of time. In this case, the maximum throughput is equal to $\left(Uf\right)/d$, where $f$ is the maximum clock frequency and $d$ the minimum latency.   [Msymbols/s]	130	63	30	13		91	82	81	74
Power consumption888Statistical power estimation at maximum clock frequency and 1.0 V supply voltage. [W]	1.13	1.97	3.43	5.74		1.04	1.70	3.17	5.80
Max. throughput/LUTs	19 529	4 733	962	177		8 413	3 756	1 853	862

Table 2. TABLE II: Implementation results for a MRT-Q-based precoder for MU-MIMO systems with U = 16 𝑈 16 U=16 UEs on a Xilinx Virtex-7 XC7VX690T FPGA

BS antennas $B$	$32$	$64$	$128$	$256$
Slices	2 543	5 097	9 444	17 630
LUTs	7 842	15 617	32 476	64 446
– LUTs as logic	7 010	13 953	29 148	57 790
– LUTs as memory	832	1 664	3 328	6 656
Flipflops	5 711	11 419	21 902	42 764
Clock freq. [MHz]	412	410	388	359
Latency [cycles]	18	18	18	18
TP999The throughput is calculated as $\left(Uf\right)/d$, where $f$ is the clock frequency and $d$ the latency. [Msymbols/s]	366	365	345	319
Power101010Statistical power estimation at max. clock freq. and 1.0 V supply voltage. [W]	0.79	1.25	1.84	3.16
Throughput/LUTs	46 665	23 356	10 621	4 945

Equations

$$E_{\mathbf{n}}\left[\Vert \mathbf{s}-\hat{}\mathbf{s}{\Vert }_2^2\right]={\Vert \mathbf{s}-\beta \mathbf{Hx}\Vert }_2^2+|\beta {|}^{2}U{N}_0,$$

$$\displaystyle\mathbf{x}=\sqrt{\frac{P}{2B}}\big{(}\operatorname{sgn}\mathopen{}\left(\Re\mathopen{}\left\{\mathbf{P}\mathbf{s}\right\}\right)+j\operatorname{sgn}\mathopen{}\left(\Im\mathopen{}\left\{\mathbf{P}\mathbf{s}\right\}\right)\!\big{)}$$

$$\text{(OPP)}\underset{\mathbf{x}\in {\mathcal{X}}^B,\beta \in \mathbb{C}}{\mathrm{minimize}}{\Vert \mathbf{s}-\beta \mathbf{Hx}\Vert }_2^2+|\beta {|}^{2}U{N}_0.$$

$$\underset{\mathbf{x}\in {\mathcal{X}}^B}{\min }\underset{\beta \in \mathbb{C}}{\min }{\Vert \mathbf{s}-\beta \mathbf{Hx}\Vert }_2^2\approx \underset{\mathbf{x}\in {\mathcal{X}}^B}{\min }\underset{\alpha \in \mathbb{C}}{\min }{\Vert \alpha \mathbf{s}-\mathbf{Hx}\Vert }_2^2.$$

$$\text{(OPP∗)}\hat{\mathbf{x}}=\underset{\mathbf{x}\in {\mathcal{X}}^B,\alpha \in \mathbb{C}}{\mathrm{argmin}}{\Vert \alpha \mathbf{s}-\mathbf{H}\mathbf{x}\Vert }_2^2.$$

$$\hat{\alpha }(\mathbf{x})=\underset{\alpha \in \mathbb{C}}{\mathrm{argmin}}{\Vert \alpha \mathbf{s}-\mathbf{H}\mathbf{x}\Vert }_2^2=\frac{{\mathbf{s}}^H\mathbf{H}\mathbf{x}}{\Vert \mathbf{s}{\Vert }_2^2}.$$

$${\Vert \hat{\alpha }(\mathbf{x})\mathbf{s}-\mathbf{H}\mathbf{x}\Vert }_2^2={\Vert \mathbf{Ax}\Vert }_2^2$$

$$\mathbf{A}=\mathbf{QH}\text{and}\mathbf{Q}={\mathbf{I}}_U-\frac{\mathbf{s}{\mathbf{s}}^H}{\Vert \mathbf{s}{\Vert }_2^2},$$

$$\text{(OPP∗∗)}\hat{\mathbf{x}}=\underset{\mathbf{x}\in {\mathcal{X}}^B}{\mathrm{argmin}}{\Vert \mathbf{Ax}\Vert }_2^2,$$

$$\hat{\mathbf{x}}=\underset{\mathbf{x}\in {\mathcal{X}}^B,\mathbf{z}\in {\mathbb{C}}^B}{\mathrm{argmin}}{\Vert \mathbf{Az}\Vert }_2^2+\gamma \Vert \mathbf{z}-\mathbf{x}{\Vert }_2^2,$$

$$\displaystyle\mathcal{B}^{B}=\bigg{\{}\mathbf{c}\in\mathbb{C}^{B}\,\bigg{|}\,$$

$$\displaystyle|\Im\{c_{b}\}|\leq\sqrt{\frac{P}{2B}},\,b=1,\ldots,B\bigg{\}}{.}$$

$$\hat{\mathbf{x}}=\underset{\mathbf{x}\in {\mathcal{B}}^B,\mathbf{z}\in {\mathbb{C}}^B}{\mathrm{argmin}}{\Vert \mathbf{Az}\Vert }_2^2+\gamma \Vert \mathbf{z}-\mathbf{x}{\Vert }_2^2,$$

$$\text{(BCR∗)}{\hat{\mathbf{x}}}^{\text{BCR}}=\underset{\mathbf{x}\in {\mathcal{B}}^B,\mathbf{z}\in {\mathbb{C}}^B}{\mathrm{argmin}}{\Vert \mathbf{Az}\Vert }_2^2+\gamma \Vert \mathbf{z}-\mathbf{x}{\Vert }_2^2-\delta \Vert \mathbf{x}{\Vert }_2^2,$$

$${\mathbf{z}}^{(t+1)}$$

$${\mathbf{x}}^{(t+1)}$$

$${\mathbf{z}}^{(t+1)}$$

$${\mathbf{x}}^{(t+1)}$$

$$\mathrm{proj}(z)=$$

$$+j\mathrm{sgn}(\Im \{z\})\min \left\{\frac{\gamma }{\gamma -\delta }|\Im \{z\}|,\sqrt{\frac{P}{2B}}\right\}$$

$$\displaystyle\hat{}\mathbf{x}=\sqrt{\frac{P}{2B}}\bigg{(}$$

$$\displaystyle+j\operatorname{sgn}\mathopen{}\left(\Im\mathopen{}\left\{\mathbf{x}^{(t_{\text{max}}+1)}\right\}\right)\!\bigg{)}.$$

$$\hat{\mathbf{x}}=\underset{\mathbf{x}\in {\mathcal{B}}^B}{\mathrm{argmin}}\frac{1}{2}{\Vert \mathbf{Ax}\Vert }_2^2-\frac{\delta }{2}\Vert \mathbf{x}{\Vert }_2^2.$$

$$\hat{}\mathbf{x}=\underset{\mathbf{x}}{\mathrm{argmin}}f(\mathbf{x})+g(\mathbf{x}),$$

$${\mathbf{z}}^{(t+1)}$$

$${\mathbf{x}}^{(t+1)}$$

$${\text{prox}}_g\left(\mathbf{z};\tau \right)=\underset{\mathbf{x}}{\mathrm{argmin}}\left\{\tau g(\mathbf{x})+\frac{1}{2}\Vert \mathbf{x}-\mathbf{z}{\Vert }_2^2\right\}.$$

$$f(\mathbf{x})=\frac{1}{2}{\Vert \mathbf{Ax}\Vert }_2^2\text{ and }g(\mathbf{x})=\chi \left(\mathbf{x}\in {\mathcal{B}}^B\right)-\frac{\delta }{2}\Vert \mathbf{x}{\Vert }_2^2,$$

$${\mathrm{prox}}_g(z)=\mathrm{sgn}(\Re \{z\})\min \left\{\frac{1}{1-\tau \delta }|\Re \{z\}|,\sqrt{\frac{P}{2B}}\right\}$$

$$+j\mathrm{sgn}(\Im \{z\})\min \left\{\frac{1}{1-\tau \delta }|\Im \{z\}|,\sqrt{\frac{P}{2B}}\right\},$$

$${\mathbf{z}}^{(t+1)}$$

$${\mathbf{x}}^{(t+1)}$$

$${\mathbf{A}}^H\mathbf{A}={\mathbf{H}}^H\mathbf{H}-\frac{{\mathbf{H}}^H\mathbf{s}{\mathbf{s}}^H\mathbf{H}}{\Vert \mathbf{s}{\Vert }_2^2}={\mathbf{H}}^H\mathbf{H}-\mathbf{v}{\mathbf{v}}^H={\overline{\mathbf{H}}}^{{\rm Y}}\overline{\mathbf{H}}.$$

$${\mathbf{z}}^{(t+1)}={\mathbf{x}}^{(t)}-\tau {\overline{\mathbf{H}}}^{{\rm Y}}\overline{\mathbf{H}}{\mathbf{x}}^{(t)},$$

Full text

I Introduction

Massive multiuser (MU) multiple-input multiple-output (MIMO) is widely believed to be a core technology in fifth-generation (5G) wireless systems as it enables substantial improvements in spectral efficiency and reliability compared to traditional, small-scale MIMO technology [2, 3, 4]. These advantages are a result of equipping the base station (BS) with hundreds or thousands of antennas, which enables fine-grained beamforming to serve tens of user equipments (UEs) in the same time-frequency resource. However, the large number of antenna elements and radio frequency (RF) chains at the BS results in a significant increase in hardware complexity, system costs, and circuit power consumption. Furthermore, massive MU-MIMO requires high interconnect and chip input/output (I/O) bandwidth between the baseband-processing unit at the BS and the radio units [5, 6]. As a consequence, a successful deployment of this technology in 5G wireless systems requires novel design approaches that jointly reduce system costs, power consumption, and interconnect bandwidth without degrading the spectral efficiency and link reliability.

I-A Massive MU-MIMO with 1-bit DACs

We consider the massive MU-MIMO downlink in which the BS is equipped with 1-bit digital-to-analog converters (DACs) and transmits data to multiple UEs in the same time-frequency resource. In traditional multi-antenna BSs, each RF port is connected to a pair of high-resolution DACs (e.g., with 10-bit precision). Scaling such architectures to massive MIMO BSs, with hundreds or thousands of antennas would result in prohibitively high power consumption and system costs. The deployment of 1-bit DACs at the BS would mitigate this problem. In addition, the use of 1-bit DACs enables one to lower the linearity and noise requirements of the surrounding RF circuitry, which has the potential to additionally reduce the circuit power consumption. Another benefit of using 1-bit DACs is the fact that lowering their resolution also reduces the interconnect bandwidth between the baseband-processing unit and the radio unit, as only one bit per sample is required by each DAC. This aspect is of practical relevance for deployment scenarios in which these two units are not co-located [5, 6].

The key challenges of 1-bit massive MU-MIMO systems are to maintain high spectral efficiency and reliability.

The work in [7] demonstrates that the performance degradation caused by 1-bit DACs in the downlink diminishes as the number of BS antennas increases. Furthermore, as shown in [7, 8, 9, 10, 1], the use of 1-bit DACs in the downlink enables reliable data transmission if sophisticated precoding algorithms that simultaneously mitigate multi-user interference (MUI) and quantization artifacts are used. While conventional linear precoding methods, such as zero-forcing (ZF) or minimum mean-squared error (MMSE) precoding followed by quantization, require low computational complexity [11, 12, 13, 14], more sophisticated, nonlinear methods are necessary to enable reliable communication at high spectral efficiency. Such precoding methods, however, typically require high computational complexity. As a consequence, a successful deployment of 1-bit massive MU-MIMO calls for the design of novel and efficient precoding algorithms that can be implemented in hardware and reliably achieve high throughput at low power consumption.

I-B Contributions

In this paper, we develop novel, computationally efficient precoding algorithms for 1-bit massive MU-MIMO systems and corresponding very-large scale integrated (VLSI) designs. Our main contributions can be summarized as follows:

• •

We use biconvex relaxation (BCR) [15] to design a nonlinear 1-bit precoding algorithm. Our algorithm, referred to as C1PO (short for biConvex 1-bit PrecOding), enables reliable, high-rate downlink transmission in 1-bit massive MU-MIMO systems for medium-sized antenna arrays.

• •

We propose a scalable and low-complexity algorithm variant, referred to as C2PO, which enables high-performance 1-bit precoding for massive MU-MIMO systems with hundreds or thousands of antenna elements.

• •

For C1PO and C2PO, which both solve nonconvex problems, we provide analytical convergence guarantees.

• •

We develop two massively parallel VLSI architectures that implement C1PO and C2PO, and achieve high throughputs in a hardware-efficient manner. Our architectures support various BS and UE antenna configurations.

• •

We present reference designs on a Xilinx Virtex-7 field-programmable gate array (FPGA) for various antenna configurations that demonstrate the efficacy of our algorithms and VLSI architectures.

• •

We compare our designs to a baseline precoder that uses maximum ratio transmission (MRT) followed by quantization (MRT-Q), a method that achieves high hardware-efficiency at the cost of poor error-rate performance.

• •

We study the trade-offs between error-rate performance and hardware efficiency (in terms of throughput per area) for the proposed FPGA designs.

Our results demonstrate the practical feasibility of 1-bit precoding in massive MU-MIMO systems, supporting reliable and high-rate downlink data transmission.

I-C Relevant Prior Art

A number of papers have studied the use of low-resolution analog-to-digital converters (ADCs) for the massive MU-MIMO uplink (UEs transmit data to the BS) with particular focus on the 1-bit case; see, e.g., [16, 17, 18, 19, 20] and the references therein. All these results have shown that the use of 1-bit ADCs is sufficient for reliable low-rate uplink transmission and that $4$ to $6$ bits are sufficient to close the gap to the infinite-precision case in most scenarios. In contrast to the uplink, the quantized downlink has gained attention only recently. Precoding in the downlink with 1-bit DACs is a more challenging problem as both MUI and quantization artifacts must be mitigated simultaneously. The results in [11, 12, 13, 14] have shown that so-called linear-quantized precoders, which perform traditional linear precoding followed by quantization, enable reliable downlink transmission for very large BS antenna arrays in the high signal-to-noise ratio (SNR) regime, even for systems that use 1-bit DACs. More sophisticated, nonlinear precoding algorithms have been proposed only recently in [7, 8, 9, 10, 1] and significantly outperform linear-quantized methods in the presence of 1-bit DACs. The computational complexity of these algorithms, however, is typically high, which prevents an efficient implementation in practical systems. In contrast to these precoding methods, we propose two novel nonlinear precoding algorithms and VLSI designs that achieve high throughput in a hardware-efficient manner.

While a large number of VLSI designs for data detection in the massive MU-MIMO uplink have been proposed in the literature (see, e.g., [21, 22, 23, 24] and references therein), only a handful of precoder designs for multi-antenna downlink systems exist [25, 26, 27, 28, 5]. Reference [25] proposes a VLSI design for vector-perturbation precoding in small-scale MIMO systems with high-precision DACs. The papers [26] and [27] discuss hardware implementations for approximate linear and ZF/MRT-based precoding, respectively, for massive MU-MIMO systems with high-precision DACs. Unfortunately, both of these publications do not provide detailed FPGA implementation results. Reference [28] describes an application specific integrated circuit (ASIC) design of a ZF precoder; reference [5] presents a decentralized ZF precoder on a graphics processing unit (GPU) cluster. Both of these precoders are, however, designed for high-precision DACs and not for 1-bit massive MU-MIMO systems. Hence, to the best of our knowledge, the VLSI designs proposed in this paper are the first hardware implementations reported in the open literature that are suitable for precoding in the 1-bit massive MU-MIMO downlink.

I-D Notation

Lowercase and uppercase boldface letters designate column vectors and matrices, respectively. For a matrix $\mathbf{A}$, we denote its transpose, Hermitian transpose, complex conjugate, and matrix $\ell_{2}$-norm by $\mathbf{A}^{T}$, $\mathbf{A}^{H}$, $\mathbf{A}^{*}$, and $\|\mathbf{A}\|_{2,2}$, respectively; the entry on the $k$th row and on the $\ell$th column of $\mathbf{A}$ is $[\mathbf{A}]_{k,\ell}$. The $M\times M$ identity matrix is denoted by $\mathbf{I}_{M}$ and the $M\times N$ all-zeros matrix is denoted by $\mathbf{0}_{M\times N}$. For a vector $\mathbf{a}$, the $k$th entry is $[\mathbf{a}]_{k}$ and we use $\mathopen{}\left\lVert\mathbf{a}\right\rVert_{2}$ to denote the $\ell_{2}$-norm of the vector $\mathbf{a}$. The real and imaginary parts of a complex vector $\mathbf{a}$ are $\Re\{\mathbf{a}\}$ and $\Im\{\mathbf{a}\}$, respectively. The signum function $\text{sgn}(\cdot)$ is defined as $\text{sgn}(a)=+1$ for $a\geq 0$ and $\text{sgn}(a)=-1$ for $a<0$ and is applied entry-wise to vectors. The multivariate complex-valued circularly-symmetric Gaussian probability density function (PDF) with covariance matrix $\mathbf{K}$ is denoted by $\mathcal{C}\mathcal{N}(\mathbf{0},\mathbf{K})$. We use $\operatorname{\mathbb{E}}_{\mathbf{x}}\mathopen{}\left[\cdot\right]$ to denote expectation with respect to the random vector $\mathbf{x}$.

I-E Paper Outline

The rest of this paper is organized as follows. In \frefsec:systemmodel, we introduce the system model and formulate the precoding problem for systems with 1-bit DACs. In \frefsec:biconvexrelaxation, we propose two new 1-bit precoding algorithms, namely C1PO and C2PO. In \frefsec:architecture1 and \frefsec:architecture2, we detail our VLSI architectures for C1PO and C2PO, respectively. In \frefsec:results, we show numerical simulations, reference FPGA implementation results, and a comparison with an MRT-based baseline precoder. We conclude the paper in \frefsec:conclusions. All proofs are relegated to Appendices A and B.

II System Model and 1-bit Precoding

We start by introducing the downlink system model and then provide the necessary details about optimal precoding in 1-bit massive MU-MIMO systems.

II-A Downlink System Model

We focus on the downlink of a single-cell, narrowband massive MU-MIMO system as illustrated in \freffig:system_overview. The system consists of a $B$-antenna BS that serves $U\leq B$ single-antenna111For simplicity, we focus on single-antenna UEs; the model can easily be expanded to support multi-antenna UEs. UEs simultaneously and in the same frequency band. We use the standard input-output relation $\mathbf{y}=\mathbf{H}\mathbf{x}+\mathbf{n}$ to model the narrowband downlink channel [2]. Here, the vector $\mathbf{y}=[y_{1},\,\dots,\,y_{U}]^{T}$ contains the received signals at all UEs, where $y_{u}\in\mathbb{C}$ is the signal received at the $u$th UE. The matrix $\mathbf{H}\in\mathbb{C}^{U\times B}$ represents the downlink channel. The so-called precoded vector is denoted by $\mathbf{x}\in\mathcal{X}^{B}$, where $\mathcal{X}$ represents the transmit alphabet; this set coincides with the set $\mathbb{C}$ of complex numbers in the case of infinite-precision DACs. In 1-bit massive MU-MIMO systems, the in-phase and quadrature components are generated separately using a pair of 1-bit DACs running at Nyquist rate and hence, the per-antenna quaternary transmit alphabet is $\mathcal{X}=\{\pm\ell\pm j\ell\}$ for a given (and fixed) $\ell>0$ that determines the transmit power. The vector $\mathbf{n}\in\mathbb{C}^{U}$ models i.i.d. circularly-symmetric complex Gaussian noise with variance $N_{0}$ per complex entry, i.e., $n_{u}\sim\mathcal{CN}(0,N_{0})$, for $u=1,\,\dots,\,U$. In what follows, we assume that the realization of the channel matrix $\mathbf{H}$ and the noise variance $N_{0}$ are perfectly known at the BS.222Knowledge of $\mathbf{H}$ is typically acquired via training in the uplink in a time-division duplexing system [2]. As discussed in [7], channel estimation errors yield only a small performance loss. Knowledge of the noise variance $N_{0}$ at the BS can be obtained by explicit feedback from the UEs to the BS.

II-B Precoding Basics

The main purpose of precoding is to transmit the constellation points $s_{u}\in\mathcal{O}$ to each UE $u=1,\ldots,U$, where $\mathcal{O}$ is the constellation set (e.g., 16-QAM). The BS uses the available channel state information (CSI) to precode the symbol vector $\mathbf{s}=[s_{1},\,\dots,\,s_{U}]^{T}$ into the precoded vector $\mathbf{x}\in\mathcal{X}^{B}$. Throughout the paper, we assume that the precoded vector $\mathbf{x}$ must satisfy an instantaneous power constraint $\mathopen{}\left\lVert\mathbf{x}\right\rVert^{2}_{2}=P$; this leads to \mathcal{X}=\big{\{}\pm\ell\pm j\ell\big{\}} with $\ell=\sqrt{P/(2B)}$.

Coherent transmission of data using multiple BS antennas leads to an array gain, which depends on the realization of the fading channel and the precoding method. As in [7, 8], we assume that the $u$th UE is able to rescale its received signals $y_{u}$ by a factor333In contrast to references [7, 8], which assumed real-valued factors $\beta_{u}$, $u=1,\ldots,U$, we allow these factors to be complex-valued. $\beta_{u}\in\mathbb{C}$ in order to compute an estimate $\hat{s}_{u}=\beta_{u}y_{u}$ for $u=1,\ldots,U$ of the transmitted symbol $s_{u}\in\mathcal{O}$.

Since the UEs cannot perform joint processing to recover the transmitted data, precoding must simultaneously reduce MUI and increase signal power at all UEs [29]. To accomplish these goals, there exist multiple formulations of this optimization problem based on different performance metrics, e.g., sum-rate throughput or error-rate (see [30] for a survey). As in [7, 8], we will focus exclusively on precoders that minimize the mean-squared error (MSE) between the estimated symbol vector \hat{\mathbf{s}}=\big{[}\hat{s}_{1},\ldots,\hat{s}_{U}\big{]}^{T}=\beta\mathbf{y} and the transmitted symbol vector $\mathbf{s}$ given by

[TABLE]

where we restrict ourselves to the case in which the precoder results in the same precoding factor $\beta$ for all UEs. Hence, in the remainder of this paper we shall assume that $\beta_{u}=\beta$ for $u=1,\ldots,U$. With this assumption, the MSE after precoding will roughly be the same for all UEs, which guarantees a certain degree of fairness among the UEs; see [7] for more details. In [8] it is shown that the UEs are able to accurately estimate the precoding factor $\beta$ using pilot-based transmission in block-fading scenarios.

In the infinite-precision case, an MSE-optimal linear precoder multiplies the symbol vector $\mathbf{s}$ with a precoding matrix $\mathbf{P}\in\mathbb{C}^{B\times U}$ so that (1) is minimized on average over all possible transmit vectors $\mathbf{s}$ subject to the power constraint. This problem, which has been studied extensively for the case of infinite-precision DACs [31, 32], enables the design of low-complexity linear precoding algorithms [2].

II-C MSE-Optimal 1-bit Precoding Problem

In the 1-bit case, linear-quantized precoders perform first linear precoding and then quantize the result to the finite transmit set $\mathcal{X}^{B}$ as

[TABLE]

for a given precoding matrix $\mathbf{P}$. Linear-quantized precoders can be analyzed theoretically and typically exhibit low complexity [7]. However, as recently shown in [7, 8, 9, 10, 1], significant performance improvements can be obtained by using sophisticated nonlinear precoding methods.

One way to design such nonlinear precoders is to solve the following MSE-optimal 1-bit precoding problem (OPP), which simultaneously finds the optimal precoding vector $\mathbf{x}^{\text{{\text{OPP}}}}$ and the associated precoding factor $\beta^{\text{{\text{OPP}}}}$:

[TABLE]

We emphasize that for a fixed value of $\beta$, the problem (OPP) is a closest vector problem that is known to be NP-hard [33, 34, 35]; this implies that there exists no known algorithm to solve it efficiently for large values of $B$. In [7, 8], approximate methods for solving (OPP) using convex relaxation have been proposed, such as the squared-infinity norm Douglas-Rachford splitting (SQUID) algorithm. Such relaxation-based methods, however, still require high computational complexity, which prevents their deployment in practical systems.

III 1-bit Precoding via Biconvex Relaxation

Since the problem (OPP) is of combinatorial nature, a brute-force search for a solution is intractable in massive MU-MIMO systems with hundreds of BS antennas. We next propose two nonlinear precoding algorithms that yield approximate but accurate solutions at low computational complexity.

III-A Approximating (OPP)

To solve (OPP) efficiently, we use the BCR framework put forward in [15], which was initially proposed for solving large semidefinite programs that appear in computer vision. In order to use this framework, we first simplify the objective function of (OPP) by assuming that $N_{0}\to 0$, i.e., we assume that the system operates in the high-SNR regime. Note that we make this assumption solely for the purpose of deriving computationally efficient algorithms; we show in \frefsec:simresults that our algorithms also work well in the low-SNR regime. As in [1, Eq. 3], we take a leap of faith with the approximation

[TABLE]

This approximation can be justified by noting that if we can find a precoded vector $\mathbf{x}\in\mathcal{X}^{B}$ for which $\mathbf{s}=\mathbf{H}\mathbf{x}$, then both problems in \frefeq:dummy are indeed equivalent. These approximations allow us to rewrite (OPP) as follows:

[TABLE]

We next get rid of the parameter $\alpha$ in (OPP∗). For a fixed $\mathbf{x}$, the optimal parameter $\hat{\alpha}(\mathbf{x})$ that minimizes the objective function of (OPP∗) is given by

[TABLE]

By inserting $\hat{\alpha}(\mathbf{x})$ into the objective function of (OPP∗), we obtain

[TABLE]

with

[TABLE]

where the matrix $\mathbf{Q}\in\mathbb{C}^{U\times U}$ is a projection onto the orthogonal complement of the space spanned by the symbol vector $\mathbf{s}$. Using \frefeq:objectivefunction, the problem (OPP∗) can be simplified to

[TABLE]

which remains to be a closest vector problem. Nevertheless, the specific form of (OPP*∗∗*) enables us to use BCR to efficiently compute approximate but accurate solutions.

III-B Biconvex Relaxation (BCR)

To solve (OPP*∗∗) using BCR, we first introduce a copy $\mathbf{z}$ of the vector $\mathbf{x}$, and replace (OPP∗∗*) with the approximation

[TABLE]

where $\gamma>0$ is a (fixed) regularization parameter. We next relax the nonconvex alphabet constraint $\mathbf{x}\in\mathcal{X}^{B}$ to its convex envelope given by

[TABLE]

This relaxation allows us to convexify the precoding problem as follows:

[TABLE]

which enables the design of algorithms that converge quickly. Unfortunately, solving this optimization problem yields, in general, the all-zeros vector, i.e., $\mathbf{x}=\mathbf{0}_{B\times 1}$. One of the key ideas of BCR is to force the solution of this new problem to satisfy the constraints in (III-B) with equality. This can be accomplished by including a nonconvex regularization term in the objective that promotes large values of $\mathbf{x}.$ As suggested in [15], we use a negative $\ell_{2}$-norm term to obtain the following biconvex relaxation optimization problem:

[TABLE]

where $\delta>0$ is a (fixed) regularization parameter. If $\delta<\gamma,$ then the formulation (BCR∗) is biconvex (i.e., the minimization with respect to $\mathbf{x}$ is convex when $\mathbf{z}$ is fixed, and vice versa). Robust parameter choices are $\gamma=\|\mathbf{A}^{H}\mathbf{A}\|_{2,2}$ and $\gamma/\delta=2$; see [15] for more details. In practice, we use numerical simulations to tune the parameters $\gamma$ and $\delta$ to further improve the empirical performance of our algorithms.

III-C C1PO: biConvex 1-bit PrecOding

We have noted above that the (BCR∗) problem is biconvex, meaning that minimization with respect to $\mathbf{x}$ alone (with $\mathbf{z}$ fixed) or $\mathbf{z}$ alone (with $\mathbf{x}$ fixed) is convex. Hence, as done in [15], we can solve the (BCR∗) problem approximately using alternating minimization. Since the problem is nonconvex, initialization critically affects the performance of our algorithm. We initialize our algorithm with the MRT precoded vector $\mathbf{x}^{(1)}=\mathbf{H}^{H}\mathbf{s}$, which yields excellent performance in practice and can be computed efficiently. Then, we solve for $\mathbf{z}$ while holding $\mathbf{x}$ fixed; afterwards, we solve for $\mathbf{x}$ while holding $\mathbf{z}$ fixed. Specifically, we repeat the following procedure:

[TABLE]

for $t=1,2,\ldots,t_{\text{max}}$, where $t_{\text{max}}$ is the maximum number of iterations. Both steps are convex optimization problems that can be solved efficiently in closed form. Hence, the above iterative procedure reduces to the following simple algorithm, which we call C1PO (short for biConvex 1-bit PrecOding).

Algorithm 1** (C1PO).**

Set $\mathbf{A}$ as in \frefeq:Amatrix, initialize $\mathbf{x}^{(1)}=\mathbf{H}^{H}\mathbf{s}$, and fix the parameters $\delta$ and $\gamma$ so that $0<\delta<\gamma$. Then, for every iteration $t=1,2,\ldots,t_{\text{max}},$ compute

[TABLE]

Here, the expansion-reprojection operator $\mathrm{proj}(\cdot)$ is

[TABLE]

and is applied element-wise to the vector $\mathbf{z}^{(t+1)}$. In the last iteration $t_{\text{max}}$, the output $\mathbf{x}^{(t_{\text{max}}+1)}$ of C1PO is quantized to the quaternary alphabet \mathcal{X}=\big{\{}\pm\ell\pm j\ell\big{\}} with $\ell=\sqrt{P/(2B)}$ as follows:

[TABLE]

Because C1PO decreases the objective function (BCR∗) on every variable update, and the objective is bounded from below, the objective values corresponding to the iterates $\{\mathbf{x}^{(t)},\mathbf{z}^{(t)}\}$ form a convergent sequence. However, by exploiting the biconvex structure of our problem, we can prove the following stronger result; the proof is given in \frefapp:c1po.

Theorem 1**.**

Any limit point of the sequence $\{\mathbf{x}^{(t)},\mathbf{z}^{(t)}\}$ generated by C1PO is a stationary point of (BCR∗).

The main computations performed by C1PO in \frefalg:C1PO are (i) the $B\times B$ matrix inversion $\mathbf{G}=(\mathbf{I}_{B}+\gamma^{-1}\mathbf{A}^{H}\mathbf{A})^{-1}$, which can be computed once during a preprocessing stage, and (ii) the per-iteration matrix-vector multiplication $\mathbf{z}^{(t+1)}=\mathbf{G}\mathbf{x}^{(t)}$ in step \frefeq:step1; the complexity of the projection in step \frefeq:step2 is negligible. Unfortunately, the complexity of the matrix inversion, evaluated in terms of operations,444For simplicity, we count the number of complex-valued multiplications to characterize the operation count. scales roughly with $B^{3}$ and the complexity of the per-iteration matrix-vector product with $B^{2}$. Both of these tasks are particularly inefficient for massive MU-MIMO systems with a large number of BS antennas. Therefore, we next propose an algorithmic variant that avoids both of these issues and whose complexity scales more favorably with the number of BS antennas.

III-D Fast Algorithm for Very-Large Systems: C2PO

To obtain our alternative algorithm, we start from the BCR formulation in (BCR∗) but rather than introducing the auxiliary variable $\mathbf{z}$, we attempt to directly solve the following nonconvex optimization problem:555To simplify notation, we have divided both terms in the objective function by a factor of two; this scaling does not affect the result.

[TABLE]

We use forward-backward splitting (FBS) [36, 37, 38], a computationally efficient method to solve large convex problems. Since the problem in \frefeq:directform is nonconvex, FBS is not guaranteed to converge to the optimal solution. Nevertheless, as shown in \frefsec:results, the proposed algorithm performs well in practice.

FBS is an efficient iterative procedure to solve convex optimization problems of the form

[TABLE]

where the function $f$ is smooth and convex, and the function $g$ is convex but not necessarily smooth or bounded. FBS consists of the following iteration [37, 36]:

[TABLE]

for $t=1,2,\ldots,t_{\text{max}}$ or until convergence. Here, $\nabla f(\mathbf{x})$ is the gradient of the smooth function $f$, and the so-called proximal operator for the function $g$ is defined as follows [39]:

[TABLE]

The sequence $\{\tau^{(t)}>0\}$ contains suitably chosen step-size parameters. For the problem (9), we show below that FBS monotonically decreases the objective (9) for any constant step size that satisfies $\tau^{(t)}=\tau<\|\mathbf{A}^{H}\mathbf{A}\|^{-1}_{2,2}.$

In order to approximately solve \frefeq:directform using FBS, we set

[TABLE]

where $\chi$ is a characteristic function that is zero if the condition $\mathbf{x}\in\mathcal{B}^{B}$ is met and infinity otherwise. For this choice of $f$ and $g$, the gradient is given by $\nabla f(\mathbf{x})=\mathbf{A}^{H}\mathbf{A}\mathbf{x}$ and the proximal operator is given by the expansion-reprojection operation

[TABLE]

which is valid for $\tau\delta<1$ and applied element-wise to vectors. By using FBS with the above-mentioned ingredients, we obtain the following simple algorithm, which we call C2PO.

Algorithm 2** (C2PO).**

Set $\mathbf{A}$ as in \frefeq:Amatrix. Initialize $\mathbf{x}^{(1)}=\mathbf{H}^{H}\mathbf{s}$ and fix the parameters $\delta$ and $\tau$ so that $\tau\delta<1$. Then, for every iteration $t=1,2,\ldots,t_{\text{max}}$ compute:

[TABLE]

Here, the $\mathrm{prox}_{g}$ operator is given in \frefeq:proxg and is applied element-wise to the vector $\mathbf{z}^{(t+1)}$. In the last iteration $t_{\text{max}}$, the output $\mathbf{x}^{(t_{\text{max}}+1)}$ of C2PO is quantized to the quaternary alphabet $\mathcal{X}$ as in \frefeq:outputquantization.

The following result shows that C2PO is well behaved, provided that the step size is chosen appropriately; the proof is given in \frefapp:c2po.

Theorem 2**.**

Suppose the step size used in C2PO satisfies $\tau<\|\mathbf{A}^{H}\mathbf{A}\|^{-1}_{2,2},$ and $\tau\delta<1.$ Then, C2PO decreases the objective (9) monotonically, and any limit point of the iterates $\{\mathbf{x}^{(t)}\}$ is a stationary point.

The most complex operation of C2PO (\frefalg:C2PO) is the matrix-vector multiplication in step \frefeq:faststep1. In contrast to C1PO (\frefalg:C1PO), however, this step requires a minimal amount of preprocessing and can be computed efficiently, especially for large BS antenna arrays. To see this, we rewrite $\mathbf{A}^{H}\mathbf{A}$, where $\mathbf{A}$ was given in \frefeq:Amatrix, as follows:

[TABLE]

Here, $\mathbf{v}=\mathbf{H}^{H}\mathbf{s}/\|\mathbf{s}\|_{2}$ is a normalized version of the MRT vector; the augmented matrices $\overline{\mathbf{H}}=[\mathbf{H};\mathbf{v}^{H}]$ and $\overline{\mathbf{H}}^{\Upsilon}=[\mathbf{H}^{H},-\mathbf{v}]$ are of dimension $(U+1)\times B$ and $B\times(U+1)$, respectively. With these definitions, we can now simplify step \frefeq:faststep1 to

[TABLE]

where we first compute $\mathbf{w}=\overline{\mathbf{H}}(\tau\mathbf{x}^{(t)})$, then $\mathbf{w}^{\prime}={\overline{\mathbf{H}}^{\Upsilon}}\mathbf{w}$, and finally $\mathbf{z}^{(t+1)}=\mathbf{x}^{(t)}-\mathbf{w}^{\prime}$. As a result, we conclude that each iteration of C2PO requires only two matrix-vector products with a cost of roughly $2B(U+1)$ operations (in contrast to $B^{2}$ operations for C1PO). In addition, the preprocessing stage of this algorithm only needs to compute the normalized MRT vector $\mathbf{v}$, which requires roughly $BU$ operations (in contrast to $B^{3}$ operations for C1PO). Hence, the complexity of C2PO can be significantly lower than that of C1PO, especially since the antenna configurations of typical massive MU-MIMO systems satisfy $U\ll B$. As we will show in \frefsec:results, the hardware efficiency of C2PO is superior to that of C1PO for large BS antenna arrays and the error-rate performance is comparable.

III-E Alternative Derivation of C2PO

It is interesting to note that there is a strong connection between \frefalg:C1PO and \frefalg:C2PO. In fact, one can obtain C2PO directly from C1PO using the following well-known series expansion. Let $\|\mathbf{A}^{H}\mathbf{A}\|_{2,2}<\gamma$. Then, we have the following Neumann series expansion [40]:

[TABLE]

As suggested in [21], we can approximate the inverse by truncating the series to the two first terms:

[TABLE]

By using this approximation in step \frefeq:step1 of \frefalg:C1PO, we immediately obtain \frefalg:C2PO after setting $\gamma^{-1}=\tau$. Note that the Neumann series expansion is only convergent for $\|\mathbf{A}^{H}\mathbf{A}\|_{2,2}<\gamma,$ which corresponds to the step size restriction $\tau<\|\mathbf{A}{{}^{H}}\mathbf{A}\|^{-1}_{2,2}$; this is exactly the same step size requirement as in \frefthm:c2po.

IV VLSI Design for C1PO

We now present a high-throughput VLSI architecture for C1PO as in \frefalg:C1PO. We then discuss the key optimizations performed in our FPGA implementation.

IV-A Architecture Overview

The proposed VLSI architecture that implements C1PO as detailed in \frefalg:C1PO is shown in \freffig:C1POarch. Our architecture consists of a linear array of $B$ identical processing elements (PEs) that share a common control unit. The PEs essentially compute the complex-valued matrix-vector product in (6), using a variant of Cannon’s algorithm [41], followed by the projection operation in (7). Each PE $b=1,2,\ldots,B$ consists of three main building blocks: (i) a $\mathbf{g}_{b}$-memory, (ii) a complex-valued multiply-accumulate (MAC) unit, and (iii) a projection unit. For the $b$th PE, the $\mathbf{g}_{b}$-memory stores the $b$th row of the matrix $\mathbf{G}=(\mathbf{I}_{B}+\gamma^{-1}\mathbf{A}^{H}\mathbf{A})^{-1}$, which we assume was computed during a separate preprocessing stage.

As mentioned in \frefsec:sub_biconvexrelaxation, simulations are used to tune the parameter $\gamma$ in order to improve the error-rate performance; the optimal value of $\gamma$ depends on the antenna configuration. The complex-valued MAC unit is used by each PE to sequentially compute an entry of the output vector $\mathbf{z}^{(t+1)}$ on line (6), while the entries of the vector $\mathbf{x}^{(t)}$ are exchanged between the PEs in a cyclic fashion; this is done to avoid an architecture with a centralized $\mathbf{x}^{(t)}$ memory that would suffer from a high fan-out because the memory’s output has to be distributed to all the PEs. The projection unit implements the expansion-reprojection operator $\mathrm{proj}(\cdot)$ on line (7) in a hardware-friendly manner. The outputs of the projection unit are also used to generate the quaternary outputs of the 1-bit precoder; to this end, each PE simply takes the sign bits of the complex-valued output vector $\mathbf{x}^{(t+1)}$.

IV-B Architecture Operation

We now detail the (rather technical) operation of the C1PO architecture illustrated in \freffig:C1POarch. In the first iteration (i.e., at $t=1$), each PE $b$ is initialized with the $b$th entry of the vector $\mathbf{x}^{(1)}$. Furthermore, the entries of the $\mathbf{g}_{b}$-memory are stored so that the first memory address corresponds to $[\mathbf{G}]_{b,b}$, the second address to $[\mathbf{G}]_{b,b+1}$, and so forth (addresses wrap around).

In the first clock cycle, each PE $b$ computes $[\mathbf{G}]_{b,b}[\mathbf{x}^{(t)}]_{b}$ and the result is stored in the accumulator. As shown on the left side of \freffig:C1POarch, in the same clock cycle the $b$th PE passes the value $[\mathbf{x}^{(t)}]_{b}$ to PE $(b-1)$, while it receives the value $[\mathbf{x}^{(t)}]_{b+1}$ from PE $(b+1)$; PE $1$ passes its value to PE $B$. In the second clock cycle, since the exchange operation made the element $[\mathbf{x}^{(t)}]_{b+1}$ available at PE $b$, each PE computes $[\mathbf{G}]_{b,b+1}\cdot[\mathbf{x}^{(t)}]_{b+1}$ and uses the accumulator to add it to the result of the previous cycle. Once again, in the same clock cycle, the $b$th PE passes the $\mathbf{x}^{(t)}$ entry that is currently being multiplied on its MAC unit to PE $(b-1)$; PE $1$ passes its value to PE $B$. Consequently, in the third clock cycle, the $b$th PE will use the values $[\mathbf{G}]_{b,b+2}$ and $[\mathbf{x}^{(t)}]_{b+2}$ to continue performing MAC operations. By repeating this procedure $B$ times, each entry of the vector $[\mathbf{x}^{(t)}]$ circulates through all the PEs exactly once, enabling each PE $b=1,2,\ldots,B$ to compute $[\mathbf{z}^{(t+1)}]_{b}$. Thus, the matrix-vector multiplication on line (6) is completed. Since the complex-valued MAC unit contains three pipeline stages, two clock cycles are required to flush the pipeline. Hence, the matrix-vector operation requires a latency of $B+2$ clock cycles. After $B+2$ clock cycles, the vector $\mathbf{z}^{(t+1)}$ is available at the outputs of the MAC units.

In the subsequent clock cycle, each PE projects their respective entry of the vector $\mathbf{z}^{(t+1)}$. According to our simulation results, the choice $\gamma/\delta=5$, which implies that $\gamma/(\gamma-\delta)=1.25$, works well for all the considered antenna configurations. Furthermore, to reduce the hardware complexity, we assume $P=2B$ so that the clipping threshold of the expansion-reprojection operator $\mathrm{proj}(\cdot)$ is $1$. As a result, the $\mathrm{proj}(\cdot)$ operator in (7) is implemented by applying the following operations independently to the real and imaginary parts of $[\mathbf{z}^{(t+1)}]_{b}$: We multiply the real (or imaginary) part of $[\mathbf{z}^{(t+1)}]_{b}$ by $1.25$; this is accomplished by adding the $[\mathbf{z}^{(t+1)}]_{b}$ value with a $2\times$ right-shifted version of itself. At the same time, the real (or imaginary) part of $[\mathbf{z}^{(t+1)}]_{b}$ is compared to $-0.8$ and $+0.8$. If the real (or imaginary) part of $[\mathbf{z}^{(t+1)}]_{b}$ is between these two numbers, then the projection unit outputs $1.25\cdot[\mathbf{z}^{(t+1)}]_{b}$. If it is smaller than $-0.8$, then the projection unit generates $-1$; if it is larger than $+0.8$, it generates $+1$. The result from this projection is stored as the next iterate $[\mathbf{x}^{(t+1)}]_{b}$ in the input register of the complex-valued MAC unit, which completes one C1PO iteration. Since the projection requires an additional clock cycle, one C1PO iteration is completed in exactly $B+3$ clock cycles.

IV-C FPGA Implementation Details

To minimize the FPGA implementation complexity of C1PO, we exclusively use fixed-point arithmetic; see \frefsec:simresults for the fixed-point error-rate performance of C1PO. To represent the entries of the vector $\mathbf{x}^{(t)}$, we use $12$-bit signed fixed-point values with $5$ fraction bits. The entries of the $\mathbf{G}$ matrix are represented using $10$ bits with $9$ fraction bits, and we use FPGA look-up tables (LUTs) as a distributed RAM to store these values. The complex-valued MAC unit uses $18$ bits with $11$ fraction bits; the projection unit uses $15$ bits with $8$ fraction bits. In our C1PO design, all adders and multipliers do not saturate, but wrap around; number resizing uses truncation. All complex-valued multipliers consist of four real-valued multipliers and two adders; we use the built-in DSP48 units for these operations.

V VLSI Design for C2PO

We now present a high-throughput VLSI architecture for C2PO as in \frefalg:C2PO. We then discuss the key optimizations used in our FPGA implementation.

V-A Architecture Overview

The proposed VLSI architecture that implements C2PO (\frefalg:C2PO) is shown in \freffig:C2POarch. In what follows, we assume that $B$ is a multiple of $U$ and $B\gg U$. Our architecture consists of $B/U$ linear arrays; each array consists of $U+1$ PEs and a control unit. The architecture divides the operation in (14) into two separate matrix-vector products: (i) $\mathbf{w}=\overline{\mathbf{H}}(\tau\mathbf{x}^{(t)})$ and (ii) $\mathbf{w}^{\prime}={\overline{\mathbf{H}}^{\Upsilon}}\mathbf{w}$; see also the discussion at the end of \frefsec:C2POalgorithm. We assume that $\overline{\mathbf{H}}$ was computed in a separate preprocessing stage. Note that for the first matrix-vector product, the matrix $\overline{\mathbf{H}}$ has more columns ($B$) than rows ($U+1$); for the second matrix-vector product, the matrix ${\overline{\mathbf{H}}^{\Upsilon}}$ has more rows than columns. Therefore, we will refer to the first matrix-vector product as the wide product, while the second one will be identified as the tall product. The final subtraction required to compute $\mathbf{z}^{(t+1)}=\mathbf{x}^{(t)}-\mathbf{w}^{\prime}$ in \frefeq:simplefaststep1 is incorporated into the tall-product operation; see \frefsec:arch2operation for more details.

To perform both the wide and tall products in a single computation unit, the matrix $\overline{\mathbf{H}}$ is divided into $B/U$ sub-matrices $\widetilde{\mathbf{H}}_{w}\in\mathbb{C}^{(U+1)\times U}$, $w=1,2,\ldots,B/U$, so that \overline{\mathbf{H}}=\big{[}\widetilde{\mathbf{H}}_{1},\widetilde{\mathbf{H}}_{2},\ldots,\widetilde{\mathbf{H}}_{B/U}\big{]}. In the same way, the matrix $\overline{\mathbf{H}}^{\Upsilon}$ is divided into $B/U$ sub-matrices $\widetilde{\mathbf{H}}^{\Upsilon}_{w}\in\mathbb{C}^{U\times(U+1)}$, where $w=1,2,\ldots,B/U$, and \overline{\mathbf{H}}^{\Upsilon}=\big{[}\widetilde{\mathbf{H}}^{\Upsilon}_{1};\widetilde{\mathbf{H}}^{\Upsilon}_{2};\ldots;\widetilde{\mathbf{H}}^{\Upsilon}_{B/U}\big{]}. Note that $\widetilde{\mathbf{H}}_{w}^{H}$ and $\widetilde{\mathbf{H}}^{\Upsilon}_{w}$ are the same matrices, except for a sign flip of the last column. Analogously to the matrices case, the vector $\mathbf{x}^{(t)}$ is divided into $B/U$ sub-vectors $\tilde{\mathbf{x}}^{(t)}_{w}\in\mathbb{C}^{U}$, $w=1,2,\ldots,B/U$, so that \mathbf{x}^{(t)}=\big{[}\tilde{\mathbf{x}}^{(t)}_{1};\tilde{\mathbf{x}}^{(t)}_{2};\ldots;\tilde{\mathbf{x}}^{(t)}_{B/U}\big{]}. We next outline the architectural principle of the wide and tall products.

(i) Wide Product: Each linear array takes one sub-matrix $\widetilde{\mathbf{H}}_{w}$ and the associated sub-vector $\tilde{\mathbf{x}}^{(t)}_{w}$ as its inputs and computes their product in a sequential, column-by-column manner. This operation is analogous to that of the C1PO architecture (cf. \frefsec:arch1operation) and within each linear array, the entries of the scaled sub-vector $\tau\tilde{\mathbf{x}}^{(t)}_{w}$ are cyclically exchanged among the PEs. The resulting vectors $\mathbf{w}_{w}=\widetilde{\mathbf{H}}_{w}(\tau\tilde{\mathbf{x}}^{(t)}_{w})$ are then added to obtain

[TABLE]

which completes the wide product. Each entry $[\mathbf{w}]_{u}$ of the resulting vector $\mathbf{w}$ is then stored in PE $u$ of all linear arrays.

(ii) Tall Product: With the $\mathbf{w}$ vector available in all the linear arrays, each array now computes $U$ entries of $\mathbf{z}^{(t+1)}$ by implementing $\mathbf{z}^{(t+1)}_{w}=\tilde{\mathbf{x}}^{(t)}_{w}-{\widetilde{\mathbf{H}}^{\Upsilon}_{w}}\mathbf{w}$. Here, however, we use a sequential procedure in which the accumulated results are exchanged between PEs of the same array. This procedure is—once again—a variant of Cannon’s algorithm [41]; see \frefsec:arch2operation for a detailed explanation. As a result, each linear array can project its computed $\mathbf{z}^{(t+1)}_{w}$ entries to generate the next iterate $\tilde{\mathbf{x}}^{(t+1)}_{w}$, which are then used by the same linear array to proceed with the next iteration. The sign bits of the new vector $\mathbf{x}^{(t+1)}$ correspond to the outputs of the C2PO architecture.

As in the C1PO architecture, each PE $u=1,2,\ldots,U+1$ is formed by three main units. The first unit is an $\tilde{\mathbf{h}}^{[w]}_{u}$-memory, which stores the $u$th row of the $\widetilde{\mathbf{H}}_{w}$ sub-matrix; $\widetilde{\mathbf{H}}^{\Upsilon}_{w}$ can be derived directly from $\widetilde{\mathbf{H}}_{w}$. The second unit is a complex-valued MAC unit, which supports (i) multiplications of $a\times b$ and $a\times b^{*}$, (ii) accumulation by addition or subtraction, and (iii) initialization of the accumulator with a non-zero value. The third unit is the projection unit, which is equivalent to the one of C1PO, although it is merged with the accumulator of the MAC unit.

V-B Architecture Operation

We now provide the (rather technical) operation details of the C2PO architecture illustrated in \freffig:C2POarch. Without loss of generality, we focus our description on the $w$th linear array of PEs, which operates on the $\widetilde{\mathbf{H}}_{w}$ sub-matrix and the $\tilde{\mathbf{x}}^{(t)}_{w}$ sub-vector. In the first iteration (i.e., at $t=1$), the entry $[\tilde{\mathbf{x}}^{(1)}_{w}]_{u}$ and its scaled version $\tau[\tilde{\mathbf{x}}^{(1)}_{w}]_{u}$ are stored in PE $u=1,2,\ldots,U$ in two different registers: The value $[\tilde{\mathbf{x}}^{(1)}_{w}]_{u}$ is stored in the register labeled with “$a$” in \freffig:C2POarch, which will later be used to initialize the accumulator in the complex-valued MAC unit; the value $\tau[\tilde{\mathbf{x}}^{(1)}_{w}]_{u}$ is stored at the input register of the MAC unit labeled with “$b$.” We restrict the stepsize $\tau$ to be of the form $2^{-\alpha}$ for some fixed $\alpha\in\mathbb{N}^{+}$, which enables us to acquire $\tau[\tilde{\mathbf{x}}^{(1)}_{w}]_{u}$ from a simple arithmetic right-shifted version of $[\tilde{\mathbf{x}}^{(1)}_{w}]_{u}$; we used numerical simulations to optimize the error-rate performance by selecting an optimal value for $\tau$. The $(U+1)$th PE stores the same value $[\tilde{\mathbf{x}}^{(1)}_{w}]_{1}$ as that in PE $1$. Similar to the C1PO architecture, the entries of the $\tilde{\mathbf{h}}^{[w]}_{u}$-memory are stored so that the first memory address contains $[\widetilde{\mathbf{H}}_{w}]_{u,u}$, the second address $[\widetilde{\mathbf{H}}_{w}]_{u,u+1}$, and so forth (addresses wrap around). For the $(U+1)$th PE, the first address of the $\tilde{\mathbf{h}}^{[w]}_{U+1}$-memory contains $[\widetilde{\mathbf{H}}_{w}]_{U+1,1}$, the second address contains $[\widetilde{\mathbf{H}}_{w}]_{U+1,2}$, etc.

(i) Wide Product: In the first clock cycle, each PE $u=1,2,\ldots,U$ computes $[\widetilde{\mathbf{H}}_{w}]_{u,u}\cdot[\tau\tilde{\mathbf{x}}^{(t)}_{w}]_{u}$ and stores the result in the accumulator. The $(U+1)$th PE computes $[\widetilde{\mathbf{H}}_{w}]_{U+1,1}\cdot[\tau\tilde{\mathbf{x}}^{(t)}_{w}]_{1}$. As shown in the upper left side of \freffig:C2POarch, in the same clock cycle, the $u$th PE passes the value $[\tau\tilde{\mathbf{x}}^{(t)}_{w}]_{u}$ to PE $(u-1)$, while it receives the value $[\tau\tilde{\mathbf{x}}^{(t)}_{w}]_{u+1}$ from PE $(u+1)$; PE $1$ passes its value to PE $U$, while PE $(U+1)$ does not pass anything. In the second clock cycle, since the cyclic exchange operation made the entry $[\tau\tilde{\mathbf{x}}^{(t)}_{w}]_{u+1}$ available at PE $u$, each PE computes $[\widetilde{\mathbf{H}}_{w}]_{u,u+1}\cdot[\tau\tilde{\mathbf{x}}^{(t)}_{w}]_{u+1}$ and uses the accumulator to add it to the result of the previous cycle. The $(U+1)$th PE uses the same value $\tau\tilde{\mathbf{x}}^{(t)}_{w}$ as PE 1; hence, it can compute $[\tau\tilde{\mathbf{x}}^{(t)}_{w}]_{2}\cdot[\widetilde{\mathbf{H}}_{w}]_{U+1,2}$. Again, in the same clock cycle, the $u$th PE passes the $\tau\tilde{\mathbf{x}}^{(t)}_{w}$ entry that is currently being multiplied on its MAC unit to PE $(u-1)$; PE $1$ passes its value to PE $B$, while PE $(U+1)$ does not pass anything. Consequently, in the third clock cycle, the $u$th PE will use the values $[\widetilde{\mathbf{H}}_{w}]_{u,u+2}$ and $[\tau\tilde{\mathbf{x}}^{(t)}_{w}]_{u+2}$ to continue performing MAC operations. During this third cycle, the $(U+1)$th PE will calculate the product $[\widetilde{\mathbf{H}}_{w}]_{U+1,3}\cdot[\tau\tilde{\mathbf{x}}^{(t)}_{w}]_{3}$. By repeating this procedure $U$ times, each entry of the sub-vector $(\tau\tilde{\mathbf{x}}^{(t)}_{w})$ cycles through all the PEs exactly once, enabling the $w$th linear array of PEs to compute $\widetilde{\mathbf{H}}_{w}(\tau\tilde{\mathbf{x}}^{(t)}_{w})$. Since the complex-valued MAC unit contains three pipeline stages, two clock cycles are required to flush the pipeline. Hence, the previous matrix-vector operation has a latency of $U+2$ clock cycles. To complete the wide product, the vectors $\widetilde{\mathbf{H}}_{w}(\tau\tilde{\mathbf{x}}^{(t)}_{w})$ must be added. We use a binary adder tree with $\log_{2}(B/U)$ pipeline stages. Hence, the vector $\mathbf{w}$ is computed after $U+\log_{2}(B/U)+2$ clock cycles. The $u$th PE in each linear array stores the entry $[\mathbf{w}]_{u}$ in the MAC unit’s input registered labeled with “$b$” in \freffig:C2POarch.

(ii) Tall Product: In the next clock cycle, the computation of the tall-product starts. During the first clock cycle of the tall product computation, the PE $u=1,2,\ldots,U$ has available $[\mathbf{w}]_{u}$, as well as $[\widetilde{\mathbf{H}}_{w}]_{u,u}$, the first entry in its memory. The PE can then compute $[\widetilde{\mathbf{H}}_{w}]^{*}_{u,u}\cdot[\mathbf{w}]_{u}=[\widetilde{\mathbf{H}}^{H}_{w}]_{u,u}\cdot[\mathbf{w}]_{u}{=[\widetilde{\mathbf{H}}^{\Upsilon}_{w}]_{u,u}\cdot[\mathbf{w}]_{u}}$. Using the accumulator, this product is then subtracted from $[\tilde{\mathbf{x}}^{(t)}_{w}]_{u}$, which was stored during the initialization phase in the register labeled with “$a$” in \freffig:C2POarch. During the same clock cycle, the $u$th PE sends its accumulated result to the $(u-1)$th PE; PE $1$ sends its accumulated result to PE $U$. Also, in the same clock cycle, the $(U+1)$th PE multiplies the conjugate of the first entry of its memory with its $\mathbf{w}$ entry. In words, the product $[\widetilde{\mathbf{H}}_{w}]^{*}_{U+1,1}\cdot[\mathbf{w}]_{U+1}=[\widetilde{\mathbf{H}}_{w}^{H}]_{1,U+1}\cdot[\mathbf{w}]_{U+1}={-[\widetilde{\mathbf{H}}^{\Upsilon}_{w}]_{1,U+1}\cdot[\mathbf{w}]_{U+1}}$ is computed. The result is sent to the $U$th PE. In the following clock cycles, this result will cycle through the linear array using the same wires and registers that were previously used to transfer the $\tau\tilde{\mathbf{x}}^{(t)}_{w}$ entries. In the second clock cycle, the $(u-1)$th PE multiplies the value $[\mathbf{w}]_{u-1}$ with $[\widetilde{\mathbf{H}}_{w}]^{*}_{u-1,u}={[\widetilde{\mathbf{H}}^{\Upsilon}_{w}]_{u,u-1}}$. The product is then subtracted from the accumulated value received from the $u$th PE during the previous cycle, and the new accumulated value is passed to the $(u-2)$th PE. In the same clock cycle, the $(U+1)$th PE multiplies the value $[\mathbf{w}]_{U+1}$ with $[\widetilde{\mathbf{H}}_{w}]^{*}_{U+1,2}={-[\widetilde{\mathbf{H}}^{\Upsilon}_{w}]_{2,U+1}}$ and sends the result to the $U$th PE, so it can cycle through the linear array. Furthermore, PE $U$ passes the ${-[\widetilde{\mathbf{H}}^{\Upsilon}_{w}]_{1,U+1}}\cdot[\mathbf{w}]_{U+1}$ (previously received from the $(U+1)$th PE) to PE $(U-1)$. In the third clock cycle, the $(u-2)$th PE calculates ${[\widetilde{\mathbf{H}}^{\Upsilon}_{w}]_{u,u-2}}\cdot[\mathbf{w}]_{u-2}$, subtracts it from the accumulated result received on the second cycle from the $(u-1)$th PE and passes the result to the $(u-3)$th PE. In the same clock cycle, the $(U+1)$th PE computes ${-[\widetilde{\mathbf{H}}^{\Upsilon}_{w}]_{3,U+1}}\cdot[\mathbf{w}]_{U+1}$ and sends it to PE $U$. Meanwhile, ${-[\widetilde{\mathbf{H}}^{\Upsilon}_{w}]_{2,U+1}}\cdot[\mathbf{w}]_{U+1}$ is passed from PE $U$ to PE $(U-1)$ and ${-[\widetilde{\mathbf{H}}^{\Upsilon}_{w}]_{1,U+1}}\cdot[\mathbf{w}]_{U+1}$ is passed from PE $(U-1)$ to PE $(U-2)$. After repeating this procedure for $U$ clock cycles, each PE $u=1,2,\ldots,U$ will contain the accumulated result for the $u$th entry of $\tilde{\mathbf{x}}^{(t)}_{w}-{\widetilde{\mathbf{H}}^{\Upsilon}_{w}}\mathbf{w}$, received from the $(u+1)$th PE during the previous cycle. However, this accumulated result is missing the product ${-[\widetilde{\mathbf{H}}^{\Upsilon}_{w}]}_{u,U+1}\cdot[\mathbf{w}]_{U+1}$, which was computed and sent by the $(U+1)$th PE. Nonetheless, in the $(U+1)$th cycle of the tall product procedure, the $u$th PE receives the missing ${-[\widetilde{\mathbf{H}}^{\Upsilon}_{w}]}_{u,U+1}\cdot[\mathbf{w}]_{U+1}$ value from the $(u+1)$th PE. The received product is accumulated with the remaining data by addition instead of subtraction. Thus, the $\mathbf{z}^{(t+1)}=\mathbf{x}^{(t)}-{\overline{\mathbf{H}}^{\Upsilon}}\mathbf{w}$ entries are calculated after $U+1$ cycles. Since the complex MAC unit is used again, two additional clock cycles are required to flush its pipeline. Hence, $U+3$ cycles are used for the tall product.

Finally, in the subsequent clock cycle after the tall product is completed, the projection operator is applied in a similar fashion as for the C1PO architecture. As in the C1PO case, our simulation results show that the choice of $\tau\lambda=0.2$ (which implies that $1/(1-\tau\lambda)=1.25$) works well for all the considered antenna configurations. Therefore, the only difference between the projection units of C1PO and C2PO is that, in the C2PO architecture, the accumulator of the MAC unit is used to multiply the real and imaginary parts of each $\mathbf{z}^{(t+1)}$ entry with $1.25$, by adding each $\mathbf{z}^{(t+1)}$ entry with a $2\times$ right-shifted version of itself. The result from this projection is stored as the next iterate $[\tilde{\mathbf{x}}^{(t+1)}_{w}]_{u}$ in the two initialization locations previously mentioned, completing one C2PO iteration. Since the projection operation requires an additional clock cycle, a full C2PO iteration is completed in exactly $2U+\log_{2}(B/U)+6$ clock cycles.

V-C FPGA Implementation Details

As for the C1PO FPGA implementation, we exclusively use fixed-point arithmetic for the C2PO FPGA design; see \frefsec:simresults for the fixed-point error-rate performance of C2PO. To represent the entries of the vector $\mathbf{x}^{(t)}$, we use $12$-bit signed fixed-point values with $5$ fraction bits. For the scaled $\tau\mathbf{x}^{(t)}$ values, we use $12$-bit signed fixed-point values with $11$ fraction bits. The entries of the $\overline{\mathbf{H}}$ matrix consist of $10$ bits with $8$ fraction bits, and we use FPGA look-up tables (LUTs) as a distributed RAM to store these values. The complex-valued MAC unit uses $18$ bits with $15$ fraction bits when doing the wide product and $11$ fraction bits when doing the tall product; the projection unit uses $18$ bits with $11$ fraction bits. The adder tree uses $21$ bits with $15$ fraction bits. Identical to the C1PO implementation, all adders and multipliers do not saturate, but wrap around; number resizing uses truncation. All complex-valued multipliers are built with four real-valued multipliers and two adders; we use DSP48 units for these operations.

VI Results

We now provide error-rate performance results for massive MU-MIMO systems and show reference FPGA implementation results for C1PO and C2PO.

VI-A Simulation Results

\fref

fig:BER_ALL shows uncoded bit error-rate (BER) curves versus the normalized transmit power $\varrho=P/N_{0}$ for massive MU-MIMO dowlink systems with $U=16$ UEs for various precoding algorithms. In \freffig:BER_32x16_BPSK we consider the case of $B=32$ BS antennas with BPSK, whereas in \freffig:BER_256x16_16QAM we consider the case of $B=128$ BS antennas with 16-QAM. For both systems, we use Gray mapping, generate i.i.d. Rayleigh fading channel matrices, and average the BER over 10 000 Monte–Carlo trials. We compare ZF followed by quantization (ZF-Q), MRT followed by quantization (MRT-Q), the nonlinear SQUID algorithm proposed in [7], as well as C1PO and C2PO, for systems with 1-bit DACs. As a baseline, we also include ZF and MRT with infinite-precision DACs (denoted by “Inf. prec. ZF” and “Inf. prec. MRT”, respectively). SQUID runs $t_{\text{max}}=50$ iterations; C1PO and C2PO both run $t_{\text{max}}=24$ iterations. For all algorithms, the curves represent MATLAB floating-point performance; for C1PO and C2PO, the markers correspond to fixed-point performance of our hardware designs. Clearly, the fixed-point implementation loss of our hardware designs is negligible, i.e., less than 0.15 dB normalized transmit power $\varrho$ at $1$% uncoded BER for both considered scenarios.

For the $16\times 32$ system (we use the notation $U\times B$ to refer to a downlink scenario with $U$ users and $B$ BS antennas) with BPSK, \freffig:BER_32x16_BPSK shows that all nonlinear precoders significantly outperform the linear-quantized precoders (ZF-Q and MRT-Q), which exhibit a high error floor. Furthermore, we see that C1PO and C2PO achieve similar performance as that of SQUID. At low values of normalized transmit power $\varrho$, SQUID is marginally better, whereas C2PO achieves the best performance at high values of $\varrho$, closely followed by C1PO and SQUID. It can also be seen that, at high values of $\varrho$, 1-bit nonlinear precoders significantly outperform the error-rate performance of MRT with infinite-precision DACs.

For the $16\times 128$ system with 16-QAM, \freffig:BER_256x16_16QAM shows a similar trend, i.e., non-linear precoders significantly outperform linear-quantized precoders. SQUID outperforms C1PO and C2PO (which perform equally well) by about 0.5 dB normalized transmit power $\varrho$ at 1% BER. However, we note that the complexity (in terms of operation counts) of SQUID is more than 2$\times$ higher than that of C1PO and C2PO, and also involves the sorting of $B$ dimensional vectors which is difficult to implement efficiently in VLSI. We also observe that non-linear precoders enable reliable transmission of higher-order modulation schemes (such as 16-QAM), which is not possible with linear-quantized methods—the error-rate performance of nonlinear 1-bit precoders for higher-order modulation schemes is studied in more detail in [8]. We also see that non-linear precoders do not exhibit an error floor in the considered BER range, which is in contrast to the linear-quantized ones. We note that a detailed theoretical analysis of the error-rate performance of non-linear 1-bit precoders is an open research problem.

Remark 1**.**

Our results are limited to a narrowband downlink channel, in which we assume that the BS has perfect knowledge of the channel matrix $\mathbf{H}$ and the noise variance $N_{0}$. We also assume that all the UEs have approximately the same large-scale fading gains, and we restricted ourselves to a single precoding factor $\beta$ for all UEs. Furthermore, we have ignored real-world hardware impairments and synchronization aspects. Hence, the provided simulation results are not necessarily representative for other, more realistic system scenarios. To enable interested readers to perform their own simulations with different system parameters, we made our MATLAB simulation framework available for download from GitHub: https://github.com/quantizedmassivemimo/1bit_precoding_VLSI.

VI-B FPGA Implementation Results

To demonstrate the efficacy of C1PO and C2PO, we implemented several FPGA designs for different antenna configurations, namely for $32$, $64$, $128$, and $256$ BS antennas; all designs support downlink transmission to $16$ UEs for modulation schemes ranging from BPSK to 16-QAM. The FPGA designs were developed on register transfer level (RTL) using Verilog, implemented using Xilinx Vivado Design Suite, and optimized for a Xilinx Virtex-7 XC7VX690T FPGA. \freftbl:implresultsCxPO shows reference FPGA implementation results for C1PO and C2PO.

We see that the logic area (in terms of slices, logic LUTs, flipflops, and DSP48 units) for all designs increases roughly linearly with the number of BS antennas; this is a result of using a linear array of PEs. The only exception is the memory requirements of C1PO (in terms of memory LUTs), which scales roughly quadratically in the number of BS antennas; this is a result of having to store the entire $B\times B$ matrix $\mathbf{G}$ in contrast to storing only the augmented $(U+1)\times B$ matrix $\overline{}\mathbf{H}$ for C2PO. We also see that the logic area for C1PO is 20% to 50% smaller than that of C2PO for all array sizes; the memory area of C1PO, however, is significantly larger for $128$ and $256$ BS antennas. This is because the architecture for C1PO is slightly simpler than that of C2PO, but the memory requirements of C1PO scale quadratically in $B$ whereas the memory requirements of C2PO only scale linearly in $B$.

The maximum clock frequency for C1PO is slightly higher than that of C2PO, which is due to the slightly simpler architecture of C1PO. As expected, the maximum clock frequency slowly decreases with $B$, since the FPGA routing overhead increases with $B$. In fact, after implementing our designs, the critical paths are typically in interconnect networks. Before mapping our designs to the FPGA, however, the critical path for the C1PO designs is in the real-valued multipliers that form part of the complex multiplier, while for the C2PO designs it is in the adders that form part of the complex multiplier. The latency of one C1PO iteration is significantly larger than that of C2PO for $64$, $128$, and $256$ antennas. This results in significantly higher throughput of C2PO for these BS antenna array sizes. In summary, C2PO is more efficient in terms of throughput per area for large BS antenna array sizes (e.g, 128 BS antennas or more), whereas C1PO is more efficient for small array sizes.

We finally note that the implementation results in \freftbl:implresultsCxPO ignore the preprocessing complexity in order to compare the complexity of the precoding stage alone. For C1PO, preprocessing requires a $B\times B$ matrix inversion, which is computationally demanding, exhibits stringent data dependencies, and requires high numerical precision, especially for large BS antenna arrays [21]. In stark contrast, preprocessing for C2PO only requires the computation of the scaled MRT output, which requires a multiplication of a $B\times U$ matrix by a $U$-dimensional vector. As a result, we consider C2PO to be the preferred 1-bit precoding method for most practical BS antenna array sizes.

VI-C Comparison with MRT-based Precoding

While the papers [27, 26, 28, 5] propose hardware designs for precoding in massive MU-MIMO systems with high-precision DACs, neither of them provide detailed FPGA implementation results. Reference [27] describes an FPGA-based testbed that uses MRT and ZF-based precoding but does not report area and clock frequency results; [26] and [28] only provide ASIC implementation results, and [5] reports implementation results on a GPU cluster. Furthermore, all of these implementations were designed for high-precision DACs. Consequently, to enable a fair comparison of conventional precoders with C1PO and C2PO, we developed a baseline design that implements MRT followed by quantization (MRT-Q).

Our MRT-Q baseline design is essentially a stripped-down and heavily optimized version of C1PO with only the necessary circuitry to implement MRT-Q. More specifically, our architecture corresponds to $B/U$ linear arrays, each one with $U$ PEs and a control unit. The arrays and PEs are organized as in \freffig:C1POarch, with the exception that the projection unit is removed from the PEs. In addition, no multipliers are required as MRT-Q computes $\mathbf{H}^{H}\mathbf{s}$ with $\mathbf{s}\in\mathcal{O}^{U},$ and hence all multiplications are with constants (given by the constellation set $\mathcal{O}$) and can be implemented with adders and shifters.

The FPGA implementation results for the MRT-Q baseline designs are reported in \freftbl:implresultsMRC. Note that these designs do not require any DSP48 units as the multiplication with constants are carried out with conventional logic. In comparison to the 1-bit precoder designs reported in \freftbl:implresultsCxPO, we see that MRT-Q-based precoding is roughly $5\times$ to $6\times$ more efficient than C2PO and up to $30\times$ more efficient than C1PO (in terms of throughput/LUTs). This efficiency advantage comes at a significant loss in terms of error-rate performance (cf. \freffig:BER_ALL). We note, however, that for massive MU-MIMO systems with significantly more BS antennas than UEs (e.g., more than $8\times$), MRT-Q is a viable low-complexity alternative—a well-known fact in the massive MU-MIMO literature [2, 3, 4].

VI-D Performance–Complexity Trade-Offs

In \freffig:performancecomplexityradeoff, we provide the performance–complexity trade-offs between C1PO (dashed lines with circle markers) and C2PO (dotted lines with square markers) for various BS antenna array sizes. This trade-off is characterized in terms of the minimum normalized transmit power $\varrho$ required to achieve 1% uncoded BER for BPSK (as in \freffig:BER_32x16_BPSK); the hardware efficiency is characterized by the throughput per area (in terms of billion symbols per second per FPGA LUT). As a reference, we also show the performance for ZF precoding with infinite-precision DACs (vertical lines). As in \freffig:BER_ALL, we consider a transmission to $U=16$ UEs.

\freffig:performancecomplexityradeoff shows that, for scenarios with a high normalized transmit power $\varrho$, only a few iterations of our algorithms are required to meet 1% uncoded BER. As the value of $\varrho$ decreases, more iterations are needed, which reduces the throughput and, hence, the hardware efficiency of the circuit. We see that for small antenna arrays (i.e., for $B=32$ and $B=64$), C1PO outperforms C2PO, while for large antenna arrays (i.e., for $B=128$ and $B=256$), C2PO significantly outperforms C1PO. We note, however, that the reported hardware efficiency does not take into account the fact that the preprocessing complexity of C1PO would be substantially higher than that of C2PO; see our discussion in \frefsec:fpgaresults. We also observe that only a small number of iterations are required (e.g., $2$ to $4$ iterations) for such large BS antenna arrays to achieve the error-rate performance limits of our algorithms.

In \freffig:performancecomplexityradeoff, we additionally show the trade-off achieved by the MRT-Q baseline design reported in \frefsec:MRCbaselinedesign. Clearly, MRT-Q achieves higher throughput per LUT than C1PO and C2PO for large BS arrays ($B=128$ and $B=256$); this gain comes, however, at the cost of rather poor error-rate performance. For small BS antenna arrays ($B=32$ and $B=64$), MRT-Q is unable to achieve the target BER of $1$%. Hence, MRT-Q is only suitable for massive MU-MIMO systems with very high BS-to-UE-antenna ratios in which best-in-class error-rate performance is not the main design objective.

Remark 2**.**

The latency of C1PO and C2PO could be reduced by modifying the architectures proposed in \frefsec:architecture1 and \frefsec:architecture2. While the proposed architectures only pass one element of the vector $\mathbf{x}^{(t)}$ (for C1PO) and one element of the sub-vector $\tilde{\mathbf{x}}_{w}^{(t)}$ (for C2PO) per clock cycle, both architectures could process two or more elements per clock cycle. Such an approach would significantly decrease the latency and improve the throughput at the cost of increased silicon area.

VII Conclusions

We have proposed two nonlinear precoding algorithms, namely C1PO and C2PO, which achieve excellent error-rate performance in 1-bit massive MU-MIMO systems at low computational complexity. To substantiate this claim, we have designed corresponding VLSI architectures—to the best of our knowledge, the first for 1-bit precoding in the downlink of massive MU-MIMO systems—and we have presented FPGA reference implementations for a variety of BS antenna array configurations. Our results demonstrate that nonlinear precoding for 1-bit massive MU-MIMO systems is feasible from a hardware implementation perspective, even for antenna arrays with hundreds of BS antennas. As a result, our hardware designs pave the way for enabling BS antenna arrays with 1-bit DACs to reliably transmit high-rate data to multiple UEs, which has the potential to keep hardware complexity, system costs, and circuit power consumption within manageable limits.

There are many avenues for future work. Besides the proposed convergence results, a theoretical error-rate performance analysis of C1PO and C2PO is a challenging open research topic. Implementing precoders for other nonlinear algorithms, such as SQUID [7], which perform better than C1PO and C2PO at low normalized transmit power $\varrho$, is left for future work. The study of 1-bit nonlinear precoders using more realistic system models and a comprehensive cost, power, and performance analysis are interesting research directions. Specifically, the design of 1-bit precoding algorithms and hardware accelerators for wideband massive MU-MIMO systems that use orthogonal-frequency division multiplexing (OFDM) is the subject of ongoing work; preliminary results are reported in [42].

Appendix A Proof of \frefthm:c1po

Let $E(\mathbf{z},\mathbf{x})=\mathopen{}\left\lVert\mathbf{A}\mathbf{z}\right\rVert^{2}_{2}+\gamma\|\mathbf{z}-\mathbf{x}\|_{2}^{2}-\delta\|\mathbf{x}\|_{2}^{2}$ denote the objective (BCR∗) minimized by C1PO. Because $\mathcal{B}^{B}$ is bounded, the sequence of iterates $\{(\mathbf{z}^{(t)},\mathbf{x}^{(t)})\}$ remains bounded and thus contains a convergent sub-sequence. Denote the limit of this sub-sequence by $(\mathbf{z}^{\star},\mathbf{x}^{\star})$ and set $E^{\star}=E(\mathbf{z}^{\star},\mathbf{x}^{\star}).$ Consider the point $\hat{}\mathbf{z}^{\star}=\operatorname*{arg\;min}_{\mathbf{z}}E(\mathbf{z},\mathbf{x}^{\star})=(\mathbf{I}_{B}+\gamma^{-1}\mathbf{A}^{H}\mathbf{A})^{-1}\mathbf{x}^{\star}.$ If $\hat{}\mathbf{z}^{\star}\neq\mathbf{z}^{\star},$ then we have the strict inequality

[TABLE]

because $E$ is strongly convex in $\mathbf{z}.$ However, this contradicts the fact that $\hat{}\mathbf{z}^{\star}=\operatorname*{arg\;min}_{\mathbf{z}}E(\mathbf{z},\mathbf{x}^{\star}),$ and so it must be the case that $\hat{}\mathbf{z}^{\star}=\mathbf{z}^{\star}.$ Because $\delta<\gamma,$ $E$ is strongly convex in $\mathbf{x},$ and a similar argument shows that $\mathbf{x}^{\star}=\operatorname*{arg\;min}_{\mathbf{x}\in\mathcal{B}^{B}}E(\mathbf{z}^{\star},\mathbf{x}).$ Hence, $(\mathbf{z}^{\star},\mathbf{x}^{\star})$ minimizes $E$ with respect to $\mathbf{z}$ and $\mathbf{x}$ separately; this, combined with the fact that $E$ is differentiable, and $\mathcal{B}$ coordinate-wise separable, guarantees that $(\mathbf{z}^{\star},\mathbf{x}^{\star})$ satisfies the first-order conditions for (BCR∗); see Theorem 2 in [43] and similar arguments in [44].

Appendix B Proof of \frefthm:c2po

Let $E(\mathbf{z},\mathbf{x})=\mathopen{}\left\lVert\mathbf{A}\mathbf{x}\right\rVert^{2}_{2}-\delta\|\mathbf{x}\|_{2}^{2}$ denote the objective (9) minimized by C2PO. Let $f$ and $g$ be defined as in (10). Using the definition of the proximal operator (11) together with (13), the second update (12) of C2PO can be written as

[TABLE]

Observe that, whenever $\tau<\|\mathbf{A}^{T}\mathbf{A}\|^{-1}_{2,2},$ the inequality

[TABLE]

holds for all $\mathbf{x}.$ Using this observation, we can write

[TABLE]

This shows that the sequence $\{E(\mathbf{x}^{(t)})\}$ is monotonically decreasing. Since the sequence is bounded below, there is some limit $L=\lim_{t\to\infty}E(\mathbf{x}^{(t)}).$ Let $\{\mathbf{x}^{(t_{k})}\}$ be a convergent sub-sequence of iterates (which must exist because the iterates are bounded) with limit point $\mathbf{x}^{\star}$. Let

[TABLE]

be the result of applying the C2PO iteration starting at $\mathbf{x}^{\star}$. Observing that $E(\mathbf{x}^{(t_{k}+1)})\leq E(\mathbf{x}^{(t_{k})})\leq E(\mathbf{x}^{(t_{k}-1)}),$ and letting $k\to\infty,$ we find that $E(\bar{}\mathbf{x}^{\star})=E(\mathbf{x}^{\star})=L,$ and so $\mathbf{x}^{\star}$ is a minimizer of (15). This is only possible if $0\in\partial g(\mathbf{x}^{\star})+\nabla f(\mathbf{x}^{\star}),$ in which case $\mathbf{x}^{\star}$ is a stationary point.

Acknowledgments

The authors would like to thank O. Tirkkonen for insightful discussions on 1-bit precoding. The authors also thank A. Burg for discussions on the hardware architecture and R. Manohar for pointing us to its connection to Cannon’s algorithm. The work of O. Castañeda and C. Studer was supported in part by Xilinx, Inc. and by the US National Science Foundation (NSF) under grants ECCS-1408006, CCF-1535897, CAREER CCF-1652065, and CNS-1717559. The work of S. Jacobsson and G. Durisi was supported by the Swedish Foundation for Strategic Research under grant ID14-0022, and by the Swedish Governmental Agency for Innovation Systems (VINNOVA) within the center ChaseOn. The work of T. Goldstein was supported in part by the US NSF under grant CCF-1535902 and by the US Office of Naval Research under grant N00014-17-1-2078.