# Simple Framework of ZF Precoding Analysis for Full-Dimensional Massive MIMO
Systems

Harsh Tataria* and Mansoor Shafi†
*Ericsson AB, Lund, Sweden
†Spark New Zealand, Wellington, New Zealand
e-mail<EMAIL_ADDRESS>and<EMAIL_ADDRESS>
###### Abstract

We provide a simple, yet general and accurate framework for analyzing zero-
forcing (ZF) precoding performance for a full-dimensional massive multiple-
input multiple-output system. Exploiting an order two Neumann series, our
framework approximates the expected (average) ZF signal-to-noise ratio and
ergodic sum spectral efficiency, while catering for a finite multipath
propagation model, as well as correlated and uncorrelated user equipment (UE)
positions in both azimuth and elevation domains. The analysis provides clear
insights on the influence of propagation and system parameters on ZF
performance. We identify how far UEs must lie in azimuth and elevation for the
ZF precoder to approach uncorrelated channels. Our framework is useful for
predicting the ZF performance degradation arising from channel correlation.
For optimal performance, UEs could be separated by around 30∘ in azimuth and
15∘ in elevation - conditions which are difficult to meet in reality.

##  I Introduction

Massive multiple-input multiple-output (MIMO) is an integral part of fifth-
generation (5G) wireless systems. Due to hundreds of electronically steerable
antenna elements at the cellular base station (BS), aggressive spatial
multiplexing to tens of user equipments (UEs) is possible in the same time-
frequency resource [1]. For sub-6 GHz massive MIMO operation, digital
precoding methods provide a means to suppress multiuser interference. Most
current commercial systems employ techniques such as zero-forcing (ZF)
precoding to _orthogonalize_ multiuser channels, such that scheduled UEs are
placed in the null space of one another [1, 2]. Naturally, geographical
location of UEs and scatterers play a vital role in conditioning the UE null
spaces. If UEs are closely spaced in either azimuth or elevation, ZF suffers
from noise enhancement resulting in lower signal-to-noise ratios (SNRs).

Nonetheless, fundamentally, “how close is close” in terms of UE positions such
that the average ZF SNR to a given UE can be maximized remains to be seen.
Naturally, this has a direct consequence on the achievable ergodic sum
spectral efficiency over all UEs in the system. Since most studies on massive
MIMO (see [3, 4] and references therein) neglect the presence of the elevation
domain, it remains to be seen whether closely separated UEs in the azimuth and
elevation domains contribute to the equal ZF performance difference or not,
relative to widely separated UEs which result in uncorrelated channels.
Important design guidelines on optimal azimuth and elevation separations of
scheduled UEs to approach the performance of uncorrelated channels remains to
be checked.

In order to answer these important questions, the focus of this paper is on
providing a simple, accurate and general framework for the performance
analysis of ZF precoding in full-dimensional [2] massive MIMO systems. Despite
the vast literature on this topic [1, 2, 3, 4, 5], an insightful framework
catering for correlated and uncorrelated finite multipath propagation channels
remains missing. To demonstrate its simplicity, our framework utilizes a
second-order Neumann series to approximate the expected (average) ZF SNR and
ergodic sum spectral efficiency. We provide clear insights on the influence of
propagation and system parameters on aforementioned performance metrics. Due
to its generality, our framework is capable of handling multiple BS array
geometries and estimating the ZF performance degradation as a result of sub-
optimal UE positions in the azimuth and elevation. As such, we provide
guidelines on the effectiveness of ZF precoding with respect to different UE
separations. To our knowledge, a study of this type has been missing from
prior literature.

Notation. Upper and lower boldface letters represent matrices and vectors. The
$M\times{}M$ identity matrix is denoted as $\mathbf{I}_{M}$, while
$M\times{}M$ diagonal matrix $\mathbf{X}$ is referred by
$\mathbf{X}=\textrm{diag}(x_{1},x_{2},\dots,x_{M})$. The $(i,j)$-th entry of
the matrix X is denoted by $(\mathbf{X})_{i,j}$, while $\mathbf{a}\hskip
1.0pt\otimes\hskip 1.0pt\mathbf{b}$ denote the Kronecker product of vectors
$\mathbf{a}$ and $\mathbf{b}$. Transpose, Hermitian transpose, inverse and
trace operators are denoted by $(\cdot)^{T}$, $(\cdot)^{H}$, $(\cdot)^{-1}$
and $\textrm{Tr}[\cdot]$. Moreover, $||\cdot||_{F}$ and $|\cdot|$ denotes the
Frobenius and scalar norms. We use
$\mathbf{h}\sim\mathcal{CN}(\mathbf{m},\mathbf{R})$ to denote a complex
Gaussian distribution for $\mathbf{h}$ with mean $\mathbf{m}$ and covariance
matrix $\mathbf{R}$. Finally, $\mathbb{E}\\{\cdot\\}$ denotes the statistical
expectation.

##  II System Model

We consider the downlink of a single-cell full-dimensional massive MIMO
system. The BS is located at the center of a circular cell with a coverage
radius of $R_{\textrm{c}}$ meters (m), and is equipped with an array of $M$
transmit antennas. The BS simultaneously serves $L$ single-antenna UEs
$\left(M\gg{}L\right)$ in the same time-frequency resource. Uniform power
allocation to the $L$ UEs is assumed, and the BS has a constraint on the total
transmit effective isotropic radiated power (EIRP) of $P_{\textrm{EIRP}}$
(quoted often in Watts or dB). Knowledge of the channel is assumed (motivated
later) with narrowband transmission. The $1\times{}M$ channel from the BS
array to UE $\ell$ is given by

$\mathbf{h}_{\ell}=\sum\limits_{p=1}^{N_{\textrm{P}}}\gamma_{\ell,p}\hskip
2.0pt\mathbf{a}_{\ell,p}\hskip 1.0pt.\vspace{-5pt}$ (1)

Here $N_{\textrm{P}}$ is the total number of contributing multipath components
(MPCs) and $\gamma_{\ell,p}\sim\mathcal{CN}\hskip 2.0pt\big{(}0\hskip
1.0pt,1/N_{\textrm{P}}\big{)}$ is the complex gain on the $p$–th MPC.
Moreover, and $\mathbf{a}_{\ell,p}$ is the $1\times{}M$ far-field array
steering vector of the BS array. In the case of a uniform linear array (ULA)
positioned on the $z$-axis at the BS,

$\displaystyle\mathbf{a}_{\ell,p}=\mathbf{a}\left(\phi_{\ell,p}\right)=$
$\displaystyle\left[\hskip 3.0pt1\hskip
8.0pte^{-j2\pi{}\frac{d_{z}}{\lambda_{f}}\cos\hskip
2.0pt\left(\phi_{\ell,p}\right)}\right.$ $\displaystyle\hskip
20.0pt\dots\hskip
4.0pt\left.e^{-j2\pi\frac{d_{z}}{\lambda_{f}}\left(M-1\right)\cos\hskip
2.0pt\left(\phi_{\ell,p}\right)}\right],$ (2)

where $d_{z}$ is the inter-element spacing along the $z$-axis of the array,
$\lambda_{f}$ is the wavelength associated with carrier frequency, $f$, and
$\phi_{\ell,p}$ is the azimuth angle-of-departure (AOD) of the $p$-th MPC.
Note that $\gamma_{\ell,p}$, and $\phi_{\ell,p}$ are assumed to be
uncorrelated over any of their indices, and with each other. Furthermore,
$\phi_{\ell,p}\in\left[\bar{\phi_{\ell}}-\Delta^{\textrm{az}}_{\ell},\bar{\phi_{\ell}}+\Delta^{\textrm{az}}_{\ell}\right]$,
where $\bar{\phi}_{\ell}$ is the azimuth line-of-sight (LOS) angle to UE
$\ell$ and $\Delta^{\textrm{az}}_{\ell}$ is the $\ell$-th UEs azimuth angular-
spread-of-departure (ASD). In the case when the BS is equipped with a uniform
planar array (UPA) on the $x$-$z$ plane, we consider $M_{x}$ antenna elements
along the $x$-axis (configured as a ULA), and $M_{z}$ antenna elements along
the $z$-axis (also configured as a ULA). The corresponding inter-element
separations in $x$ and $z$ directions are denoted by $d_{x}$ and $d_{z}$,
respectively. The total number of BS antennas at the UPA is retained as
$M=M_{x}M_{z}$. Different from the ULA, the UPA is able to exploit the
elevation domain of the channel. To this end, we denote the elevation AOD
belonging to the $p$-th MPC from the BS to the $\ell$-th UE as
$\theta_{\ell,p}$, and note that
$\theta_{\ell,p}\in\left[\bar{\theta}_{\ell}-\Delta^{\textrm{ze}}_{\ell},\bar{\theta}_{\ell}+\Delta^{\textrm{ze}}_{\ell}\right]$.
Here $\Delta^{\textrm{ze}}_{\ell}$ is the $\ell$-th UE specific elevation ASD,
and $\bar{\theta}_{\ell}$ is the elevation LOS angle to UE $\ell$. In light of
the above, the $1\times{}M$ array steering response vector for the $p$-th MPC
from the BS to UE $\ell$ is given by

$\mathbf{a}_{\ell,p}=\mathbf{a}\left(\phi_{\ell,p},\theta_{\ell,p}\right)=\mathbf{a}_{x}\left(\phi_{\ell,p}\hskip
3.0pt,\theta_{\ell,p}\right)\otimes\mathbf{a}_{z}\left(\phi_{\ell,p}\hskip
3.0pt,\theta_{\ell,p}\right),$ (3)

where

$\displaystyle\mathbf{a}_{x}\left(\phi_{\ell,p}\hskip
3.0pt,\theta_{\ell,p}\right)=$ $\displaystyle\left[\hskip 2.0pt1\hskip
8.0pte^{-j2\pi\frac{d_{x}}{\lambda_{f}}\sin\hskip
1.0pt\left(\theta_{\ell,p}\right)\cos\hskip
1.0pt\left(\phi_{\ell,p}\right)}\hskip 3.0pt\right.$ $\displaystyle\dots\hskip
3.0pt\left.e^{-j2\pi\frac{d_{x}}{\lambda_{f}}\left(M_{x}-1\right)\sin\hskip
1.0pt\left(\theta_{\ell,p}\right)\cos\hskip
1.0pt\left(\phi_{\ell,p}\right)}\right],$ (4)

and

$\displaystyle\mathbf{a}_{z}\left(\phi_{\ell,p}\hskip
3.0pt,\theta_{\ell,p}\right)=$ $\displaystyle\left[\hskip 1.0pt1\hskip
8.0pte^{-j2\pi\frac{d_{z}}{\lambda_{f}}\sin\hskip
1.0pt\left(\theta_{\ell,p}\right)\sin\hskip
1.0pt\left(\phi_{\ell,p}\right)}\right.$ $\displaystyle\left.\hskip
3.0pt\dots\hskip
3.0pte^{-j2\pi\frac{d_{z}}{\lambda_{f}}\left(M_{z}-1\right)\sin\hskip
1.0pt\left(\theta_{\ell,p}\right)\sin\hskip
1.0pt\left(\phi_{\ell,p}\right)}\right].$ (5)

For maximum clarity, we defer the discussion of the aforementioned parameters
till Sec. IV.

Remark 1. Two reasons can be pointed out to justify the rather impractical
assumption of perfect channel knowledge. Firstly, unlike previous studies, the
central focus of the work is to devise a simple, yet accurate, analysis
framework to gain insights into the behavior of expected ZF SNR and ergodic
sum spectral efficiency when considering full-dimensional massive MIMO
systems. In contrast to prior studies, spatial channel parameters in both
azimuth and zenith domains are utilized to capture their variations across
multiple terminals. Under such a heterogeneous scenario, it is extremely
difficult to make analytical progress without perfect channel knowledge.
Secondly, it is worth noting that the results obtained from the subsequent
analysis and evaluations can be treated as a useful upper bound on the
performance which may be seen in practice with estimated channels.

The received signal at UE $\ell$ is given by

$y_{\ell}=\sqrt{\frac{P_{\textrm{EIRP}}\hskip 1.0pt\beta_{\ell}}{\eta}}\hskip
1.0pt\mathbf{h}_{\ell}\hskip 1.0pt\mathbf{g}_{\ell}{}\hskip
1.0pts_{\ell}+\sum\limits_{\begin{subarray}{c}i=1\\\
i\neq{}\ell\end{subarray}}^{L}\sqrt{\frac{P_{\textrm{EIRP}}\hskip
1.0pt\beta_{\ell}}{\eta}}\hskip 1.0pt\mathbf{h}_{\ell}\hskip
1.0pt\mathbf{g}_{i}\hskip 1.0pts_{i}+n_{\ell},\vspace{-2pt}$ (6)

where $\beta_{\ell}$ denotes the link gain of UE $\ell$ composing of the
large-scale propagation effects of pathloss and shadow fading. Parameters
involved in $\beta_{\ell}$ are discussed in Sec. IV. Additionally,
$\mathbf{g}_{\ell}$ is the $M\times{}1$ un-normalized precoding vector from
the BS to UE $\ell$, obtained from the $\ell$-th column of $\mathbf{G}$, the
composite un-normalized $M\times{}L$ precoding matrix of all $L$ UEs. The data
symbol for UE $\ell$ is given by $s_{\ell}$, s.t.
$\mathbb{E}\\{|s_{\ell}|^{2}\\}=1$, $\forall\ell=\\{1,2,\dots,L\\}$. Finally,
the additive white Gaussian noise term at UE $\ell$ is denoted by $n_{\ell}$,
such that $n_{\ell}\sim\mathcal{CN}(0,\sigma^{2})$. For simplicity, we fix
$\sigma^{2}$ $\forall{}\ell=\\{1,2,\dots,L\\}$. The precoder normalization
parameter is given by $\eta=||\mathbf{G}||^{2}_{F}/L$, s.t.
$\mathbb{E}\\{||\mathbf{g}_{\ell}||^{2}\\}=1,\forall\ell=\\{1,2,\dots,L\\}$.

In the case of ZF precoding, $\mathbf{g}_{\ell}$ forms the $\ell$-th column of
$\mathbf{G}=\mathbf{H}^{H}(\mathbf{HH}^{H})^{-1}$, where
$\mathbf{H}=[\mathbf{h}_{1}^{T},\mathbf{h}_{2}^{T},\dots,\mathbf{h}_{L}^{T}]^{T}$
is the $L\times{}M$ matrix of all $L$ UE channels. Noting that
$\mathbf{HG}=\mathbf{HH}^{H}(\mathbf{HH}^{H})^{-1}=\mathbf{I}_{L}$, the ZF SNR
for UE $\ell$ is given by [4]

$\textrm{SNR}_{\ell}=\frac{P_{\textrm{EIRP}}\hskip
1.0pt\beta_{\ell}}{\sigma^{2}\hskip 1.0pt\eta}=\frac{P_{\textrm{EIRP}}\hskip
1.0pt\beta_{\ell}}{\sigma^{2}\left\\{\frac{1}{L}\left\\{\textrm{Tr}\left[\left(\mathbf{HH}^{H}\right)^{-1}\hskip
1.0pt\right]\hskip 1.0pt\right\\}\hskip 1.0pt\right\\}},$ (7)

with $\eta=\|\mathbf{G}\|^{2}_{F}/L=\textrm{Tr}\hskip 2.0pt[\hskip
1.0pt(\hskip 1.0pt\mathbf{HH}^{H})^{-1}\hskip 1.0pt]/L$. The expression in (7)
can be readily used to derive the ergodic sum spectral efficiency (in
bits/s/Hz) over all $L$ UEs. This is given by

$\textrm{R}=\mathbb{E}\left\\{\hskip 1.0pt\sum\limits_{\ell=1}^{L}\hskip
2.0pt\log_{2}\left(1+\textrm{SNR}_{\ell}\right)\right\\},\vspace{-2pt}$ (8)

with the expectation performed over small-scale fading in $\mathbf{H}$.

Remark 2. Analyzing the ergodic sum spectral efficiency performance with ZF
precoding in full-dimensional massive MIMO systems is a difficult task. This
is since finding _exact_ moments of the underlaying ZF SNR in (7) is
challenging due to the matrix trace in the denominator of (7) being a random
function of the inverse of the composite channel correlation matrix. Moreover,
the structure of $\mathbf{H}$ is fully heterogeneous, involving both azimuth
and elevation components. Below we devise a simple, yet accurate, analysis
framework of the ZF ergodic sum spectral efficiency via the expected ZF SNR.

##  III Analysis Framework and Insights

As exact analysis of ZF SNR moments is difficult, we _approximate_ the inverse
in (7) with a finite order Neumann series expansion [4, 5]. We first express
$\mathbf{HH}^{H}$ as

$\mathbf{HH}^{H}=\mathbb{E}\left\\{\mathbf{HH}^{H}\right\\}+\mathbf{\Xi}\hskip
1.0pt=\mathbf{\Delta}+\mathbf{\Xi},\vspace{-5pt}$ (9)

where $\mathbf{\Xi}=\mathbf{HH}^{H}-\mathbf{\Delta}$ and
$\mathbb{E}\left\\{\mathbf{\Xi}\right\\}$ may not be zero. This is in sharp
contrast to the analysis carried out in [4] for simple Rayleigh fading
channels, where $\mathbb{E}\\{\mathbf{\Delta}\\}=\mathbb{\mathbf{\Xi}}=0$.
With an order $N$ Neumann series, the inverse of (9) can be written by

$\displaystyle\left(\mathbf{HH}^{H}\right)^{-1}=$
$\displaystyle\left(\mathbf{\Delta}+\mathbf{\Xi}\hskip
1.0pt\right)^{-1}=\left[\mathbf{\Delta}\left(\mathbf{I}_{L}+\mathbf{\Delta}^{-1}\hskip
2.0pt\mathbf{\Xi}\hskip 1.0pt\right)\hskip 1.0pt\right]^{-1}$
$\displaystyle\approx$ $\displaystyle\sum\limits_{n=0}^{N}\hskip
1.0pt\left(\hskip 1.0pt-\mathbf{\Delta}^{-1}\hskip 2.0pt\mathbf{\Xi}\hskip
1.0pt\right)^{n}\mathbf{\Delta}^{-1}.$ (10)

Substituting the definition of $\mathbf{\Xi}$, we can write (10) as

$\left(\mathbf{HH}^{H}\right)^{-1}\approx\sum\limits_{n=0}^{N}\hskip
3.0pt\left[\hskip 1.0pt-\hskip
1.0pt\mathbf{\Delta}^{-1}\left(\mathbf{HH}^{H}-\mathbf{\Delta}\right)\hskip
1.0pt\right]^{n}\mathbf{\Delta}^{-1}.\vspace{-3pt}$ (11)

Analyzing the Neumann series in (11), when $N=2$ gives the expressions of (12)
(see top of the following page for space reasons)111Knowing the series order
apriori is a rather difficult task, which in general does not have a
prescribed methodology. For most common channel descriptions (see e.g.[5]),
expansion with $N=2$ has shown to converge to the true inverse of
$\mathbf{HH}^{H}$. Using this as our motivation, we also employ $N=2$ in this
paper to develop the analysis framework..

$\displaystyle\left(\mathbf{HH}^{H}\right)^{-1}\approx$ $\displaystyle\hskip
4.0pt\mathbf{\Delta}^{-1}-\left[\hskip
1.0pt\mathbf{\Delta}^{-1}\left(\mathbf{HH}^{H}-\mathbf{\Delta}\right)\mathbf{\Delta}^{-1}\right]+\left[\hskip
1.0pt\mathbf{\Delta}^{-1}\left(\mathbf{HH}^{H}-\mathbf{\Delta}\right)\mathbf{\Delta}^{-1}\left(\mathbf{HH}^{H}-\mathbf{\Delta}\right)\mathbf{\Delta}^{-1}\right]$
$\displaystyle=$ $\displaystyle\hskip 3.0pt3\hskip
1.0pt\mathbf{\Delta}^{-1}-3\hskip 1.0pt\mathbf{\Delta}^{-1}\hskip
1.0pt\mathbf{HH}^{H}\hskip
1.0pt\mathbf{\Delta}^{-1}+\mathbf{\Delta}^{-1}\hskip
1.0pt\mathbf{HH}^{H}\mathbf{\Delta}^{-1}\hskip
1.0pt\mathbf{HH}^{H}\mathbf{\Delta}^{-1}.$ (12)

Defining $\mathbf{Z}=\mathbf{\Delta}^{-1}\mathbf{HH}^{H}\mathbf{\Delta}^{-1}$
and substituting the resulting expression from (12) into (7) yields

$\textrm{SNR}_{\ell}\approx\frac{P_{\textrm{EIRP}}\hskip
2.0pt\beta_{\ell}}{\sigma^{2}\bigg{\\{}\frac{1}{L}\Big{\\{}\hskip
1.0pt\textrm{Tr}\hskip 1.0pt\big{[}\hskip 1.0pt3\mathbf{\Delta}^{-1}-3\hskip
1.0pt\mathbf{Z}+\mathbf{Z}\hskip
1.0pt\mathbf{HH}^{H}\mathbf{\Delta}^{-1}\hskip 1.0pt\big{]}\hskip
2.0pt\Big{\\}}\hskip 1.0pt\bigg{\\}}}\hskip 1.0pt.\vspace{-2pt}$ (13)

After further algebraic manipulations, the matrix trace in the denominator of
(13) can be broken down into the trace of the involved terms [6]. This is
given by

$\displaystyle\textrm{Tr}\left[\hskip 1.0pt\left(\hskip
1.0pt\mathbf{HH}^{H}\hskip 1.0pt\right)^{-1}\right]\approx$
$\displaystyle\hskip 3.0pt3\hskip 2.0pt\textrm{Tr}\left[\hskip
1.0pt\mathbf{\Delta}^{-1}\hskip 1.0pt\right]-3\hskip
2.0pt\textrm{Tr}\left[\hskip 1.0pt\mathbf{HH}^{H}\mathbf{\Delta}^{-2}\hskip
1.0pt\right]$ $\displaystyle+\textrm{Tr}\left[\hskip
1.0pt\mathbf{HH}^{H}\mathbf{\Delta}^{-1}\hskip
1.0pt\mathbf{HH}^{H}\mathbf{\Delta}^{-2}\hskip 1.0pt\right].$ (14)

Remark 3. In what follows, we evaluate the expected value of (13), using the
result in (14). To overcome the cumbersome nature of the expectation over
small-scale fading, we employ the univariate special case of the commonly used
Laplace approximation [4]. This allows us to write (15) (on top of the
following page due to space constraints).

$\mathbb{E}\left\\{\textrm{SNR}_{\ell}\right\\}\approx\frac{P_{\textrm{EIRP}}\hskip
2.0pt\beta_{\ell}}{\sigma^{2}\hskip 2.0pt\bigg{\\{}\hskip
2.0pt\frac{1}{L}\hskip 2.0pt\bigg{\\{}\mathbb{E}\left\\{3\hskip
2.0pt\textrm{Tr}\left[\hskip 1.0pt\mathbf{\Delta}^{-1}\right]-3\hskip
2.0pt\textrm{Tr}\left[\mathbf{HH}^{H}\mathbf{\Delta}^{-2}\right]+\textrm{Tr}\left[\mathbf{HH}^{H}\mathbf{\Delta}^{-1}\hskip
1.0pt\mathbf{HH}^{H}\mathbf{\Delta}^{-2}\right]\hskip 2.0pt\bigg{\\}}\hskip
3.0pt\bigg{\\}}\hskip 2.0pt\right\\}}.$ (15)

The approximation introduced in (15) is of the form
$\mathbb{E}\\{\delta/X\\}=\delta/\mathbb{E}\\{X\\}$, where $\delta$ is a
scalar value. The accuracy of such approximations rely on $X$ having a small
variance compared to its mean value - an effect which can be observed by
applying a Taylor series expansion to $\delta/X$ around
$\delta/\mathbb{E}\\{X\\}$. The terms in (15) are well suited to this
approximation, especially when $M$ and $L$ start to grow (common for a massive
MIMO system), since the implicit averaging in the denominator of (15) gives
rise to the variance reduction required for convergence (see e.g., [4]).

We further analyze (15) to derive an analytical solution with an order two
Neumann series. Taking the expectation through the trace in the denominator of
(15) allows us to write its three individual terms as

$T_{1}=3\hskip 2.0pt\textrm{Tr}\left[\hskip 1.0pt\mathbb{E}\hskip
1.0pt\\{\mathbf{\Delta}^{-1}\\}\hskip 1.0pt\right],\vspace{-3pt}$ (16)
$T_{2}=3\hskip 2.0pt\textrm{Tr}\left[\hskip
1.0pt\mathbb{E}\left\\{\mathbf{HH}^{H}\mathbf{\Delta}^{-2}\right\\}\hskip
1.0pt\right],\vspace{-2pt}$ (17)

and

$T_{3}=\textrm{Tr}\hskip 1.0pt\left[\hskip
1.0pt\mathbb{E}\left\\{\mathbf{HH}^{H}\mathbf{\Delta}^{-1}\hskip
2.0pt\mathbf{HH}^{H}\mathbf{\Delta}^{-2}\right\\}\hskip
1.0pt\right],\vspace{-2pt}$ (18)

respectively. As $\mathbb{E}\\{\mathbf{HH}^{H}\\}=\mathbf{\Delta}$, the first
two terms result in a cancellation, allowing us to express (15) as

$\displaystyle\mathbb{E}\left\\{\textrm{SNR}_{\ell}\right\\}\approx$
$\displaystyle\frac{P_{\textrm{EIRP}}\hskip
2.0pt\beta_{\ell}}{\sigma^{2}\hskip 1.0pt\bigg{\\{}\frac{1}{L}\hskip
1.0pt\bigg{\\{}\hskip 1.0pt\textrm{Tr}\left[\hskip
1.0pt\mathbb{E}\left\\{\mathbf{HH}^{H}\mathbf{\Delta}^{-1}\hskip
1.0pt\mathbf{HH}^{H}\right\\}\mathbf{\Delta}^{-2}\hskip
2.0pt\right]\bigg{\\}}\hskip 3.0pt\bigg{\\}}}\hskip 1.0pt.$ (19)

We now evaluate the remaining expectation in the denominator of (19).
Separating the expectation into its respective terms, for any
$\ell\neq{}j\in\left\\{1,2,\dots,L\right\\}$ we can write

$\displaystyle\mathbb{E}\left\\{\hskip 1.0pt\left(\hskip
1.0pt\mathbf{HH}^{H}\mathbf{\Delta}^{-1}\hskip 1.0pt\mathbf{HH}^{H}\hskip
1.0pt\right)_{\ell,j}\right\\}\overset{(a)}{=}$ $\displaystyle\hskip
2.0pt\mathbb{E}\left\\{\hskip
2.0pt\mathbf{h}_{\ell}\left(\mathbf{\Psi}\right)\mathbf{\Delta}^{-1}\left(\mathbf{\Psi}\right)^{H}\mathbf{h}_{j}^{H}\right\\}$
$\displaystyle=\mathbb{E}\left\\{\sum\limits_{r=1}^{L}\sum\limits_{s=1}^{L}\hskip
2.0pt\mathbf{h}_{\ell}\hskip
2.0pt\mathbf{h}_{r}^{H}\left(\mathbf{\Delta}^{-1}\right)_{r,s}\mathbf{h}_{s}\hskip
2.0pt\mathbf{h}_{j}^{H}\right\\}$ $\displaystyle\overset{(b)}{=}\hskip
2.0pt\sum\limits_{r=1}^{L}\sum\limits_{s=1}^{L}\mathbb{E}\left\\{\mathbf{h}_{\ell}\hskip
2.0pt\mathbf{h}_{r}^{H}\hskip 1.0pt\mathbf{h}_{s}\hskip
2.0pt\mathbf{h}_{j}^{H}\right\\}\left(\mathbf{\Delta}^{-1}\right)_{r,s}$
$\displaystyle=M^{2}\left(\mathbf{\Delta}^{-1}\right)_{\ell,j},$ (20)

where $(a)$ contains $\mathbf{\Psi}=(\mathbf{h}_{1}^{H}\hskip
2.0pt\mathbf{h}_{2}^{H}\hskip 2.0pt\dots\hskip 2.0pt\mathbf{h}_{L}^{H})$,
while $(b)$ contains straightforward manipulations to arrive at the result in
(20). Likewise, when $\ell=j\in\\{1,2,\dots,L\\}$, we can express

$\displaystyle\mathbb{E}\left\\{\left(\mathbf{HH}^{H}\mathbf{\Delta}^{-1}\hskip
1.0pt\mathbf{HH}^{H}\right)_{\ell,\ell}\right\\}=\mathbb{E}\left\\{\mathbf{h}_{\ell}\hskip
1.0pt\mathbf{h}_{\ell}^{H}\hskip 1.0pt\mathbf{h}_{\ell}\hskip
1.0pt\mathbf{h}_{\ell}^{H}\right\\}\left(\mathbf{\Delta}^{-1}\right)_{\ell,\ell}$
$\displaystyle+\sum\limits_{\begin{subarray}{c}r=1\\\
r\neq{}\ell\end{subarray}}^{L}\mathbb{E}\left\\{\mathbf{h}_{\ell}\hskip
2.0pt\mathbf{h}_{r}^{H}\mathbf{h}_{r}\hskip
1.0pt\mathbf{h}_{\ell}^{H}\right\\}\left(\mathbf{\Delta}^{-1}\right)_{r,r}$
$\displaystyle=\left(\mathbf{\Delta}^{-1}\right)_{\ell,\ell}\hskip
1.0pt\mathbb{E}\left\\{\mathbf{h}_{\ell}\hskip
2.0pt\mathbf{h}_{\ell}^{H}\hskip 2.0pt\mathbf{h}_{\ell}\hskip
2.0pt\mathbf{h}_{\ell}^{H}\right\\}$
$\displaystyle+\sum\limits_{\begin{subarray}{c}r=1\\\
r\neq{}\ell\end{subarray}}^{L}\left(\mathbf{\Delta}^{-1}\right)_{r,r}\textrm{Tr}\left[\hskip
1.0pt\mathbf{R}_{r}\mathbf{R}_{\ell}\hskip 1.0pt\right],$ (21)

where $\mathbf{R}_{\ell}=\mathbb{E}\\{\mathbf{h}_{\ell}\hskip
1.0pt\mathbf{h}_{\ell}^{H}\\}$. Recognizing that the expectation in (21) can
be evaluated by stating $\mathbb{E}\left\\{\mathbf{h}_{\ell}\hskip
1.0pt\mathbf{h}_{\ell}^{H}\hskip 1.0pt\mathbf{h}_{\ell}\hskip
1.0pt\mathbf{h}_{\ell}^{H}\right\\}=\mathbb{E}\left\\{|\hskip
1.0pt\mathbf{h}_{\ell}\hskip 1.0pt\mathbf{h}_{\ell}^{H}|^{2}\right\\}$,
following lengthy but straightforward calculations result in
$\mathbb{E}\left\\{|\mathbf{h}_{\ell}\hskip
1.0pt\mathbf{h}_{\ell}^{H}|^{2}\right\\}=M^{2}+\textrm{Tr}\hskip
2.0pt[\mathbf{R}_{\ell}^{2}]$. We only resort to only a brief presentation due
to space reasons. As such,

$\displaystyle\mathbb{E}\left\\{\hskip 2.0pt\left(\hskip
1.0pt\mathbf{HH}^{H}\mathbf{\Delta}^{-1}\mathbf{HH}^{H}\right)_{\ell,\ell}\right\\}=$
$\displaystyle\left(\mathbf{\Delta}^{-1}\right)_{\ell,\ell}\Big{\\{}M^{2}+\textrm{Tr}\left[\mathbf{R}_{\ell}^{2}\right]\Big{\\}}\hskip
2.0pt$ $\displaystyle+\sum\limits_{\begin{subarray}{c}r=1\\\
r\neq{}\ell\end{subarray}}^{L}\left(\mathbf{\Delta}^{-1}\right)_{r,r}\Big{\\{}\textrm{Tr}\left[\mathbf{R}_{r}\hskip
1.0pt\mathbf{R}_{\ell}\hskip 1.0pt\right]\Big{\\}}.$ (22)

Combining the results in (22) and (20), the overall expectation

$\displaystyle\mathbb{E}\left\\{\mathbf{HH}^{H}\mathbf{\Delta}^{-1}\mathbf{HH}^{H}\right\\}$
$\displaystyle=M^{2}\mathbf{\Delta}^{-1}+\textrm{diag}\left(d_{1},d_{2},\dots,d_{L}\right)$
$\displaystyle=M^{2}\mathbf{\Delta}^{-1}+\mathbf{D},$ (23)

where $d_{\ell}=\sum\nolimits_{r=1}^{L}(\mathbf{\Delta}^{-1})_{\hskip
2.0ptr,r}\hskip 2.0pt\textrm{Tr}\hskip
2.0pt[\mathbf{R}_{r}\mathbf{R}_{\ell}]$. Substituting (23) into (19), we
obtain
$\mathbb{E}\\{\textrm{SNR}_{\ell}\\}=\mathbb{E}\\{\textrm{SNR}_{\ell,\textrm{approx}}\\}$
as

$\displaystyle\mathbb{E}\left\\{\textrm{SNR}_{\ell,\textrm{approx}}\right\\}\approx$
$\displaystyle\frac{P_{\textrm{EIRP}}\hskip
2.0pt\beta_{\ell}}{\sigma^{2}\hskip 1.0pt\bigg{\\{}\frac{1}{L}\hskip
1.0pt\bigg{\\{}\hskip 1.0pt\textrm{Tr}\left[\hskip
2.0pt\left(M^{2}\mathbf{\Delta}^{-1}+\mathbf{D}\right)\mathbf{\Delta}^{-2}\hskip
1.0pt\right]\bigg{\\}}\bigg{\\}}}.$ (24)

This can be used to directly approximate the ergodic sum spectral efficiency
by

$\textrm{R}_{\textrm{approx}}\approx{}L\log_{2}\left(1+\mathbb{E}\left\\{\textrm{SNR}_{\ell,\textrm{approx}}\right\\}\right).\vspace{-3pt}$
(25)

While the tightness of (24) and (25) is assessed and discussed in Sec. IV
relative to their simulated counterparts, below we give some remarks to draw
insights on the derived results.

Remark 4. The results in (24) and (25) present a simple solution to a
difficult problem, lending itself to some useful insights. The denominator of
(24) (expected noise power of UE $\ell$) is primarily governed by
$\textrm{Tr}[\mathbf{\Delta}]=\textrm{Tr}[\hskip
1.0pt\mathbb{E}\\{\mathbf{HH}^{H}\\}\hskip 1.0pt]$ and
$\textrm{Tr}[\mathbf{R}_{\ell}\mathbf{R}_{r}]$. Maximization of the both trace
terms would lead to maximum expected noise power. This happens when
$\mathbf{R}_{\ell}$ starts to become similar to $\mathbf{R}_{r}$, i.e., the
_maximum eigenvectors_ of $\mathbf{R}_{\ell}$ and $\mathbf{R}_{r}$ are
aligned/in-phase. In general, such an alignment is caused by similarity of
channels to UEs $\ell$ and $r$ where $r\neq{}\ell$. This is in turn controlled
by the propagation channel’s spatial parameters such as the LOS angles and
angular spreads in both azimuth and elevation domains. In general, non-
overlapping angular support sets would result in independent channels, meaning
that the UE and scatterer locations do not correlate the individual channel
responses. Typically, when the UEs are closely spaced, their impulse responses
are highly correlated, where the $\textrm{Tr}[\hskip
1.0pt\mathbb{E}\\{\mathbf{HH}^{H}\\}\hskip 1.0pt]$ and
$\textrm{Tr}[\mathbf{R}_{\ell}\hskip 1.0pt\mathbf{R}_{r}]$ approach their
maximum value. Later we demonstrate this effect qualitatively by experimenting
with the relative UE angular separations in both azimuth and elevation
domains. Furthermore, another prominent factor influencing the maximization of
expected noise power is the ratio of $M^{2}/L$, implying that for a given
number of BS antennas, increasing the number of UEs will _linearly_ increase
the expected noise power and thus decrease the overall ZF SNR. Such insights
are difficult to obtain from more complex solutions derived in the literature
(see e.g., [7]), for simpler propagation channels which require a linked set
of equations, even in the large system regime. In sharp contrast, our analysis
poses no such constraints and is general to the operational system dimension,
channel model, channel correlation structure and antenna array geometries.

##  IV Numerical Results

We consider an urban non-LOS macrocellular environment from the propagation
channel model in [8]. Operating at 3.7 GHz ($\lambda_{f}$=8.1 cm), we have
$R_{c}$ = 500 m with $\sigma^{2}=1;\hskip
3.0pt\forall{}\ell=\\{1,2,\dots,L\\}$. As such we refer to the term “operating
SNR” to denote $P_{\textrm{EIRP}}/\sigma^{2}$. We employ inter-element spacing
of 0.5$\lambda_{f}$ in the azimuth and 0.7$\lambda_{f}$ in elevation. The link
gain for UE $\ell$, $\beta_{\ell}$, follows the model quoted in Table 7.4.1-1
of [8] with pathloss exponent of 3.67 and zero-mean log-normally distributed
shadowing having a standard deviation of 4 dB. The total number of MPCs,
$N_{P}$, and other parameters such as the ASDs in azimuth and elevation are
employed from Table 7.5-6 Part I in [8]. Unless otherwise stated, we assume
that the BS is equipped with UPA having $M=128$ antennas (16 rows, 8 columns)
serving $L=4$ UEs. This is since the UPA better captures the full-dimensional
channel, relative to the ULA. Though this is our default system configuration,
we confirm that the overall trends do not change with varying values of $M$ or
$L$. Inter-UE angular separations in the azimuth and elevation are defined in
relation to a given UE’s azimuth and elevation LOS angles. Relative to these,
a constant offset in the LOS angles determines the angular separation of other
UEs. This allows us to study the impact of closely spaced UEs relative to UEs
which may be far apart. We first evaluate the accuracy of the proposed Neumann
series approximation, followed by expected SNR and ergodic sum spectral
efficiency estimates for a range of UE placements.

Figure 1 depicts the Neumann series approximation accuracy assessed by
comparing the true $(\mathbf{HH}^{H})^{-1}$ relative to that derived in (12).
In order to quantify the total approximation accuracy, with $N=2$, we define
the following metric: $\textrm{error magnitude},\alpha=\textrm{Tr}[\hskip
1.0pt\mathbf{I}_{L}-[(\mathbf{HH}^{H})^{-1}/\hskip
1.0pt\textrm{approx.NS}\hskip 1.0pt]\hskip 1.0pt]$, where “approx.NS”
represents (12). The figure shows the cumulative distribution functions (CDFs)
of the error magnitude with variation in the number of BS antennas, $M$ with
$L=4$ at operating SNR = 0 dB. Two interesting trends can be observed:
Firstly, increasing $M$ significantly reduces the error magnitude (comparing
when $M$=10 and 32 to $M$=64). This difference is due to the improvement
brought by scaling up the degrees-of-freedom induced by higher values of $M$
relative to $L$. At $M=64$ the error magnitude CDF seems to have largely
converged to zero, while the top end of the CDF still exhibits very small
error magnitudes. Secondly, errors at the top end of the CDF start to
disappear as $M$ is increased from 64 to 128 and 256, as seen by the zoomed
part of the figure, where the inverse becomes well conditioned. Though not
shown here due to space reasons, we confirm that if the _order_ of the Neumann
series is reduced to $N=1$, the approximation suffers from 40% loss in
accuracy around CDF = 0.5.

Figure 1: CDFs of Neumann series approximation error magnitude as a function
of BS antenna numbers when $N=2$, $L=4$ and operating SNR = 0 dB.

Figure 2 shows the expected per-UE ZF SNR as a function of operating SNR with
$M$=128, $L$=4 for various inter-UE separations in both azimuth and elevation
domains. Some important conclusions can be drawn from the figure. Firstly, in
the top sub-figure, the impact of UE separation in the azimuth domain is
examined. For consistency in comparison, elevation UE separation of 15∘ and
azimuth UE separation of 30∘ denotes uncorrelated channels. As the UE
separation is reduced, one can notice a significant reduction (up to 6 dB) in
the expected ZF SNR. This is seen by comparing azimuth UE separations of 30∘
relative to 7.5∘ or 5∘, while keeping elevation separation at 15∘. As UEs
become more closely spaced, the expected ZF SNR suffers from noise
enhancements due to similar channels and maximization of the trace terms
referred to in Remark 4. Even separations of 10∘ are insufficient to approach
the performance of uncorrelated channels. The separations needed for achieving
uncorrelated channels are unlikely to be seen in reality, due to the width of
the required angular separations. Naturally, this is sensitive to the type of
multipath environment. Secondly, we can observe that our derived expression of
(24) approximates the simulated expected ZF SNR extremely well, as it remains
tight for both correlated and uncorrelated scenarios and across all operating
SNRs. The bottom sub-figure depicts the impact of reduced elevation
separations, where one can notice a smaller difference in the expected ZF SNR,
as the elevation UE separation is reduced from 15∘ to 10∘. To this end, UE
separation reduction in the elevation does not have as prominent impact as in
the azimuth. This is a subject of further study and investigations.

Figure 2: Expected per-UE ZF SNR vs. operating SNRs with $M$=128, $L$=4, and
varying azimuth and elevation UE separations. Figure 3: Ergodic sum spectral
efficiency vs. operating SNR with $M$=128, $L$=4 and different azimuth UE
separations. Elevation separations are 15∘.

Figure 3 shows the ergodic sum spectral efficiency performance as a function
of the operating SNRs. The exact same simulation scenarios and parameters are
employed as for the top sub-figure of Fig. 2. We can observe similar trends to
those for top sub-figure of Fig. 2, where the ergodic sum spectral efficiency
exhibits an almost linear growth with increasing operating SNR levels. Here
also for each scenario, our derived expression in (25) closely matches the
simulated performance.

##  V Conclusions

We articulate a simple, tight and widely applicable analytical framework for
understanding the performance of ZF precoding in full-dimensional massive MIMO
systems. Via an order two Neumann series, we approximate the expected per-UE
ZF SNR and ergodic sum spectral efficiency. Our framework caters for finite
multipath propagation channels, correlated and uncorrelated UE positions, and
multiple BS array configurations. Our analysis reveals several important
insights and describes the interaction of propagation and system parameters on
ZF performance. To achieve optimal performance, we show that UEs could be
separated in the azimuth by 30∘ in azimuth and $15^{\circ}$ in elevation, yet
this depends on the underlying multipath conditions. Nevertheless, this will
be difficult to achieve in practice.

## References

  * [1] F. Rusek,  _et al.,_ “Scaling up MIMO: Opportunities and challenges with very large arrays,” _IEEE Signal Process. Mag._ , vol. 30, no. 1, pp. 40-60, Jan. 2013.
  * [2] G. Xu,  _et al.,_ “Full dimension MIMO (FD-MIMO): Demonstrating commercial feasibility,” _IEEE J. Sel. Areas in Commun._ , vol. 35, no. 8, pp. 1876-1886, Aug. 2017.
  * [3] T. L. Marzetta,  _et al._ , _Fundamentals of Massive MIMO_ , Cambridge University Press, 2016.
  * [4] H. Tataria,  _et al.,_ “On the general analysis of coordinated regularized zero-forcing precoding: An application to two-tier small-cell networks,” _IEEE Trans. Commun._ , vol. 65, no. 7, pp. 3133-3150, Jul. 2017.
  * [5] H. Wang,  _et al._ , “An efficient detector for massive MIMO based on improved matrix partition,” _IEEE Trans. Signal Process._ , Mar. 2021.
  * [6] R. A. Horn and C. R. Johnson, _Matrix Analysis_ , Cambridge University Press, 1985.
  * [7] A. Adhikary,  _et al._ , “Joint spatial division and multiplexing for mm-Wave channels,” _IEEE J. Sel. Areas Commun._ , vol. 32, no. 6, pp. 1239-1255, Jun. 2014.
  * [8] 3GPP TR 38.901, “Study on channel model for frequencies from 0.5 to 100 GHz (v16.1.0)”, 3GPP Release 14, Nov. 2020.