Vector Precoding for Gaussian MIMO Broadcast Channels: Impact of Replica Symmetry Breaking

Introduction

The multiple-input multiple-output (MIMO) Gaussian broadcast channel (GBC) is the focus of many research activities, addressing the growing demand for higher throughput wireless systems, and in particular the increasing use of multiple-antenna systems in essentially all modern wireless standards (see, e.g., ). The capacity region of the MIMO GBC is the dirty paper coding (DPC) capacity region , and several attempts have been made in recent years to propose practically oriented approaches for implementing DPC, as e.g., . DPC still remains, however, a difficult, computationally demanding task, which motivates the search for more practical (suboptimum) precoding alternatives.

Since linear precoding, such as zero-forcing (ZF), leads to reduced performance (especially when the channel is ill-conditioned), much attention has been given to nonlinear precoding schemes. In particular, lattice-based precoding approaches have often been investigated, as for example the vector perturbation approach suggested in (see also for a general framework). The vector perturbation approach was inspired by the idea of Tomlinson-Harashima precoding (THP) . In this scheme, a scaled complex integer vector is added to each data vector, chosen to minimize the energy penalty imposed by a linear zero-forcing (ZF) front-end. A modulo function is employed at the receivers, uniquely determining the transmitted symbols in the absence of noise. An analogous precoding scheme based on a linear minimum-mean-squared-error (MMSE) front-end was considered in . An approach based on optimizing mutual information was taken in . Vector perturbation is however still complex as it involves the solution of an NP-hard integer-lattice least squares problem (commonly implemented using the sphere-decoding algorithm ). Addressing the complexity aspect of the method, related approaches can also be found, e.g., in (see references therein for additional literature in this framework), where lattice-basis reduction techniques are employed.

The analytical performance analysis of such nonlinear precoding schemes is not at all trivial. It is common to consider, therefore, uncoded symbol error probabilities (via simulations), asymptotic capacity scaling laws and diversity orders (the asymptotic slope of the error probability in the high signal-to-noise ratio (SNR) regime), or to employ Monte-Carlo simulations to obtain information-theoretically achievable rates (see e.g., , and also for a semi-tutorial review in this respect). The energy penalty induced by the linear front-end is another commonly addressed performance measure. A lower bound on the energy penalty based on lattice theoretic arguments can be found in . The optimum constellation shaping for a ZF front-end (in terms of the energy penalty), allowing for data to be independently decoded by the users, is investigated in , where a selective mapping technique is introduced based on random coding arguments, implementable using nested lattice coding in a trellis precoding framework (see also for a more recent study on selective mapping).

The energy penalty minimization was also investigated in where another nonlinear precoding approach in this framework was recently proposed. The transmitter comprises a linear front-end combined with nonlinear precoding. The nonlinear part relies on relaxation of the transmitted symbols’ alphabets to larger alphabet sets. The idea is to optimize the vector of transmitted symbols over the extended alphabet sets, so as to minimize the energy penalty imposed by the linear front-end, which is essentially the idea behind vector perturbation. However, a notable feature of this precoding scheme is that it can also be combined with convex extended alphabet sets (in contrast to ), lending themselves to efficient practical energy minimization algorithms. It can be considered in this sense as a generalization of the vector perturbation scheme (see also in this respect).

Another interesting contribution of is the harnessing of statistical physics tools for the analysis of the nonlinear precoding scheme, while considering the large system limit in which both the number of users and the number of transmit antennas grow large, while their ratio goes to some finite constant. One of the main objectives of statistical physics is the quantitative description of macroscopic properties of many-body systems while starting from the fundamental interactions between microscopic elements. In this framework, a general tool for the analysis of random (“disordered”) systems, referred to as the “replica method”, was originally invented for the analysis of spin glasses. The latter term describes a spin orientation that has similarity to the type of location of atoms in glasses, which are random in space but frozen in time . However, the replica method turns out to have a much wider range of applications (see, e.g., for recent tutorial manuscripts). In recent years, in particular following Tanaka’s pioneering work , the method has been successfully applied to various problems in wireless communications. The replica method has also been recognized by now as an important tool for information-theoretic analyses in cases where “conventional” random matrix theory does not apply. Although the replica method is heuristic in nature, extensive simulations and exact analytical results in the literature suggest that the replica analysis generally yields excellent approximations in many cases of interest (see again , and also, e.g., and references therein).

The replica analysis usually employs a number of underlying assumptions regarding the behavior of the quantities in concern in the large-system limit. One such fundamental assumption is the “self-averaging” property, which relies on the expectation that macroscopic properties of large random systems converge to deterministic values as the system dimensions grow large. Self-averaging is a property of most physical systems with large (or infinite) degrees of freedom, and is the result of the high probability of occurrence of typical events or samples. Nevertheless, in the case of glassy systems this property has been not trivial to prove, since the underlying randomness of the interactions makes the system inherently non-ergodic. The self-averaging property for the conventional so-called Sherrington-Kirkpatrick (SK) spin glass model was first proven in . More recently, generalized it to a more general class of spin glass models and, using an ingenious method, showed that the averages over the disorder actually do converge in the large system limit. Even more recently, code-division multiple-access (CDMA) systems were shown to be self-averaging . Although the particular system we study is not explicitly covered by the above analysis, it can readily be proved to be self-averaging using the same method. For the sake of space, we will not cover the proof here, however, we will refer to self-averaging as a fundamental property of the large system limit rather than an assumption. Another common assumption in replica analyses is that of replica symmetry (RS) (see, e.g., ), according to which it is assumed that the crosscorrelations between replicated microscopic system configurations are independent of the replica indices. The RS assumption, however, is known to produce incorrect conclusions for certain physical quantities such as, e.g., the minimum energy configuration. This led to the development of the replica symmetry breaking (RSB) theory . Recently, the full RSB solution of the SK spin glass model, first proposed in , was shown to be an upper bound to the minimum energy configuration and later to be the exact solution of the model . Apart from its general seminal importance, it is profoundly relevant in the context of vector precoding because the SK-model is a particular case of the more general models discussed in and in the sequel.

In this paper we consider a communication system setting in which RSB indeed occurs, and demonstrate the significant impact of the RSB treatment on the validity of the approximations produced by the replica analysis. We focus here on a wireless MIMO broadcast channel (BC) setting, where the transmitter has $N$ transmit antennas and the $K$ users have single receive antennas. Full channel state information (CSI) is assumed available at the transmitter, while the receivers are cognizant of their own channels only (more on this later). No user-cooperation of any kind is assumed. The received signals are embedded in additive white Gaussian noise (AWGN). The precoding approach considered in is revisited. Note that the focus in is mainly on presenting the method, and on the derivation of the energy penalty in the asymptotic regime, in which both the number of transmit antennas $N$ and the number of users $K$ go to infinity, while $K/N\rightarrow\alpha<\infty$ (commonly referred to as the system load). Furthermore, the analysis in is based on the RS assumption. It turns out however that the RS assumption can only produce valid asymptotic approximations in this setting when the extended alphabets are convex sets (see, e.g., supportive simulation results in ). In contrast, for the non-convex alphabets considered in these approximations turn out to be rather loose, and produce overoptimistic results, especially as the system load gets close to unity. This behavior can be readily observed by comparing the RS based energy penalty to the asymptotic lower bound of .

Here, an alternative analysis is provided based on what is referred to in the statistical physics literature as the one-step RSB (1RSB) ansatz, which allows one to search for more general solutions than the RS ansatz, but does not cover the full complexity of solutions of full RSB. In addition to an energy penalty analysis, analogous to the one in , we complement the results by providing an information-theoretic perspective of the proposed precoding approach. Coded transmissions and achievable throughputs are considered. The employed performance measure is the normalized spectral efficiency, defined as the total number of bits/sec/Hz per transmit antenna that can be transmitted arbitrarily reliably through the broadcast channel. The limiting marginal conditional distribution of the nonlinear precoder’s output which is required for the calculation of spectral efficiency, as well as the limiting energy penalty, are analytically formulated. Focusing on a ZF front-end, the spectral efficiency is expressed via the input-output mutual information of the equivalent single-user channel observed by each of the receivers. The analysis is applied next to a particular family of discrete extended alphabet sets (following ), focusing on a QPSK input, demonstrating the RSB phenomenon. To complete the analysis we repeat the derivations while employing the RS ansatz, which is, as said, adequate for convex relaxation schemes, and the results are then applied to a convex alphabet example . For both extended alphabets, numerical spectral efficiency results indicate significant performance enhancement over linear ZF preprocessing for medium to high SNRs. Furthermore, performance enhancement is also revealed compared to a generalized THP approach (which is a popular practical nonlinear precoding alternative for such settings). Comparison of the two types of extended alphabet examples leads to interesting conclusions regarding the performance vs. complexity tradeoff of precoding schemes of the kind considered here.

The remainder of this paper is organized as follows. Section 2 describes the system model. Section 3 provides an outline of the replica analysis and includes some general results. In particular, it clarifies the concept of RSB which later results are based upon. In order to analyze the mutual information and later the trade-off between spectral and power efficiency of various precoding schemes, we need to characterize the limiting conditional distribution of the precoder output. This task is solved in Section 4 providing a set of nonlinear equations whose solutions characterize the desired distributions. Section 5 particularizes to the ZF front-end and shows that the channel model can be represented as an equivalent concatenated single-user channel. Then, it derives the spectral efficiency of this equivalent concatenated channel. Section 6 particularizes the results of the previous sections to a discrete lattice-based alphabet relaxation of QPSK. Numerical solutions of the analytical results are provided. Those based on RSB are shown to match simulation results while those based on RS are demonstrated to fail. Section 7 is the corresponding counterpart to Section 6 for convex relaxation. Unlike Section 6, it finds the RS ansatz to provide accurate approximations. Section 8 presents a comparative analysis of the spectral efficiency of the two alphabet relaxation schemes against some other precoding approaches. Finally, Section 9 ends this paper with some concluding remarks.

System Model

Consider the following Gaussian MIMO broadcast channel

where $\boldsymbol{r}_{[K\times 1]}$ is the vector of received signals, $\boldsymbol{H}_{[K\times N]}$ is the (random) complex channel transfer matrix, assumed to be of unit expected row norm, $\boldsymbol{t}_{[N\times 1]}$ is the vector of transmitted signals, and $\boldsymbol{n}_{[K\times 1]}$ is the vector of i.i.d. zero mean proper complex AWGNs at the users’ receivers. We denote the noises’ spectral levels by $\sigma^{2}$ so that $\boldsymbol{n}\sim\mathcal{N}_{c}(\boldsymbol{0},\sigma^{2}\boldsymbol{I})$.

Let $\boldsymbol{u}_{[K\times 1]}$ denote the vector of the encoders’ outputs, i.e., $\boldsymbol{u}=[u_{1},\dots,u_{K}]^{T}\in\mathscr{U}^{K}$. The vector $\boldsymbol{u}$ is the input to a nonlinear precoding block that minimizes the energy penalty of the precoder through input alphabet relaxation (see below), and outputs a $K\times 1$ vector $\boldsymbol{x}$. The vector $\boldsymbol{x}$ is then taken as input to the linear front-end block where it is multiplied by the linear front-end matrix $\boldsymbol{T}_{[N\times K]}$, which is, in general, a function of the channel transfer matrix $\boldsymbol{H}$ (note that $\boldsymbol{x}$ depends on $\boldsymbol{T}$ and, hence, we can use the functional notation $\boldsymbol{x}(\boldsymbol{T},\boldsymbol{u})$). The result is then normalized so that the actually transmitted vector $\boldsymbol{t}$ satisfies an instantaneous total power (energy per symbol) constraint $P_{\textsf{tot}}$, i.e.,

where $\mathscr{E}^{\textsf{tot}}(\boldsymbol{T},\boldsymbol{x})$ denotes the energy penalty induced by the precoding matrix $\boldsymbol{T}$, and the particular choice of $\boldsymbol{x}$, as well as the average symbol energy of the underlying alphabet $\mathscr{U}$ (the explicit dependence on the arguments is omitted henceforth for simplicity). Denoting by $P$ the individual power constraint per user (taken as equal for all), so that $P_{\textsf{tot}}=KP$, we define the transmit SNR as

We note at this point that as an alternative to the normalization taken in (2-3), ensuring an instantaneous transmit power constraint, a weaker average transmit power constraint can be applied, by simply replacing $\mathscr{E}^{\textsf{tot}}$ with $E\left\{\mathscr{E}^{\textsf{tot}}\right\}$, where $E\left\{\cdot\right\}$ denotes expectation. However, since we later concentrate on the energy penalty per symbol,

and in view of the self-averaging property of the large system limit (as shall be made clear in the following), the two types of energy constraints yield the same asymptotic results. We thus focus for convenience throughout this paper on the instantaneous power constraint (as implied by (2-3)). Note also that in order to differentiate between the energy penalty induced by the precoding scheme, and the effect of the underlying symbol energy of the input alphabet $\mathscr{U}$, one can alternatively represent the results in terms of what we refer to here as the precoding efficiency, defined through

where $\sigma_{u}^{2}=E\{\left|u\right|^{2}\}$ (with the expectation taken with respect to (2-2)).

Outline of the Replica Analysis

In the following we describe the main ideas behind the replica analysis of the problem in hand, and provide a heuristic outline of the approach taken to derive the main results of this paper. The reader is referred to tutorial manuscripts such as for an elaborated background on the replica analysis. The fully detailed proofs are deferred to the appendices.

We start here by focusing on the energy penalty, and note that the task of the nonlinear precoding block at the transmitter (see Figure 1) can be described as follows. Its task is equivalent to the minimization of an objective function (called the Hamiltonian in physics literature) having the quadratic form

with $(\cdot)^{\dagger}$ denoting transpose conjugation and $\boldsymbol{J}$ being a random matrix of dimensions $K\times K$. Thus, the minimum energy penalty per symbol can be expressed as

where we use the shortened notation $\mathscr{B}_{\boldsymbol{u}}\triangleq\mathscr{B}_{u_{1}}\times\cdots\times\mathscr{B}_{u_{K}}$. Note also that to comply with (2-5) one should take $\boldsymbol{J}=\boldsymbol{T}^{\dagger}\boldsymbol{T}$, however since the results derived in the sequel hold, at least in part, for a more general class of matrices, we retain the formulation as in (3-1).

To calculate the minimum of the objective function as defined in (3-1), it is convenient to introduce some notions from statistical physics (see, e.g., ). In particular, we define a discrete probability distribution on the set of state vectors $\left\{\boldsymbol{x}\right\}$, namely the Boltzmann distribution, as

where the parameter $\beta>0$ is referred to as the inverse temperature $\beta=1/T$, while the normalization factor $\mathcal{Z}$ is the so-called partition function, which is defined as

The definitions above hold for both discrete and continuous alphabets $\mathscr{B}_{\boldsymbol{u}}$. The only difference is that for continuous alphabets the sums over $\boldsymbol{x}\in\mathscr{B}_{\boldsymbol{u}}$ are replaced by integrals.

At thermal equilibrium, the energy of the system is preserved, while the second law of thermodynamics states that the entropy is the maximum possible. This is equivalent to minimizing the free energy of the system

where $\beta$, the inverse temperature, is in fact the Lagrange multiplier in the maximization of (3-6), subject to the mean energy constraint. At equilibrium, the free energy can be expressed as

Note that from Lagrangian duality the Boltzmann distribution (3-3) is also the solution to the problem of minimizing the energy for a given entropy.

All mean thermodynamic quantities can now be derived directly from the free energy. In particular, the energy of the system is

while its thermodynamic entropy (disorder) is

In addition to the above quantities we can use the free energy and the partition function to obtain the empirical joint distribution of the precoder input $u$ and output $x$, which is defined for general $\beta$ as

Eqs. (3-3) to (3-11) will be useful in deriving some of the results presented in the sequel.

The rationale behind the introduction of the Boltzmann distribution is that as $\beta\rightarrow\infty$, the partition function becomes dominated by the terms corresponding to the minimum energy. Hence, taking the logarithm and further normalizing with respect to $\beta$, one gets the desired limiting quantity (energy, entropy or empirical distribution) at the minimum energy subspace of $\mathscr{B}_{\boldsymbol{u}}$. Note that even if the energy minimizing vector is not unique, or in fact even if the number of such vectors is exponential in $K$, one still gets the desired quantity when taking the limit $\beta\rightarrow\infty$.

It is crucial to point out that in the above summation over the set of state-vectors $\mathscr{B}_{\boldsymbol{u}}$, both the input vector $\boldsymbol{u}$ and the matrix $\boldsymbol{J}$ are fixed. These random variables are called quenched. Therefore, all the above manipulations still do not alleviate the difficulty of calculating the desired quantities. In particular, the main difficulty comes from the free energy being a random variable itself, which depends on the particular realizations of $\boldsymbol{J}$ and $\boldsymbol{u}$. As discussed in Section 1, the proofs of the self-averaging property of the SK-model in can be generalized to apply to the form of $\mathcal{H}(\boldsymbol{x})$ analyzed here. This means that the free energy converges in probability at the asymptotic limit to a non-random quantity, i.e.,

where the expectation $E\{\cdot\}$ is over all realizations of $\boldsymbol{J}$ and $\boldsymbol{u}$. As a result, all quantities that can be obtained from the free energy in an analytic manner, e.g., by differentiation of a parameter, are also self-averaging. The empirical joint distribution of the precoder input and output converges to a non-random distribution which is expressed by (3-16). This self-averaging property makes the problem more straightforward to tackle, since we may now hope to get analytic results for the average of the free energy and its derivatives.

With that in mind, the limiting energy penalty (per symbol) can be represented as

If we add this term to $\mathcal{H}(\boldsymbol{x})$ in the exponent, the partition function gets modified to

where we have dropped the explicit dependence of $\mathcal{Z}(h)$ and $V(h)$ on $\upsilon$ and $\xi$ (as well as $\boldsymbol{u}$ and $\boldsymbol{x}$) for the sake of notational compactness. In the sequel, any dependence on $h$ shall implicitly also indicate a dependence on $\upsilon$ and $\xi$. Using the above partition function we obtain a modified free energy using (3-8). Upon differentiation with respect to $h$, setting $h=0$, and letting $\beta\rightarrow\infty$ we get

where $\mathcal{F}(\beta,h)$ denotes the free energy for the modified partition function $\mathcal{Z}(h)$. An alternative method for deriving the limiting empirical distribution, which relies on the limiting moments, can be found in , albeit with more restrictive assumptions on the limiting distribution.

The next step in the analysis is to invoke some underlying assumptions. The first assumption is that the random matrix $\boldsymbol{J}$ can be decomposed as

where $\boldsymbol{D}$ is a diagonal matrix with diagonal elements being the eigenvalues of $\boldsymbol{J}$, and $\boldsymbol{U}$ is a unitary Haar distributed matrix . It is further assumed that the empirical distribution of the diagonal elements of $\boldsymbol{D}$ converges to a nonrandom distribution uniquely characterized by its $R$-transform For a definition of the $R$-transform, see Appendix E. $R(\cdot)$, which is assumed to exist.

Going back to the original communication system model, note that we are in fact interested in the normalized averages of most of the quantities described above, at the limit as $K\rightarrow\infty$. Therefore, to make a distinction, while retaining the relation between the quantities, we shall use henceforth the following notational convention

Calculating the expectation of a logarithm of a sum of exponents (see (3-13)) is a formidable task. The standard approach in statistical physics is to invoke the so-called replica ‘‘trick’’. The latter is based on the following identity An equivalent representation often encountered in the literature is $E\left\{\log\mathcal{Z}\right\}=\lim_{n\rightarrow 0^{+}}\frac{E\left\{\mathcal{Z}^{n}\right\}-1}{n}$.

which holds in general for real $n$. The “trick” here relies on the assumption that the right hand side (RHS) of (3-22) can be evaluated for integer $n$, and that the desired quantity can be found by analytic continuation in the vicinity of $n=0^{+}$. Although this “trick” does not a priori have any justified validity, its success in statistical physics, and more recently in communications theory, makes it a reasonable approach. Further assuming that the limits with respect to $K$ and $n$ can be interchanged (which is the common practice in replica analyses), (3-13) can be rewritten as

where we use the notation $\sum_{\left\{\boldsymbol{x}_{a}\right\}}=\sum_{\boldsymbol{x}_{1}\in\mathscr{B}_{\boldsymbol{u}}}\cdots\sum_{\boldsymbol{x}_{n}\in\mathscr{B}_{\boldsymbol{u}}}$, and $\Tr(\cdot)$ denotes the trace operator.

The summation over the replicated precoder output vectors $\left\{\boldsymbol{x}_{a}\right\}_{a=1}^{n}$ in (3-23) is performed by splitting the replicas into subshells, defined through an $n\times n$ matrix $\boldsymbol{Q}$

The limit $K\to\infty$ allows us to perform the following derivations by saddle point integration. This first yields the following general result.

For any inverse temperature $\beta$, any structure of $\boldsymbol{Q}$ consistent with (3-24), and any $R$-transform $R(\cdot)$ such that $R(\boldsymbol{Q})$ is well-defined Note that if $R(\cdot)$ has a series expansion, $R(\boldsymbol{Q})$ is well-defined. Since $R(\cdot)$ is the free cumulant generating function, $R(\boldsymbol{Q})$ is well-defined, if all moments of the asymptotic eigenvalue distribution of $\boldsymbol{J}$ exist. , the energy is given by

where $\boldsymbol{Q}$ is the solution to the saddle point equation

with $\mathscr{B}_{u}^{n}$ denoting the $n$-fold Cartesian product of $\mathscr{B}_{u}$.

With the help of Proposition 3.1, the energy can be written as

with $\boldsymbol{Q}$ given by (3-26). In (3-27), $\boldsymbol{x}$ is a $K$-dimensional vector and its components represent users. The contributions of the users to the energy arise due to the inner product $\boldsymbol{x}^{\dagger}\boldsymbol{Jx}$ and are coupled, unless $\boldsymbol{J}$ is diagonal. In (3-28), $\bf x$ is an $n$-dimensional vector and its components represent replicas of the same user. The contributions of the users to the energy arise due to integration over the distribution $F_{U}(u)$, and are decoupled and additive. This is just another incarnation of the decoupling principle that, under the assumption of replica symmetry, was addressed in . Here, we find that it holds for the energy of general (also replica symmetry breaking) spin glass systems and their equivalents in communication theory.

Another interesting observation is the following. In , an analogy between the $R$-transform and effective interference in linear MMSE detection was discovered, and the additivity of the effective interference of coupled users was explained based on the additivity of the $R$-transforms of free random variables. Relying on the code symbols of different users $\left\{u_{k}\right\}$ being i.i.d., we can rewrite (3-28) as

and interpret $\mathscr{E}_{k}(\beta)$ as the effective energy of user $k$. Like the effective interference in , it depends only on the signal constellation of user $k$ and the $R$-transform, and it is additive among users. In contrast to , (3-29) is more general and neither constrained to linear detectors nor to Gaussian symbol alphabets.

for some constants $\left\{q_{0},\chi_{0}\right\}$. The 1RSB assumption leads to a more involved structure, formulated as

using the constants $\left\{q_{1},p_{1},\chi_{1},\mu_{1}\right\}$. The above constants (i.e., $\left\{q_{0},\chi_{0}\right\}$ for RS, and $\left\{q_{1},p_{1},\chi_{1},\mu_{1}\right\}$ for 1RSB) are referred to as macroscopic parameters, and obtained from the corresponding saddle point equations. The limiting energy penalty can then be expressed in terms of these macroscopic parameters, as shown in the following sections. An analogous procedure can be employed to obtain the limiting empirical joint distribution of the precoder input and output using (3-16).

Replica symmetry breaking is not limited to one step, and in fact in order to exactly characterize the limiting energy penalty and precoder output statistics, we would eventually need to consider full RSB, as discussed in Section 1. However, we will only present here precoding results up to the accuracy of 1RSB for purposes of analytical tractability. For the interested reader and the sake of completeness, we include general results on multiple-step RSB in Appendix A.

Limiting Characterization of the Precoder Output

We restrict ourselves in the following to 1RSB analysis of the limiting characteristics of the precoder output. As demonstrated in the sequel, when compared to simulation results at finite numbers of antennas, 1RSB gives quite accurate approximations for the quantities of interest, while the RS ansatz does so only in special cases.

Applying the 1RSB ansatz, as outlined in Section 3 (see in particular (3-31)), the limiting properties of the precoder output are characterized by means of four macroscopic parameters $q_{1},p_{1},\chi_{1},\mu_{1}\in(0,\infty)$, which are determined as specified below. Let $\boldsymbol{J}$ be a $K\times K$ random matrix satisfying the decomposability property (3-18), and let $R(\cdot)$ denote the $R$-transform of its limiting eigenvalue distribution. Consider now the following function of complex arguments

where $\Re\left\{\cdot\right\}$ takes the real part of the argument, and the parameters $\varepsilon_{1}$, $g_{1}$ and $f_{1}$ are defined as

Furthermore, denote its normalized version by

the parameters $\left\{q_{1},p_{1},\chi_{1},\mu_{1}\right\}$ are given by the solutions to the four coupled equations In general these coupled equations have multiple solutions and one needs to choose the solution that minimizes the energy penalty.

The limiting properties of the precoder outputs can now be summarized by means of the following two propositions. The detailed proofs are provided in Appendices B.1 and B.2, respectively.

Suppose the random matrix $\boldsymbol{J}$ satisfies the decomposability property (3-18). Then under some technical assumptions, including in particular one-step replica symmetry breaking, the effective energy penalty per symbol $\mathscr{E}^{\textsf{tot}}/K$ converges in probability as $K,N\rightarrow\infty$, $K/N\rightarrow\alpha<\infty$, to

The conditional limiting empirical distribution of the precoder’s outputs is specified next.

With the same underlying assumptions as in Proposition 4.1, the limiting conditional empirical distribution of the nonlinear precoder’s outputs given an input symbol $u$ satisfies

where $1\!\left\{\cdot\right\}$ denotes the indicator function.

2 A Replica Symmetric Reduction

Although the 1RSB solution of the replica analysis leads in principle to a more accurate description of the large system limit, corresponding results can also be derived using the simplifying assumption that the system exhibits a replica symmetric behavior (see (3-30)). These results shall be used in the sequel to demonstrate the impact of replica symmetry breaking. However they can also be extremely useful for more conveniently analyzing settings that do exhibit replica symmetric properties, such as the case of convex extended alphabet sets addressed in . A convex alphabet example is discussed in Section 7.

The limiting energy penalty under the RS assumption was in fact already derived in , and the result is recalled in the following proposition. The result is given in terms of the two macroscopic parameters $q_{0},\chi_{0}\in(0,\infty)$, which are obtained through the solution of the two coupled equations

Suppose the random matrix $\boldsymbol{J}$ satisfies the decomposability property (3-18). Then under some technical assumptions, including in particular replica symmetry, the effective energy penalty per symbol $\mathscr{E}^{\textsf{tot}}/K$ converges in probability as $K,N\rightarrow\infty$, $K/N\rightarrow\alpha<\infty$, to In , the self-averaging property was stated as an assumption, since the authors were not aware of .

The limiting conditional distribution of the precoder outputs can also be characterized under the RS assumption, in an analogous manner to Proposition 4.2.

With the same underlying assumptions as in Proposition 4.3, the limiting conditional empirical distribution of the nonlinear precoder’s outputs given an input symbol $u$ satisfies

This is the measure of the corresponding Voronoi region in the scaled conditional signal constellation $\mathscr{B}_{\upsilon}$, with respect to the (complex) Gaussian probability measure.

Proof: The proof follows the same steps as in the proof of Proposition 4.2, while replacing (3-31) with (3-30).

3 Zero-Temperature Entropy

One way to demonstrate the degree of consistency of the RS and 1RSB solutions is to look at their limiting (thermodynamic) zero-temperature entropy defined as $\bar{\mathscr{S}}=\lim_{\beta\to\infty}{\mathscr{S}}(\beta)$. It can also be obtained in a manner similar to Propositions 4.1 and 4.3. In Appendix B.3, we show:

With the same underlying assumptions as in Propositions 4.1 or 4.3, the limiting entropy per symbol converges to

with $\chi$ denoting $\chi_{1}$ and $\chi_{0}$ for 1RSB and RS, respectively.

In any stable thermodynamic system the entropy is non-negative for all temperatures. However, one of the main pitfalls of the RS solution of the original SK-model is that its zero-temperature entropy is negative, indicating an instability . For all $R$-transforms that are strictly increasing functions of negative real arguments, Proposition 4.5 clearly implies that the entropy is always negative, becoming zero only when the zero temperature value of $\chi_{1}$, respectively $\chi_{0}$, approaches zero. While the full RSB solution has been shown to have vanishing entropy at zero temperature and corresponds to the correct solution, the following lemma proven in Appendix E, indicates that negative entropy is a rather common effect for finite RSB steps.

The R-transform, wherever its derivative with respect to a real argument exists, is an increasing function. If the probability distribution is different from a single mass point, the R-transform is strictly increasing.

Note that the above argument for the entropy holds only for discrete state variables. In the case of continuous alphabets, the (then differential) entropy of a system can in fact be negative. Therefore, a negative zero-temperature entropy is not an alarm bell per se. For discrete state variables, the zero-temperature entropy serves as a measure of accuracy: the closer it is to zero, the better the approximation.

Zero-Forcing Front-End

To gain more insight into the impact on system performance of the nonlinear precoding scheme under investigation, we now particularize to a specific linear front-end, namely the ZF front-end. The precoding matrix $\boldsymbol{T}$ in this case is given by the pseudo-inverse of the channel transfer matrix, which we write here as

The underlying assumptions are that $N\geq K$ and that the matrix $\boldsymbol{H}\boldsymbol{H}^{\dagger}$ is almost surely (a.s.) positive definite In Section 8, we will also allow for $N<K$ following the treatment in .. Focusing on the asymptotic regime for $K/N\rightarrow\alpha\leq 1$, then using (2-1), (2-3), and Proposition 4.1, the equivalent single-user channel observed by user $i$ is

where $\check{n}_{i}$ is a zero mean circularly symmetric complex Gaussian noise with variance $\frac{1}{\rho}$,

denotes the effective received SNR, and $\bar{\mathscr{E}}_{\textsf{rsb1}}$ is given by (4-12).

with $P_{X|U}(\mathsf{x}|u)$ given by (4-13), and

is the (complex) Gaussian density with mean $x$ and variance $1/\rho$.

Proof: The Proposition follows straightforwardly from Proposition 4.2 and (5-2).

The achievable throughput of the nonlinear precoding scheme can be derived from the equivalent single-user channel model using Proposition 5.1. Accordingly, the achievable rate of a randomly chosen user is given by the mutual information Note that in the large-system limit the receivers only need information about the state of their own channel, but not about the states of the other channels due to the self-averaging property which makes the impact of the other users’ channels and data deterministic. between the input $u$ and received signal $y$, i.e.,

where ${\rm h}(\cdot)$ and ${\rm h}(\cdot|\cdot)$ denote differential entropy and conditional differential entropy, respectively (which can be readily calculated using Proposition 5.1). The normalized spectral efficiency is then given by

and it is functionally dependent on the system average $E_{b}/N_{0}$ through the relation

To get a better insight into the impact of the nonlinear precoding scheme, it is useful to compare the results to the spectral efficiency of DPC with Gaussian input (specifying the ultimate performance), as well as to the spectral efficiency of linear ZF (for both Gaussian and discrete alphabet input). Another interesting comparison is to the spectral efficiency of generalized THP (GTHP), which is a popular practical nonlinear precoding alternative to the scheme considered here (see, e.g., ). For the sake of comparison we further particularize henceforth to the case in which the entries of the channel transfer matrix $\boldsymbol{H}$ are i.i.d. zero-mean circularly symmetric complex Gaussian random variables, with variance $1/N$ (“a Gaussian $\boldsymbol{H}$”). Note that in this case the $R$-transform of the limiting eigenvalue distribution of the random matrix $\boldsymbol{J}=(\boldsymbol{H}\boldsymbol{H}^{\dagger})^{-1}$, and its derivative, simplifiy to

Starting with DPC, the limiting spectral efficiency in this setting coincides with the corresponding spectral efficiency of the dual uplink channel with uniform power distribution . This follows from the limiting conclusion in , and by observing that the optimization problem over diagonal input covariance matrices, that specifies the maximum achievable sum-rate (see and references therein), is solved by a uniform power distribution . The spectral efficiency of DPC is hence given by

where ${\mathcal{F}}({\textsf{snr}},\alpha)$ is defined as

Regarding linear ZF, we restrict the discussion to the case in which the active user population can only be controlled through the system load $\alpha$, as is in fact assumed for the nonlinear precoding scheme (see also the discussion in Section 8). In this setting, as shown, e.g., in , the induced precoding efficiency (2-7) (equivalent here to the inverse multiuser efficiency) converges in the large system limit to

and again for Gaussian input the spectral efficiency coincides with the corresponding result in (see also )

The corresponding spectral efficiency with discrete input alphabet can be derived, e.g., following the guidelines in . Considering the particular case of binary phase shift keying (BPSK) input, one obtains

The spectral efficiency of linear ZF precoding combined with QPSK input is obtained via the relation

yielding (through (5-9)) $C^{{\textsf{zf}},{\textsf{qpsk}}}(\frac{E_{b}}{N_{0}})=2C^{{\textsf{zf}},{\textsf{bpsk}}}\left(\frac{E_{b}}{N_{0}}\right)$. The spectral efficiency of GTHP for the corresponding setting is derived in Appendix F.

Lattice Precoding: An RSB Example

Adhering to , we consider in the following a particular example of a discrete relaxed alphabet set for QPSK signaling, which exhibits replica symmetry breaking. The original QPSK constellation alphabet is represented by the set

and quadrature symmetric transmissions are assumed (note that the above definition induces $\sigma_{u}^{2}=2$). The relaxed alphabets in this particular example can be represented as points from the extended lattice

where it is assumed that $-\infty=c_{0}<c_{1}<\cdots<c_{L}<c_{L+1}=\infty$. The parameter $L$ thus specifies the number of lattice points used in the extended alphabet in each dimension, and we particularize here to the set $\left\{+1,-3,+5,-7,+9,\dots\right\}$. The alphabet relaxation scheme is depicted in Figure 3. Due to the complete quadrature symmetry of this setting, all QPSK constellation points and their corresponding relaxed alphabet subsets are completely equivalent, and we focus in the following, for notational convenience, on the QPSK constellation point represented by $u=1+j$, and $\mathscr{B}_{1+j}$.

where we introduced the real argument function

Applying this observation to (4-8)–(4-11), and exploiting the quadrature symmetry property, the derivation simplifies considerably by noticing that the inner integrals therein can be represented as sums of separate integrals over the regions specified by (6-5). Accordingly, consider the two real argument functions

Then, following some tedious algebra, it can be shown from (4-8)–(4-11) that the parameters $\left\{q_{1},p_{1},\chi_{1},\mu_{1}\right\}$ are the solutions to the coupled equations

The corresponding energy penalty is obtained by plugging the four solutions into (4-12). Applying the same approach to Proposition 4.2, the limiting conditional probability of the precoder output being $x=c_{m}+jc_{n}\in\mathscr{B}_{1+j}$ is given by

The limiting conditional probabilities that correspond to the rest of the QPSK constellation points are readily obtained from (6-13) by symmetry considerations. Note also that (6-13) implies that the real part and the imaginary part of the precoder output $x$ behave as independent random variables.

Numerical results for the limiting energy penalty of the discrete lattice relaxation scheme are plotted in Figure 4.

The figure shows the limiting energy penalty (in dB) as a function of the system load $\alpha$, for the particular case of a Gaussian $\boldsymbol{H}$ and a ZF front-end. Since $\sigma_{u}^{2}=2$, the corresponding precoding efficiency (2-7) can be immediately obtained by subtracting $3\rm dB$ from the energy penalty shown in the figure. The results in Figure 4 correspond to alphabet relaxations with $L=2$ and $L=3$. Note that the two curves are essentially indistinguishable and the energy penalty with $L=2$ becomes only negligibly larger as $\alpha$ gets close to unity. This implies that increasing $L$ beyond $2$ in this setting provides diminishing returns. Empirical energy penalties obtained through Monte Carlo simulations are also included in the figure. The results are for systems in which the number of users is fixed to $K=8$, $K=16$, and $K=32$ (averaged over $10^{4}$, $10^{3}$, and $10^{2}$ channel realizations, respectively). The energy penalty is shown to decrease with the system size, and the simulation results exhibit a good match to the limiting energy penalty predicted by the 1RSB replica analysis. The lower bound for the limiting energy penalty obtained in is also plotted in this figure which, with appropriate scaling to match the current setting, is given by

Figure 4 shows that the 1RSB prediction approaches the lower bound as the load approaches unity. Note however that the 1RSB result stays strictly higher than the lower bound. In fact, a careful numerical examination of the limiting 1RSB energy penalty at $\alpha=1$ shows that it hits the value of $7.0744$ dB for $L\geq 4$, while the lower bound in this case is $\frac{16}{\pi}\approx 7.0697$ dB. The numerical analysis of the limiting energy penalty is considerably simplified in this region of $\alpha$ by the (numerical) observation that the macroscopic parameter $\chi_{1}$ approaches $0$ as $\alpha\rightarrow 1$ (although it stays strictly positive). The small $\chi_{1}$ approximation of the equations employed to calculate the limiting 1RSB energy penalty is shortly discussed in Appendix D. The RSB phenomena is demonstrated by considering the limiting energy penalty obtained via the RS approximation, as stated by Proposition 4.3 (the explicit expression for the current example is given in [24, Eq. (26)]). As shown in Figure 4, the RS approximation fails to predict the limiting energy penalty for $\alpha>0.3$, and in fact it even violates the lower bound (6-14) for $\alpha>0.55$.

The better accuracy of 1RSB is also visible looking at the zero-temperature entropy. We can analytically evaluate Proposition 4.5 in the case of a Gaussian $\boldsymbol{H}$, which becomes

The entropy for both the RS and 1RSB approximations for a relaxation level $L=2$ are shown in Figure 5. Although the 1RSB solution of the above model also has negative zero-temperature entropy, it is much closer to zero, corresponding to a much weaker instability, and approaches zero as $\alpha\rightarrow 1$. In contrast, the RS entropy drifts away from zero as $\alpha\rightarrow 1$.

The limiting conditional probabilities of (6-13) are plotted in Figure 6,

as well as the empirical conditional probabilities, based on the Monte Carlo simulations employed to produce the energy penalties of Figure 4. The results correspond to a relaxation level of $L=2$, and focus on the real part of the extended alphabet points, given that the real part of the original QPSK constellation point satisfies $\Re\left\{u\right\}=1$ (recall the decoupling of the real and imaginary parts implied by (6-13)). The simulation results exhibit again a good match to the limiting analytical 1RSB prediction. It is also clearly demonstrated that, when the system load $\alpha$ is low, hardly any relaxation is required, while the probability of using symbols from the extended alphabet set increases as $\alpha$ approaches unity.

Convex Precoding: An RS Example

This section is devoted to another alphabet relaxation scheme, also introduced in for QPSK signaling. The key feature of this relaxation scheme is that the extended alphabet set is continuous and convex, allowing for an efficient solution to the corresponding quadratic programming problem of minimizing the energy penalty. Convex optimization problems are generally believed not to exhibit replica symmetry breaking . In certain special cases this has been shown explicitly . Furthermore, as will be demonstrated in the sequel, the replica symmetric solution for this alphabet relaxation scheme agrees well with numerical simulations and thus considerably simplifies the analysis of the limiting regime.

the relaxed alphabet subsets are defined by

The alphabet relaxation scheme is depicted in Figure 7, and it is referred to henceforth as convex relaxation for QPSK (CR-QPSK).

The RS approximation for the limiting energy penalty with the CR-QPSK relaxation scheme is obtained through Proposition 4.3, and it is given by the solution to the following fixed point equation [24, Eq. (30)]

Note that (7-3) yields finite energy penalties for all loads $0\leq\alpha<2$. Although loads greater than unity imply that the matrix $\boldsymbol{HH}^{\dagger}$ in (5-1) is singular, this does not lead to interference at the receivers in the large system limit, as shown rigorously in .

Numerical results for the limiting energy penalty of CR-QPSK are plotted in Figure 8.

Empirical results based on Monte Carlo simulations are provided as well. These results were obtained by fixing the number of users to $K=32$, and averaging over 1000 channel realizations. The results exhibit an excellent match to the limiting RS analytical results, thus supporting the validity of the RS approximation. The corresponding results for the discrete lattice-based alphabet relaxation scheme of Section 6 are also provided for the sake of comparison, and it is clearly observed that in terms of the limiting energy penalty, the discrete scheme is superior to the CR-QPSK scheme for all $\alpha\in(0,1]$. The limiting energy penalty difference approaches its maximum value of $2.41\,$dB at $\alpha=1$. As will be shown in Section 8, however, the comparison becomes more subtle when spectral efficiency is investigated.

The RS approximation of the limiting conditional distribution of the precoder outputs is obtained using Proposition 4.4. The idea here is to start from a discretized version of the continuous CR-QPSK relaxed alphabet set, and obtain the limiting conditional distribution of each relaxed alphabet point using (4-17). The final step is then to take the limit as the areas of the Voronoi cells corresponding to each such point vanish. Using this approach, while restricting the discussion to a Gaussian $\boldsymbol{H}$ and focusing for convenience on the QPSK constellation point $u=1+j$, one gets the corresponding conditional probability density function (pdf)

where we decompose the complex argument as $x\triangleq x_{\textsf{re}}+jx_{\textsf{im}}$, $\mathcal{U}(x)$ denotes the unit step function, $\bar{\mathscr{E}}^{{\textsf{cr-qpsk}}}$ denotes the limiting energy penalty of the CR-QPSK scheme obtained from (7-3), and the constant $Q_{1}$ is defined as

The conditional pdf given the rest of the QPSK constellation points (i.e., $u\in\{-1+j,$ $-1-j,1-j\}$) is obtained in an analogous manner, while considering the full symmetry of the extended constellation.

Returning to (7-4), note that this pdf contains masses on the boundaries of $\mathscr{B}_{1+j}$, and in particular a mass point at the original QPSK constellation point (i.e., $x=1+j$). Plots that demonstrate this behavior of the pdf as a function of $\alpha$ are provided in Figure 9. The upper left plot shows the weight of the mass point at $x=1+j$, as a function of $\alpha$ (corresponding to $Q_{1}^{2}$). The lower left plot shows the pdf mass on the lower boundary of the extended alphabet subset $\mathscr{B}_{1+j}$ (i.e., when the imaginary part of the precoder’s output is fixed to $x_{\textsf{im}}=j$). The plots on the right show the pdf on the interior of $\mathscr{B}_{1+j}$, for $\alpha=0.1$ (upper right) and $\alpha=0.9$ (lower right). The increase in probability of using extended alphabet points as the system load increases, is clearly demonstrated in the figure.

Additional numerical results comparing the analytical RS approximation for the pdf to empirical simulation results are shown in the upper left plot of Figure 9 and in Figure 10.

The upper left plot of Figure 9 compares the probability mass at $x=1+j$ to corresponding simulation results for $K=32$ (averaged over 1000 channel realizations). The corresponding comparison for the cumulative distribution function (CDF) of $x_{\textsf{re}}$, given that $\Re\left\{u\right\}=1$, is shown in Figure 10. The left plot shows the CDF for $\alpha=0.7442$ (i.e., for $N=43$), while the right plot shows the results for unit load. As observed, all empirical results exhibit a very good match to the analytical RS approximation, further supporting the validity of the RS analysis for the CR-QPSK scheme.

Spectral Efficiency Comparison

The two previous sections focused on the transmitting end of the system, and investigated the limiting behavior of the precoder output while employing two particular alphabet relaxation schemes. In the following we turn to investigate the limiting behavior of the system as a whole, by considering the normalized spectral efficiency in view of the analysis of Section 5. Accordingly, we restrict the discussion to a ZF front-end and a Gaussian $\boldsymbol{H}$, and apply Proposition 5.1 to obtain the spectral efficiencies of the discrete lattice-based alphabet relaxation scheme and of CR-QPSK.

Starting with the discrete scheme, the spectral efficiency is obtained by incorporating (4-12) and (6-13) into (5-4)–(5-8). The observation made in Section 6 regarding the independence of the real and imaginary parts of the precoder’s output, leads to the following conclusion. The achievable rate in (5-7) for QPSK input can be obtained by treating QPSK signaling as two independent corresponding BPSK signaling settings. Accordingly, the conditional precoder output probabilities, given a real BPSK input of $u=1$, are given by (cf. (6-13))

where we set $\mathscr{B}_{1}=\left\{c_{k}\right\}_{k=1}^{L}$, and the conditional probabilities given $u=-1$ can be immediately obtained from symmetry considerations. It is then straightforward to show that the corresponding spectral efficiency is given by

The spectral efficiency with QPSK input is then obtained through the relation $C^{\textsf{qpsk}}(\frac{E_{b}}{N_{0}})=2C^{{\textsf{bpsk}}}(\frac{E_{b}}{N_{0}})$, while substituting

where we decomposed the complex argument as $y\triangleq y_{\textsf{re}}+jy_{\textsf{im}}$, and the real argument function $Q_{2}(\xi)$ is defined as

The marginal distribution of the equivalent single user channel output $y$ is given by

Finally, following (5-8) and accounting for the inherent symmetry in (8-4) and (8-6), the spectral efficiency of the CR-QPSK scheme is given by

Comparative numerical spectral efficiency results are plotted in Figure 11. The figure shows the spectral efficiencies of the discrete extended alphabet relaxation scheme (while taking $L=2$), and of the CR-QPSK scheme, as well as the spectral efficiency of linear ZF precoding for Gaussian and QPSK input (see (5-15) and (5-17), respectively), and the spectral efficiency of GTHP with QPSK input (following (F-32)–(F-33)). The spectral efficiencies were evaluated for the optimum choice of the system load $\alpha$. The optimum load is a function of $\frac{E_{b}}{N_{0}}$ and shown in Figure 12.

In Figure 11, the DPC spectral efficiency (5-12) is also provided for comparison, evaluated both for $\alpha=1$, and for $\alpha\rightarrow\infty$ (specifying the ultimate performance). The optimization with respect to $\alpha$ emphasizes its role as a crucial system design parameter, facilitating the proper working point for each transmission scheme, per each $\frac{E_{b}}{N_{0}}$. It also naturally translates to a practical scheduling scheme, specifying the desired number of simultaneously active scheduled users per transmit antenna (see, e.g., ).

The results indicate that nonlinear precoding can provide significant performance enhancement for medium to high $\frac{E_{b}}{N_{0}}$ values. The discrete lattice-based relaxation scheme is shown to outperform linear ZF with QPSK input for $\frac{E_{b}}{N_{0}}>3.43\,\text{dB}$. The beneficial effect of the lattice relaxation scheme becomes more pronounced, the more the spectral efficiency approaches the upper limit of $2$ bits/sec/Hz per transmit antenna. For example, a spectral efficiency of $1.75$ bits/sec/Hz can be obtained with lattice relaxation already at $\frac{E_{b}}{N_{0}}\approx 7.66\,\text{dB}$, whereas linear ZF requires additional $7.26\,\text{dB}$ for the same spectral efficiency. In fact, the QPSK-based lattice precoding scheme is shown to marginally outperform linear ZF with Gaussian input for $4.19\,\text{dB}<\frac{E_{b}}{N_{0}}<7.26\,\text{dB}$. The lattice relaxation scheme also outperforms GTHP for all $\frac{E_{b}}{N_{0}}$ values, becoming more effective for medium to high $\frac{E_{b}}{N_{0}}$ (for example, GTHP needs $2.93\,\text{dB}$ more energy per bit to achieve $1.75$ bits/sec/Hz) Note that in general the modulo-receiver employed by GTHP induces poor performance in the low spectral efficiency region (see Appendix F).. The gap from the DPC upper bound is, however, still essentially retained ($4.49\,\text{dB}$ at $1.75$ bits/sec/Hz, considering DPC with $\alpha=1$, to make a fairer comparison).

As for the CR-QPSK scheme, Figure 11 shows that it also provides a considerable performance enhancement over linear ZF with QPSK input. It is outperformed by the lattice relaxation scheme for 4.38 dB $<\frac{E_{b}}{N_{0}}<$ 9.40 dB. It performs better at low values of $\frac{E_{b}}{N_{0}}$, and in fact it even negligibly outperforms linear ZF with Gaussian input in the low $\frac{E_{b}}{N_{0}}$ region. Moreover, unlike the discrete scheme, CR-QPSK outperforms linear ZF precoding (with QPSK input) for all $\frac{E_{b}}{N_{0}}$ values. Furthermore, it outperforms lattice relaxation in the high $\frac{E_{b}}{N_{0}}$ region, since it allows for loads up to $\alpha<2$ and therefore its spectral efficiency is no longer upper bounded by 2 bits/sec/Hz, but rather by 4 bits/sec/Hz per transmit antenna. Though, the convergence to the limiting spectral efficiency of 4 bits/s/Hz at high $\frac{E_{b}}{N_{0}}$ is rather slow. CR-QPSK also outperforms GTHP for all $\frac{E_{b}}{N_{0}}$ values, but the advantage is more significant for high $\frac{E_{b}}{N_{0}}$, where overloading is employed. These results are of particular interest since the CR-QPSK scheme lends itself to efficient implementation, whereas the discrete relaxation scheme involves the solution of an NP-hard optimization problem. It is also important to note that, as shown in Figure 8, the CR-QPSK scheme is always inferior to the lattice relaxation scheme in terms of the limiting energy penalty. Hence, in view of the observations made here, one can conclude that restricting the analysis to the energy penalty alone provides only limited insight into the behavior of large coded systems, as it essentially focuses only on the transmitter, while ignoring the impact of the nonlinear precoding scheme on the receiver.

Concluding Remarks

The replica symmetry breaking ansatz of statistical physics was employed in this paper to investigate the large system limit behavior of nonlinear precoding for the MIMO Gaussian broadcast channel based on linear zero-forcing and alphabet relaxation. For lattice relaxations, the replica symmetric ansatz was shown to yield misleading results for system loads greater than approximately 0.3 while the one-step replica symmetry breaking ansatz provides sensible results for any load. For exact results, however, multiple-step replica symmetry breaking must be considered.

Introducing a nonlinear superchannel comprising the actual channel and the precoder, allows for a Markov chain description of an individual user’s channel. This enables the calculation of mutual information and spectral efficiency in the large system limit. While convex QPSK relaxations are significantly outperformed by lattice relaxations in terms of transmitted energy per bit, they are very competitive when combined with strong error-correction coding as shown by the spectral efficiency analysis. Except for medium signal-to-noise ratios, they are superior to lattice relaxations. Both schemes were shown to outperform Tomlinson-Harashima precoding with QPSK input for all signal-to-noise ratios.

The combination of polynomial complexity and high spectral efficiency makes convex alphabet relaxation schemes, as introduced in , a promising alternative to the NP-hard lattice relaxations due to their polynomial complexity, and to Tomlinson-Harashima precoding due to their superior performance. The results motivate the search for convex schemes amenable to efficient implementation. Additional examples for extended alphabets are currently investigated, see for preliminary results. Note, however, that the problem of finding the optimum precoding scheme that maximizes the spectral efficiency in this framework is not at all trivial, as the corresponding equivalent channel statistics depend, in this setting, on the choice of input distribution and extended alphabet sets.

Appendix A Higher RSB Orders

using the constants $\left\{q_{r},p_{r}^{(1)},\dots,p_{r}^{(r)},\chi_{r},\mu_{r}^{(1)},\dots,\mu_{r}^{(r)}\right\}$. The limit as $r\to\infty$ is called full RSB and gives the exact solution to the problem . The particular temperature-dependent scaling of some parts of $\boldsymbol{Q}$ is used to evaluate the free energy at zero temperature without getting divergent terms. If a finite temperature is of interest a different scaling may be considered. Plugging (A-1) into (3-25), while exploiting the particular structure of $\boldsymbol{Q}$, we find:

For any temperature, the energy for $r$-step RSB is

where $R^{\prime}(\cdot)$ denotes the derivative of the function $R(\cdot)$.

In order to proceed to full RSB, the limit $r\to\infty$ must be taken. Naively, one might think this would make the sums in (A-2) diverge. However, the macroscopic parameters are determined by the saddle point equations, which guarantee that the sums stay finite through decreasing the macroscopic parameters. Thus, we introduce a continuum of macroscopic parameters $\mu(x)$ and $p(x)$, taken over $0\leq x\leq 1$, such that

Accordingly, we find for the energy in the limit $r\to\infty$

Using integration by parts, (A-10) simplifies to

The functions $p(x)$ and $\mu(x)$ must be determined by the respective saddle point equations.

Appendix B Proofs for 1RSB

The joint distribution of the entries of the vector $\boldsymbol{x}$, conditioned on both the input vector $\boldsymbol{u}$ and the channel transfer matrix $\boldsymbol{H}$, is given for a non-zero temperature by the Boltzmann distribution

where $\mathcal{Z}$ is the partition function defined in (3-4). Taking the limit $\beta\rightarrow\infty$ (zero temperature), the denominator in (B-1) is dominated by its maximum value term, and the limiting joint distribution of the entries of $\boldsymbol{x}$, conditioned on all inputs, converges to the Dirac measure at $\argmin_{\boldsymbol{x}\in\mathscr{B}_{\boldsymbol{u}}}\boldsymbol{x}^{\dagger}\boldsymbol{J}\boldsymbol{x}$, corresponding to the minimum normalized energy penalty, as given by Proposition 4.1.

To prove Proposition 4.1, we will need to evaluate the free energy averaged over all realizations of $\boldsymbol{u}$ and $\boldsymbol{H}$. For future convenience, we also include the dummy variable $h$ and the function $V(\cdot)$ defined in (3-14) and rewrite the free energy as

where $\mathcal{Z}({h};\boldsymbol{u},\boldsymbol{H})$ is given by (3-15). The second equality is a manifestation of the underlying assumption that the coded symbols of all users are drawn randomly and independently of the channel transfer matrix $\boldsymbol{H}$. In view of this formulation, we consider now the limit of the term in the parentheses above

As shown later on, this inner limit is a deterministic quantity, for almost every realization of the input vector $\boldsymbol{u}$. It will hence be concluded that in fact

As indicated earlier, the key tool in the derivation of the above quantity is the replica method of statistical physics, using the identity

and following the outline in Section 3. With that in mind, the quantity $[\mathcal{Z}(h;\boldsymbol{u},\boldsymbol{H})]^{n}$ is regarded as consisting of $n$ identical replicas of the original (unnormalized) probability model in the following way

Interchanging the limits of $K\rightarrow\infty$ and $n\rightarrow 0$, the focus is first on the derivation of the limit

Since the first exponential term within the expectation is independent of the channel transfer matrix $\boldsymbol{H}$, $\Xi_{n}$ can be rewritten as

The inner expectation in (B-8) is the Harish-Chandra-Itzykson-Zuber integral (see and references therein), and the objective here is its evaluation for fixed-rank matrices $\sum_{a=1}^{n}\boldsymbol{x}_{a}\boldsymbol{x}_{a}^{\dagger}$, in the large $K$ limit. This problem was recently considered in , and invoking Theorem 1.7 therein, (B-8) can be represented for large $K$ as $o(K)$ is used here to denote quantities that satisfy $\lim_{K\to\infty}o(K)/K=0$.

where $R(w)$ is the $R$-transform of the limiting eigenvalue distribution of the matrix $\boldsymbol{J}$, and $\left\{\lambda_{a}\right\}$ denote the eigenvalues of the $n\times n$ matrix $\beta\boldsymbol{Q}$ with $\boldsymbol{Q}$ defined through Here [53, Theorem 1.7] is applied individually for all given vectors $\left\{\boldsymbol{x}_{a}\right\}$.

Since additive exponential terms of order $o(K)$ have no effect on the results in the limiting regime as $K\rightarrow\infty$, due to the $\frac{1}{K}$ factor outside the logarithm in (B-9) (this shall become clear in the derivation to follow), any such terms are dropped henceforth for notational simplicity.

In order to calculate the summation in (B-9), the procedure employed in is repeated here, and the $Kn$-dimensional space spanned by the replicas is split into subshells by means of (3-24). Assuming $\boldsymbol{Q}^{\dagger}=\boldsymbol{Q}$, $\Xi_{n}$ can be represented as

since the trace is the sum of the eigenvalues,

is the probability weight of the subshell.

Denoting by $\boldsymbol{P}$ the Hermitian matrix with elements $P_{ab}=\boldsymbol{x}_{a}^{\dagger}\boldsymbol{x}_{b}-KQ_{ab}$, this yields

Considering the inner summation in (B-24), then rearranging terms and using (3-14) the expression can be rewritten as

Now, using the underlying assumption that the coded symbols transmitted by different users are i.i.d., one can apply the strong law of large numbers for $K\rightarrow\infty$ to get

where the convergence is in the almost sure sense, for any extended alphabets such that the expectation in (B-30) exists. Note that this observation implies that any randomness due to $\boldsymbol{u}$ in the RHS of (B-11) effectively vanishes at the large system limit, due to the normalization with respect to $K$ outside the logarithm.

The next step in the evaluation of (B-11) is the observation that in the limit as $K\rightarrow\infty$, the integrand therein is dominated by the exponential term with maximal exponent. Therefore, only the subshell that corresponds to this extremal value of the correlation between the vectors $\left\{\boldsymbol{x}_{a}\right\}$ is relevant for the calculation of the integral. Thus, we have at the saddle point

Since the trace is the sum of the eigenvalues, we can write (B-13) as

Furthermore, we observe that also the integrand in (B-28) is dominated by the exponential term with maximal exponent in the limit $K\to\infty$. Thus, at the saddle point we have

We now invoke the 1RSB assumption (3-31) regarding the structure of the matrices $\boldsymbol{Q}$ at the saddle-point that dominate the integral. In a similar manner we set

introducing the macroscopic parameters $f_{1}$, $g_{1}$, and $\varepsilon_{1}$.

With these assumptions one can explicitly obtain the eigenvalues of the matrix $\beta\boldsymbol{Q}\,$ The eigenvalue $(\beta nq_{1}+\mu_{1}p+\chi_{1})$ occurs with multiplicity $1$, the eigenvalue $(\mu_{1}p_{1}+\chi_{1})$ occurs with multiplicity $(\frac{n\beta}{\mu_{1}}-1)$, and the eigenvalue $\chi_{1}$ occurs with multiplicity $(n-\frac{n\beta}{\mu_{1}})$., and $\mathcal{G}(\boldsymbol{Q})$ can be rewritten as

It also follows from the 1RSB assumption that

Due to (B-32), the partial derivatives of

with respect to $q_{1}$, $p_{1}$, and $\chi_{1}$ must vanish as $K\rightarrow\infty$ by definition of the saddle point. Using (B-38) and (B-39) this yields the following set of equations

Solving for $\varepsilon_{1}$, $g_{1}$, and $f_{1}$, while focusing on the limit as $n\rightarrow 0$, one gets

We now rewrite the expression for $M_{k}(f_{1},g_{1},\varepsilon_{1},\mu_{1})$ in (B-40) using the Hubbard-Stratonovich transform and the shortened notation of (4-6)

Due to (B-35), the partial derivatives of

with respect to $f_{1}$, $g_{1}$ and $\varepsilon_{1}$, must also vanish as $K\rightarrow\infty$. This produces the following set of equations (while taking the limit as $n\rightarrow 0$)

The parameter $\mu_{1}$ should also be chosen such that the partial derivative of

with respect to $\mu_{1}$ vanishes. This yields at the limit as $n\rightarrow 0$

Incorporating all previous results, we get that the quantity $\Xi_{n}$ of (B-7) is equal to

where the macroscopic parameters $\left\{f_{1},g_{1},\varepsilon_{1},q_{1},p_{1},\chi_{1},\mu_{1}\right\}$ are obtained from the saddle point fixed-point equations (B-45)–(B-47), (B-53)–(B-55), and (B-57). Now in view of (B-5), the next step in the derivation is to take the limit

which is justified by the observation that $\Xi_{n}$ converges to the same limit for almost every realization of $\boldsymbol{u}$, applying the law of large numbers in (B-30) (see also (B-3) and the discussion that follows).

We note at this point that the energy penalty $\bar{\mathscr{E}}_{\textsf{rsb1}}$ satisfies

and (4-12) can be readily expressed from Proposition A.1 as a function of the macroscopic parameters $\left\{q_{1},p_{1},\chi_{1},\mu_{1}\right\}$. In order to evaluate the energy penalty, it is thus left to derive the fixed point equations that determine these parameters, as given by (4-8)–(4-11), which is obtained by substituting $h=0$ and taking the limit as $\beta\rightarrow\infty$ in (B-53)–(B-55), and (B-57) after back-substitution of (B-51). This completes the proof of Proposition 4.1.

B.2 Proposition 4.2

We derive the limiting conditional distribution of the precoder output $x$ given the input $u$ starting from (3-17). We therefore need to evaluate the derivative of the free energy with respect to $h$. This can be done directly given (B-60). Taking the partial derivative in (B-59) and using (B-51), we get

After taking the limit $\beta\rightarrow\infty$, while applying the saddle point integration rule, we finally get (4-13). This completes the proof of Proposition 4.2.

B.3 Proposition 4.5

We start from (3-10), (3-19), and (3-20) which yield

The partial derivative with respect to $\beta$ above reflects the fact that all implicit dependencies of $\mathscr{F}(\beta)$ on $\beta$ through its dependence on other parameters, e.g., $f_{1},g_{1},\chi_{1},\mu_{1}$, have vanishing derivatives since $\mathscr{F}(\beta)$ is evaluated at a saddle point. Making use of the saddle point equations, we find that

Using the fixed point equations we may re-express the first line as follows:

Plugging this into the above equation and using (B-45)–(B-47), we eventually get the following equation for the zero-temperature entropy

Remarkably the above equation for the entropy holds also for the RS case. To recover the RS structure of the equations above we start with $\mu_{1}/\beta=1$ and $\chi_{0}=\chi_{1}+\mu_{1}p_{1}$. Then, we find that $q_{0}=q_{1}$, $\epsilon_{0}=\epsilon_{1}-\beta g_{1}^{2}$, and $f_{0}=f_{1}$. After that we find the equations to reduce to the RS case analyzed in .

Appendix C Proof of Proposition 3.1

We will now apply (3-9) to express the energy in a compact fashion. We start by considering the representation of the normalized average free energy in terms of $\boldsymbol{Q}$ (see (3-19) and (3-24)), and let us denote this representation, for the sake of clarity, as $\mathscr{F}(\boldsymbol{Q},\beta)$. In general, the replica crosscorrelation matrix $\boldsymbol{Q}$ depends on $\beta$. However, at the saddle point we have (by definition)

Thus, the total derivative in (3-9) becomes a partial derivative at the saddle point, i.e.

Referring to the proof in Appendix B, then with (B-2), (B-5), (B-7), (B-11), and (B-15), while substituting $h=0$, this gives

which is easily shown to be equivalent to (3-25). Furthermore, we get (3-26) by plugging (B-34) into (B-36) while substituting $h=0$.

Appendix D Discrete Lattice Relaxation: Small χ1\chi_{1} Approximation Near Unit Load (1RSB)

This appendix provides an approximate derivation of the 1RSB equations for the discrete lattice-based alphabet relaxation scheme of Section 6, while assuming a Gaussian $\boldsymbol{H}$, and a ZF front-end. The approximation is based on the numerical observation that the macroscopic parameter $\chi_{1}$, employed in the 1RSB ansatz for this setting, approaches zero as the system load gets close to unity. This approximation considerably simplifies the numerical solution of the 1RSB equations in this region of the system load.

For $\alpha=1$, the $R$-transform of $\boldsymbol{J}=\boldsymbol{T}^{\dagger}\boldsymbol{T}$ (see Proposition 4.1) satisfies

Considering the small $\chi_{1}$ regime, we get from (4-2)–(4-4)

Particularizing to the two-dimensional discrete lattice-based alphabet relaxation scheme in concern, one gets from (6-6)

We now rewrite the function $\Theta_{k}(\xi)$ of (6-7) as

and observe the following. Starting with exponential argument, we get

while the arguments of the $Q(\cdot)$ functions satisfiy

Now recall that from the underlying definition of the extended alphabet set, it follows that $c_{k}\geq c_{k-1}\ \forall k$, and it can hence be concluded that

In a similar manner one can observe that the exponential terms in the RHS of (6-8) vanish as $\chi_{1}\rightarrow 0$, and conclude that

In view of the above we can now restate the coupled equations that determine the macroscopic parameters $q_{1}$, $p_{1}$, and $\mu_{1}$ in the following way (cf. (6-9)–(6-12), and note that the equation for determining $\chi_{1}$ can be ignored):

The energy penalty in this case is given by (cf. (4-12))

D.2 Case II: α<1\alpha<1, α→1\alpha\rightarrow 1

In a similar manner to the previous section, we start with the $R$-transform of $\boldsymbol{J}=\boldsymbol{T}^{\dagger}\boldsymbol{T}$, and rewrite it for small $w$, using the Taylor expansion around $w=0$, as

We focus in the following on the regime in which $\alpha\rightarrow 1$, so that $\frac{1}{1-\alpha}\gg 1$, but still $\frac{1}{1-\alpha}\ll\frac{1}{\chi_{1}}$. It hence follows that

Particularizing again to the two-dimensional discrete lattice alphabet relaxation scheme for QPSK signaling, it follows from (6-6) that

Next, the arguments of the $Q(\cdot)$ functions in (6-7) satisfy

Finally, note that $\int_{0}^{\mu_{1}p_{1}}R(-w)\,dw$ exists for $\alpha<1$, and the approximation

is employed to derive the three coupled equation that determine the macroscopic parameters $q_{1}$, $p_{1}$, and $\mu_{1}$. The three equations are thus

The expression for the energy penalty is given by

We also note that the exact expressions for the $R$-transform and its derivative were employed for the purpose of producing more accurate numerical results, while using this small $\chi_{1}$ approximation for $\alpha<1$.

Appendix E Proof of Lemma 4.6

The Stieltjes transform of the probability distribution $F(x)$ is defined by

In terms of the Stieltjes transform, the $R$-transform is defined as

where $m^{-1}(s)$ denotes the inverse function of $m(s)$ with respect to composition, i.e., $m(m^{-1}(s))=s$.

We start with the observation that the derivative of the Stieltjes transform is lower bounded by its square

by means of Jensen’s inequality, with equality if and only if the distribution $F(x)$ is a single mass point. Next, we consider the derivative of the $R$-transform. Letting $w=m(s)$, it follows that

with equality if and only if the distribution $F(x)$ is a single mass point. Lemma 4.6 then follows immediately.

Appendix F Spectral Efficiency of Generalized Tomlinson-Harashima Precoding

For the sake of comparison, we review here the derivation of the spectral efficiency of generalized Tomlinson-Harashima precoding (GTHP), which is another practical alternative to the capacity achieving DPC. The approach is based on inflated lattice strategies, and borrows ideas from the recent analysis of pulse amplitude modulation (PAM) in . The spectral efficiency is derived following (see also ), while employing successive encoding using the inflated lattice strategy at each stage, where the signals of previously encoded users are treated as causally known interference. We consider here the “canonical” channel model as in (2-1), and note that a comparative analysis of other variants of GTHP can be found, e.g., in .

The underlying idea of the scheme considered here is first to induce a “triangular” channel structure using the $LQ$-factorization of the channel transfer matrix. Assuming $\boldsymbol{H}$ is full rank, we denote

where $\boldsymbol{L}_{[K\times K]}$ is lower triangular with positive diagonal entries and $\bar{\textbf{Q}}_{K\times N}$ has orthonormal rows. The transmitted signal is then given by

where $\boldsymbol{x}$ is the nonlinear precoder’s output (cf. (2-3)). The signal received by the $k$th user is thus given by

where $\left\{L_{ij}\right\}$ denote the entries of $\boldsymbol{L}$ and $x_{i}$ is the nonlinear precoder’s output that corresponds to user $i$. Normalizing both sides of the equation by $L_{kk}$, we get the following equivalent channel

where we denote the multiuser interference experienced by user $k$ as $s_{k}\triangleq\sum_{j=1}^{k-1}\frac{L_{kj}}{L_{kk}}x_{j}$, and $\breve{n}_{k}$ is a zero-mean circularly symmetric complex Gaussian noise with variance $\frac{\sigma^{2}}{L_{kk}^{2}}$.

In the GTHP setting, instead of DPC (as employed at this point, e.g., by the “zero-forcing dirty-paper” scheme of ), we follow for each user the THP-type strategy described in for canceling the interference due to previously encoded users. This strategy, which applies for canceling causally known interference, leads in the broadcast setting to a considerably reduced complexity as it involves only scalar quantizations (as opposed to vector quantizations in the noncausal case, see therein). To make a fair comparison to the precoding schemes discussed in Sections 6-7, we particularize here to the case in which the information bearing signal takes on binary values per each dimension (so that the total spectral efficiency for quadrature modulation, as a function of $\frac{E_{b}}{N_{0}}$, is twice as much as the one obtained for binary input). The spectral efficiency for continuous input is derived as well for completeness. The basic transmission scheme is reviewed first, while considering real channels.

The underlying real channel model is given by

where $x$ is subject to an average power constraint $P_{x}$, $n$ is a zero-mean AWGN with variance $P_{n}$, and $s$ is an interference signal which is known causally at the transmitter (i.e., at the current time instance), but not at the receiver. This channel model is also referred to in the literature as the “dirty-tape” model . Consider the one-dimensional lattice

Let $\mathcal{V}=[-\Delta,\Delta)$ denote the basic Voronoi region of $\Lambda$. Let $d$ be a dither signal uniformly distributed over $\mathcal{V}$. Under a common randomness assumption, this dither signal is assumed to be available at the receiver as well.

Starting with continuous information bearing signals, then by the GTHP scheme the transmitter sends the signal

Effectively, the induced channel is equivalent to (see , Lemma 6)

and we note that the dither signal ensures that $x$ is uniformly distributed over the Voronoi region, and is independent of either the information bearing signal $v$, or the noise $n$.

The capacity of this channel is achieved by a uniform input distribution over the Voronoi region, $v\sim\Unif\left\{\mathcal{V}\right\}$, for which the relation between the lattice constant $\Delta$ and the transmit power $P_{x}$ is given by $\Delta=\sqrt{3P_{x}}$. The corresponding achievable rate is equal to the input-output mutual information of the equivalent channel (F-9)

The entropy of the effective noise is derived via the following observation. Denoting the “self-noise” term by

The pdf of the effective noise (F-10) is thus given by

Turning to discrete input with M-pulse amplitude modulation (M-PAM) (representing the information bearing signals), the setting is equivalent to the case in which the continuous information bearing signal considered above is quantized (cf. ). Instead of (F-7), the channel input is now given by

where $\mathcal{Q}(\cdot)$ denotes the nearest-neighbor uniform quantizer with step size $\Delta$ , and $v$ is assumed to be uniformly distributed over the Voronoi region. We note here that this transmission scheme differs from the one considered in , where the channel input is quantized to comply with an M-PAM constellation (see therein). Note also that as in the continuous setting, due to the dither signal, the channel input $x$ is still uniformly distributed over the Voronoi region. The effective channel can now be represented in the form (cf. (F-9))

where the effective noise is still given by (F-10). Restricting this review to the case of binary information bearing signals per dimension, the channel input signal is limited to the interval $\mathcal{V}=[-\Delta,\Delta)$, while the quantized information bearing signal is obtained using

For consistency we retain the relation $\Delta=\sqrt{3P_{x_{\mathcal{Q}}}}$ .

The achievable rate for binary input is given again by the mutual information

Note that the pdf of the random quantity inside the modulo function in (F-18) is given by

Hence, the pdf of the equivalent channel output $y^{\prime}_{\mathcal{Q}}$ is equal to

and the achievable rate of (F-20) can be rewritten as

The above principles can now be applied to the channel in (F-3), where the transmitter pre-cancells using the GTHP scheme, per each transmitted symbol, the interference due to the corresponding symbols of previously encoded users. Using (F-11) and (F-20), the achievable rate of the $k$th user can be obtained by substituting $P_{x}=L_{kk}^{2}\,{\textsf{snr}}$ for continuous input, and $P_{x_{\mathcal{Q}}}=L_{kk}^{2}\,{\textsf{snr}}$ for the binary setting, yielding, respectively, for real channels

To complete the analysis, it is left to derive the normalized spectral efficiency of GTHP in the large system limit. This is obtained using the following observation (see [49, Lemma 3]).

Let $\boldsymbol{H}$ be a $K\times N$ random matrix, having i.i.d. circularly symmetric zero-mean entries with variance $\frac{1}{N}$ and finite fourth moment, and let $\boldsymbol{H}^{(k)}$, $k<K$, denote the matrix constructed by striking out the last $K-k$ rows of $\boldsymbol{H}$. Then

and for $k,K,N\to\infty$, s.t. $\frac{K}{N}\to\alpha<\infty$ and $\frac{k}{K}\to\nu\in[0,1)$, it follows that

where for the continuous case we substitute

and for the case of discrete binary information bearing signals we substitue

The spectral efficiency for QPSK modulation satisfies (following the convention in ):

where $C^{{\textsf{gthp}}}_{\textsf{bpsk}}({\textsf{snr}})$ is given by (F-30) and (F-32), and it can be expressed as a function of $\frac{E_{b}}{N_{0}}$ through (5-9). An analogous result for the case of continuous input can be readily obtained using (F-31). Both spectral efficiencies can be optimized with respect to the choice of the system load $\alpha$.

References