Energy-Efficient Power Control in Cell-Free and User-Centric Massive MIMO at Millimeter Wave

TL;DR

This paper explores energy-efficient power control in cell-free and user-centric massive MIMO systems at millimeter wave frequencies, introducing novel channel models, low-complexity beamforming, and power allocation strategies that enhance performance without requiring channel estimation at mobile stations.

Contribution

It introduces a multiuser clustered millimeter wave channel model and proposes channel-independent beamforming and power control schemes for energy efficiency in cell-free and user-centric massive MIMO.

Findings

Power control significantly improves energy efficiency.

Hybrid beamforming incurs manageable performance loss.

Proposed schemes do not require channel estimation at mobile stations.

Abstract

In a cell-free massive MIMO architecture a very large number of distributed access points simultaneously and jointly serves a much smaller number of mobile stations; a variant of the cell-free technique is the user-centric approach, wherein each access point just serves a reduced set of mobile stations. This paper introduces and analyzes the cell-free and user-centric architectures at millimeter wave frequencies, considering a training-based channel estimation phase, and the downlink and uplink data transmission phases. First of all, a multiuser clustered millimeter wave channel model is introduced in order to account for the correlation among the channels of nearby users; second, an uplink multiuser channel estimation scheme is described along with low-complexity hybrid analog/digital beamforming architectures. Third, the non-convex problem of power allocation for downlink global energy efficiency maximization is addressed. Interestingly, in the proposed schemes no channel estimation is needed at the mobile stations, and the beamforming schemes used at the mobile stations are channel-independent and have a very simple structure. Numerical results show the benefits granted by the power control procedure, that the considered architectures are effective, and permit assessing the loss incurred by the use of the hybrid beamformers and by the channel estimation errors.

Tables

Table 1. Table I: List of acronyms

5G	Fifth-Generation
AP	Access-Point
BCD-SD	Block Coordinate Descent for Subspace Decomposition
CF	Cell-Free
CoMP	Coordinated Multipoint
CPU	Central Processing Unit
CSI	Channel State Information
FD	Fully-Digital
HY	Hybrid
LMMSE	Linear Minimum Mean Square Error
LTE	Long-Term Evolution
MIMO	Multiple-Input-Multiple-Output
mmWave	Millimeter Wave
MS	Mobile Station
TDD	Time Division Duplex
ULA	Uniform Linear Array
ZF	Zero-Forcing

Table 2. Table II: Parameters for Pathloss model

Scenario	Model Parameters
UMi Street Canyon LOS	n=1.98, $\sigma$=3.1 dB
UMi Street Canyon NLOS	n=3.19, $\sigma$=8.2 dB
UMi Open Square LOS	n=2.89, $\sigma$=7.1 dB
UMi Open Square NLOS	n=1.73, $\sigma$=3.02 dB

Equations

$${\text{Y}}_m=\sum _{k=1}^K\sqrt{p_k}{\text{S}}_{k,m}{\mathbf{\Phi }}_k+{\text{W}}_m,$$

$${\text{y}}_m=\sum _{k=1}^K\sqrt{p_k}{\text{A}}_k{\text{s}}_{k,m}+{\text{w}}_m$$

$$\begin{split}&\mathbb{E}[\parallel\textbf{V}_{k,m}^{H}\textbf{y}_{m}-\textbf{s}_{k,m}\parallel^{2}]=\text{tr}(\textbf{V}_{k,m}^{H}\mathbb{E}[\textbf{y}_{m}\textbf{y}_{m}^{H}]\textbf{V}_{k,m})+\\ &+\mathbb{E}[\parallel\textbf{s}_{k,m}\parallel^{2}]-\mathbb{E}[2\Re\{\text{tr}(\textbf{s}_{k,m}^{H}\textbf{V}_{k,m}^{H}\textbf{y}_{m})\}]=\\ &=\text{tr}\bigg{(}\textbf{V}_{k,m}^{H}\bigg{(}\sum\limits_{l=1}^{K}p_{k}\textbf{A}_{l}\textbf{A}_{l}^{H}+\sigma^{2}\textbf{I}\bigg{)}\textbf{V}_{k,m}\bigg{)}+\\ &+\mathbb{E}[\parallel\textbf{s}_{k,m}\parallel^{2}]-\sqrt{p_{k}}\text{tr}(\textbf{V}_{k,m}^{H}\textbf{A}_{k}+\textbf{V}_{k,m}^{T}\textbf{A}_{k}^{*})\end{split}$$

$$\textbf{V}_{k,m}^{LMMSE}=\sqrt{p_{k}}\big{(}\sum\limits_{l=1}^{K}p_{l}\textbf{A}_{l}\textbf{A}_{l}^{H}+\sigma^{2}\textbf{I}\big{)}^{-1}\textbf{A}_{k}\;.$$

$${\text{Q}}_{k,m}=(\hat{\mathcal{G}}{\hat{\mathcal{G}}}^H{)}^{-1}{\hat{\text{S}}}_{k,m}.$$

$${\text{Q}}_{k,m}=\frac{{\text{Q}}_{k,m}}{\sqrt{\text{mbox}tr({\text{Q}}_{k,m}{\text{Q}}_{k,m}^H)}},$$

$${\mathbf{Q}}_m^{opt}=\left[{\mathbf{Q}}_{1,m},\dots ,{\mathbf{Q}}_{K,m}\right],$$

$${\mathbf{Q}}_m^{BB}=\left[{\mathbf{Q}}_{1,m}^{BB},\dots ,{\mathbf{Q}}_{K,m}^{BB}\right],$$

$${\text{s}}_m^{CF}(n)=\sum _{k=1}^K\sqrt{{\eta }_{m,k}}{\text{Q}}_{k,m}{\text{x}}_k^{DL}(n),$$

$$\begin{array}{rl}{\text{y}}_k^{CF}(n)& =\sum _{m=1}^M{\text{H}}_{k,m}^H{\text{s}}_m^{CF}(n)+{\text{z}}_k(n)=\\ & =\sum _{m=1}^M\sqrt{{\eta }_{m,k}}{\text{H}}_{k,m}^H{\text{Q}}_{k,m}{\text{x}}_k^{DL}(n)+\\ & +\sum _{l=1,l\ne k}^K\sum _{m=1}^M\sqrt{{\eta }_{m,l}}{\text{H}}_{k,m}^H{\text{Q}}_{l,m}{\text{x}}_l^{DL}(n)+{\text{z}}_k(n),\end{array}$$

$${\hat{\text{x}}}_k^{DL,CF}(n)={\text{L}}_k^H{\text{y}}_k^{CF}(n).$$

$$\mathcal{M}(k)=\{m:k\in \mathcal{K}(m)\}.$$

$${\text{s}}_m^{UC}(n)=\sum _{k\in \mathcal{K}(m)}\sqrt{{\eta }_{m,k}}{\text{Q}}_{k,m}{\text{x}}_k^{DL}(n).$$

$$\begin{array}{rl}{\text{y}}_k^{UC}(n)& =\sum _{m=1}^M{\text{H}}_{k,m}^H{\text{s}}_m^{UC}(n)+{\text{z}}_k(n)=\\ & =\sum _{m\in \mathcal{M}(k)}\sqrt{{\eta }_{m,k}}{\text{H}}_{k,m}^H{\text{Q}}_{k,m}{\text{x}}_k^{DL}(n)+\\ & +\sum _{l=1,l\ne k}^K\sum _{m\in \mathcal{M}(l)}\sqrt{{\eta }_{m,l}}{\text{H}}_{k,m}^H{\text{Q}}_{l,m}{\text{x}}_l^{DL}(n)+{\text{z}}_k(n),\end{array}$$

$${\hat{\text{x}}}_k^{DL,UC}(n)={\text{L}}_k^H{\text{y}}_k^{UC}(n).$$

$${\mathcal{R}}_k=\text{B}{\text{log}}_2\mid \text{I}+{\text{R}}_k^{-1}{\text{A}}_{k,k}{\text{A}}_{k,k}^H\mid ,$$

$${\text{A}}_{k,l}=\sum _{m\in \mathcal{M}(k)}\sqrt{{\eta }_{m,l}}{\text{L}}_k^H{\text{H}}_{k,m}^H{\text{Q}}_{l,m}.$$

$${\text{R}}_k=\sum _{l\ne k}{\text{A}}_{k,l}{\text{A}}_{k,l}^H+{\sigma }_z^2{\text{L}}_k^H{\text{L}}_k,$$

$${\text{B}}_{k,l,m}={\text{L}}_k^H{\text{H}}_{k,m}^H{\text{Q}}_{l,m}⟹{\text{A}}_{k,l}=\sum _{m\in \mathcal{M}(k)}\sqrt{{\eta }_{m,l}}{\text{B}}_{k,l,m}$$

$$\begin{array}{rl}{\text{R}}_k& =\sum _{l\ne k}\sum _{m\in \mathcal{M}(l)}\sum _{m^{\prime }\in \mathcal{M}(l)}\sqrt{{\eta }_{m,l}{\eta }_{m^{\prime },l}}{\text{B}}_{k,l,m}{\text{B}}_{k,l,m^{\prime }}^H+\\ & +{\sigma }_z^2{\text{L}}_k^H{\text{L}}_k,\end{array}$$

$${\mathcal{R}}_k=\text{B}{\text{log}}_2\left|\begin{array}{c}\text{I}+{\text{R}}_k^{-1}\sum _{m,m^{\prime }}\sqrt{{\eta }_{m,k}{\eta }_{m^{\prime },k}}{\text{B}}_{k,k,m}{\text{B}}_{k,k,m^{\prime }}^H\end{array}\right|$$

$$\begin{array}{rl}{\mathcal{R}}_k& =\underset{g_1(\mathbf{\eta })}{\underbrace{\text{B}{\text{log}}_2\left|\begin{array}{c}{\sigma }_z^2{\text{L}}_k^H{\text{L}}_k+\sum _l^K\sum _m\sum _{m^{\prime }}\sqrt{{\eta }_{m,l}{\eta }_{m^{\prime },l}}{\text{B}}_{k,l,m}{\text{B}}_{k,l,m^{\prime }}^H\end{array}\right|}}\\ & -\underset{g_2(\mathbf{\eta })}{\underbrace{\text{B}{\text{log}}_2\left|\begin{array}{c}{\sigma }_z^2{\text{L}}_k^H{\text{L}}_k+\sum _{l\ne k}^K\sum _m\sum _{m^{\prime }}\sqrt{{\eta }_{m,l}{\eta }_{m^{\prime },l}}{\text{B}}_{k,l,m}{\text{B}}_{k,l,m^{\prime }}^H\end{array}\right|}},\end{array}$$

$${\text{y}}_m(n)=\sum _{k=1}^K\sqrt{{\tilde{\eta }}_k}{\text{H}}_{k,m}{\text{L}}_k{\text{x}}_k^{UL}(n)+{\text{w}}_m(n)$$

$${\tilde{\text{y}}}_{m,k}(n)={\text{Q}}_{k,m}^H{\text{y}}_m(n),$$

$${\hat{\text{x}}}_k^{UL}(n)=\sum _{m=1}^M{\tilde{\text{y}}}_{m,k}(n),k=1,\dots ,K.$$

$${\hat{\text{x}}}_k^{UL,UC}=\sum _{m\in \mathcal{M}(k)}{\tilde{\text{y}}}_{m,k}(n),k=1,\dots ,K,$$

$$\begin{array}{rl}{\hat{\text{x}}}_k^{UL,UC}& =\sum _{m\in \mathcal{M}(k)}{\text{Q}}_{k,m}^H\sqrt{{\tilde{\eta }}_k}{\text{H}}_{k,m}{\text{L}}_k{\text{x}}_k^{UL}(n)+\\ & +\sum _{j=1,j\ne k}^K\sqrt{{\tilde{\eta }}_j}\sum _{m\in \mathcal{M}(k)}{\text{Q}}_{k,m}^H{\text{H}}_{j,m}{\text{L}}_j{\text{x}}_j^{UL}(n)+\\ & +\sum _{m\in \mathcal{M}(k)}{\text{Q}}_{k,m}^H{\text{w}}_m(n).\end{array}$$

$${\tilde{\mathcal{R}}}_k=B{\log }_2\left|\mathbf{I}+{\tilde{\eta }}_k{\mathbf{R}}_k^{-1}\left(\sum _{m\in \mathcal{M}(k)}{\mathbf{B}}_{k,k,m}^H\right)\left(\sum _{m\in \mathcal{M}(k)}{\mathbf{B}}_{k,k,m}\right)\right|,$$

$$\begin{array}[]{lll}\mathbf{R}_{k}=&\displaystyle\sum_{j\neq k}\widetilde{\eta}_{j}\!\!\left(\sum_{m\in{\cal M}(j)}\mathbf{B}_{j,k,m}^{H}\right)\!\!\!\left(\sum_{m\in{\cal M}(j)}\mathbf{B}_{j,k,m}\right)\\ \\ &+\sigma^{2}_{w}\mbox{card}({\cal M}(k))\;,\end{array}$$

$$\underset{\mathbf{\eta }}{\max }\frac{\sum _{k=1}^K{\mathcal{R}}_k(\mathbf{\eta })}{\sum _{m=1}^M\left[\sum _{k\in \mathcal{K}(m)}\delta {\eta }_{m,k}+P_{c,m}\right]}$$

Full text

I Introduction

Future 5G wireless systems will heavily rely on the use of large-scale antenna arrays, a.k.a. massive MIMO, and of carrier frequencies above $10\,\text{GHz}$, the so called mmWave frequencies [1]. Indeed, on one hand, the use of massive MIMO permits serving several users on the same time-frequency resource, while, on the other one, mmWave carrier frequencies will enable the use of much larger bandwidths, at least on short-distances (up to about one hundred meters). The combined use of massive MIMO systems along with mmWave frequencies is indeed one of the key technological enablers of the envisioned wireless Gbit/s experience [2]. For conventional sub-$6\,\text{GHz}$ frequencies, a new communications architecture, named CF massive MIMO, has been recently introduced in [3, 4], in order to alleviate the cell-edge problem and thus increase the system performance of unlucky users that happen to be located very far from their serving AP. In the CF architecture, instead of having few base stations with massive antenna arrays, a very large number of simple APs randomly and densely deployed serve a much smaller number of MSs. This approach has some similarities with the CoMP [5, 6] and network MIMO [7]; these technologies exploit coordinate beamforming, scheduling, and joint transmission using multiple distributed antennas and/or base stations, with the aim of mitigating the interference and achieving diversity gains. There are however some key differences between these strategies and the CF massive MIMO approach discussed in this paper. In particular, we have that: (a) the CF massive MIMO strategy fully leverages the advantages and unique features of massive MIMO, such as the channel hardening effect and the use of the TDD protocol;and (b) the CF massive MIMO also makes a more limited use of the backhaul link since the channel coefficients are locally estimated at the APs and are mutually shared, and, moreover, beamformers are computed at the APs using locally available information. The CF massive MIMO architecture can be thus thought as the scalable way of implementing distributed network MIMO deployments. In the CF architecture described in [3, 4], single-antenna APs and MSs are considered, all the APs serve all the MSs, all the APs are connected through a backhaul link to a CPU, but every AP performs channel estimation locally, and channel estimates are not sent to the CPU, but are locally exploited. In particular, for the downlink communication phase, the CPU sends to the APs the data symbols to be sent to all the MSs, while for the uplink communication phase, the APs use the backhaul to send the sufficient statistics for all the MSs to the CPU, which combines them and performs uplink data decoding. In [8] a UC variant of the CF approach is introduced, wherein each APs, instead of serving all the MSs in the considered area, just serves the ones that he receives best; the results of [8] show that the UC approach provides savings on the required backhaul capacity and, also, provides better data-rates to the vast majority of the users.

The paper [9] analyzes a CF massive MIMO system from the point of view of its energy consumption, taking into account the power consumed by the backhaul links, the number of active APs and the number of antennas at each AP. The paper also embraces the UC philosophy by considering AP-MS selection schemes aimed at maximizing the system energy efficiency. While previous papers deal with the case in which the APs and the MSs are both equipped with one antenna (with the exception of [9] which considers multiple antennas at the APs), in [10] the CF and UC architectures are generalized to the case in which both the APs and the MSs are equipped with multiple antennas. This generalization is not trivial and indeed the system proposed in [10], despite the use of multiple antennas at the MSs and the use of a multiplexing order larger than one, develops a beamforming scheme at the APs that does not require channel estimation at the MSs, which adopt a channel-independent beamformer. While [10] considers a simple uniform power allocation, in [11, 12] power control procedures aimed at rate maximization are proposed. The paper [13] investigates the use of pilot signals on the downlink, so as to enable channel estimation at the MSs; even though such an approach provides some performance improvement, it contradicts the TDD-based "Massive MIMO philosophy" wherein no channel estimation is required at the MSs. A compute-and-forward approach to CF massive MIMO was then proposed in [14] in order to reduce the load on the backhaul links.

All of the above cited papers consider the case in which a conventional sub-6 GHz carrier frequency is used. On the other hand, as already stated, mmWave frequencies will play a big role in future wireless cellular systems due to availability of wide amounts of unused bandwidths [15]. The design of a wireless cellular system operating at mmWave poses different and new challenges with respect to conventional sub-6 GHz frequencies, due to different propagation mechanisms [15], to the need of using large antenna arrays at both side of the links to counteract the increased path-loss [16], and to the presence of hardware constraints that prevent the realization of FD precoding and postcoding beamforming structures [17]. Despite these tough challenges, the attractiveness of the mmWave frequencies for cellular communications has led to intense research efforts in the last few years. Due to the difficulty of realizing FD beamforming architectures at mmWave with antenna arrays of large size, the CF massive MIMO architecture at these high frequencies appears particularly attractive, mainly for two reasons; first of all, it substitutes large co-located antennas with several APs with antenna arrays of smaller dimension, and, thus, with less hardware complexity; second, due to the limited range of mmWave communications, the distributed AP dense deployment alleviates the cell-edge problem and creates path-diversity, which is precious since signal blockages may frequently happen at these frequencies.

Energy efficiency is another key topic that should be considered when designing modern wireless communication systems. Indeed, increasing the bit-per-Joule energy efficiency is regarded as a key requirement of future 5G networks [18]. A recent survey on the most promising energy-efficient techniques for 5G has recently appeared in [19]; in that paper, radio resource allocation is identified, among other techniques, as a relevant approach to improve the energy efficiency of future wireless networks. Recent contributions on energy-efficient radio resource allocation for 5G are [20, 21, 22, 23].

Motivated by the above discussion, this paper considers details a CF and UC massive MIMO wireless system operating at mmWave and focuses on resource allocation strategies maximizing the global energy efficiency. To the best of authors’ knowledge, this is the first paper that considers the CF and UC massive-MIMO deployments at mmWave; preliminary results on this issue have been reported by the authors in the conference papers [24] and [25], wherein it has also been shown that the UC approach generally outperforms the CF approach.

The contribution of the paper can be summarized as follows. First of all, we introduce a multiuser mmWave channel model that permits taking into account channel correlation for close users. Building upon the clustered channel model [26] widely used at mmWave frequencies, we extend this model to take into account the fact that if several APs and MSs are in the same area, their channels must be built using the same set of scatterers; adopting this model, users that are very close will receive beams with very close direction of arrival and so channel correlation for nearby users is intrinsically taken into account. Then, we study the UC and CF approaches at mmWave frequencies; we assume that both the APs and MSs are equipped with multiple antennas, use HY analog-digital partial ZF beamforming at the APs, while a very simple 0-1 beamforming architecture, independent of the channel estimate, is used at the MSs. We present simulation results for scenarios that can be representative of a lightly-loaded system, with $M=80$ APs and $K=6$ MSs transmitting on the same time-frequency slot, and of a highly-loaded system, with $M=80$ APs and $K=16$ MSs. Results show that over a bandwidth of $200\,\textrm{MHz}$ and with a maximum transmit power of $0.1\,\textrm{W}$ at each AP, taking into account channel estimation errors and using low-complexity beamforming structures, the downlink average rate per user is $1.5\,\textrm{Gbit/s}$ in the lightly-loaded scenario and $400\,\textrm{Mbit/s}$ in the heavily-loaded situation. Likewise, in the lightly-loaded scenario, the average uplink rate-per-user is about $1\,\textrm{Gbit/s}$.

This paper is organized as follows. Next section is devoted to the discussion of the system and to the description of the processing (it involves uplink training, downlink data transmission and uplink data transmission); Section III concerns global energy efficiency maximization; Section IV concerns the sum-rate optimization under (theoretical) perfect ZF beamforming; Section V contains details about the used multiuser channel model and the discussion of the numerical results. Finally, concluding remarks are given in Section VI.

*Notation. *We use lower case boldface characters to denote column vectors and capital boldface characters to denote matrices. We denote $N\times N$ identity matrix with $\textbf{I}_{N}$; $(\cdot)^{T}$, $(\cdot)^{H}$ and $(\cdot)^{*}$ denote the transpose, conjugate transpose and complex conjugate of a matrix. Then we denote with $\text{tr}(\textbf{A})$, $\parallel\textbf{A}\parallel$, $\mid\textbf{A}\mid$ the trace, the norm and the determinant of a matrix A. $\textbf{x}\sim\mathcal{CN}(\mu,\sigma^{2})$ indicates a Complex Gaussian random variable with mean $\mu$ and variance $\sigma^{2}$. $\otimes$ indicate the Kronecker product between two matrices. We denote with $\text{card}(\cdot)$ a set’s cardinality. A list of used acronyms is reported in Table I.

II System model and communication protocol

We consider an area where $K$ MSs and $M$ APs are randomly located. The APs are connected by means of a backhaul network to a CPU wherein data-decoding is performed (see fig. 1). Communications take place on the same frequency band; downlink and uplink are separated through TDD111In TDD, using calibrated hardware, the uplink channel is the reciprocal of downlink channel.. The communication protocol is made of three different phases: uplink training, downlink data transmission and uplink data transmission. During the uplink training phase, the MSs send pilot sequences to the APs and each AP estimates the channels; during the second phase the APs use the channel estimates to perform pre-coding and transmit the data symbols; finally, in the third phase the MSs send uplink data symbols to the APs. While in the CF approach all the APs simultaneously serve all the MSs (a fully-cooperative scenario), in the UC approach each AP serves a pre-determined number of MSs, say $N$, and in particular the ones that it receives best.

We assume that each AP (MS) is equipped with an antenna array with $N_{AP}$ ($N_{MS}$) elements. The ($N_{AP}\times N_{MS}$)-dimensional matrix $\textbf{H}_{k,m}$ denotes the channel matrix between the $k$-th user and the $m$-th AP. Details about the channel model will be reported in Section V

Additionally, we assume that each MS employs a very simple 0-1 beamforming structure; in particular, denoting by $P$ the multiplexing order, namely the number of parallel streams sent to a given receiver, the $(N_{MS}\times P)$-dimensional beamformer used at the $k$-th MS receiver is denoted by $\textbf{L}_{k}$ and is defined as $\textbf{L}_{k}=\textbf{I}_{P}\otimes\textbf{1}_{N_{MS}/P}$, denoting with $\textbf{1}_{N_{MS}/P}$ an all-1 vector of length $N_{MS}/P$. Otherwise stated, we assume that the MS receive antennas are divided in $P$ disjoint groups of $N_{MS}/P$ elements, and the received data collected at the antennas of each group are simply summed together. It is APs’ task, based on the uplink channel estimates and exploiting the TDD channel reciprocity, to ensure that the summed samples are, at least approximately, aligned in phase. Similarly when considering uplink transmission, the antennas in each group send the same signal with the same phase. We describe now the three phases of the communication protocol.

II-A Uplink training

During the uplink training the MSs transmit pilot sequences in order to enable channel estimation at the APs. Let $\tau_{c}$ be the length of the channel coherence time and $\tau_{p}$ be the length of uplink training phase, both in discrete time samples. Of course we must have $\tau_{p}<\tau_{c}$. We define by $\bm{\Phi}_{k}\in\mathcal{C}^{P\times\tau_{p}}$ the matrix containing on its rows the pilot sequences sent by the $k$-th MS. We assume that $\bm{\Phi}_{k}\bm{\Phi}_{k}^{H}=\textbf{I}_{P}$, i.e. the rows of $\bm{\Phi}_{k}$ are orthogonal, but no orthogonality is required for the pilot sequences assigned to other MSs222Of course, when $KP\leq\tau_{p}$ it would be possible to assign to all the MSs mutually orthogonal pilot sequences. In general, the pilot sequences contain symbols from the same constellation as the data symbols, e.g. QAM symbols or similar. In this paper, however, for the sake of simplicity, we assume that the pilot sequences are binary random sequences, and we just require that each matrix $\bm{\Phi}_{k}$ has orthogonal rows.. Obviously, using orthogonal pilots tout court would lead to a system immune to pilot contamination, but this would put a limit on the maximum value of the product $KP$ that could be accommodated in the channel coherence time. The received signal at the $m$-th AP in the $\tau_{p}$ signaling intervals devoted to uplink training can be cast in the following $N_{AP}\times\tau_{p}$-dimensional matrix $\textbf{Y}_{m}$:

[TABLE]

where $\textbf{S}_{k,m}=\textbf{H}_{k,m}\textbf{L}_{k}$, $\textbf{W}_{m}$ is the matrix of thermal noise samples, whose entries are assumed to be i.i.d. $\mathcal{CN}(0,\sigma_{w}^{2})$ RVs. In the following we briefly outline the structure of the LMMSE channel estimator. Defining $\textbf{y}_{m}=\text{vec}(\textbf{Y}_{m})$, $\textbf{w}_{m}=\text{vec}(\textbf{W}_{m})$, $\textbf{s}_{k,m}=\text{vec}(\textbf{S}_{k,m})$, we obtain the vectorized model:

[TABLE]

with $\textbf{A}_{k}=\Phi_{k}^{T}\otimes\textbf{I}$. Since we need to estimate a vector, we process $\textbf{y}_{m}$ by a matrix $\textbf{V}_{k,m}^{H}$, i.e. $\hat{\textbf{s}}_{k,m}=\textbf{{}V}_{k,m}^{H}\textbf{y}_{m}$. Then, the MSE is:

[TABLE]

Recalling that ([27, Ch. 4]) $\nabla_{\textbf{Z}^{*}}\text{tr}(\textbf{Z}^{T}\textbf{M})=\textbf{0}$, $\nabla_{\textbf{Z}^{*}}\text{tr}(\textbf{Z}^{H}\textbf{M})=\textbf{M}$, and $\nabla_{\textbf{Z}^{*}}\text{tr}(\textbf{Z}^{H}\textbf{M}\textbf{Z})=\textbf{MZ}$, setting the gradient of the MSE with respect to the complex matrix $\textbf{V}_{k,m}^{*}$ equal to zero and solving for $\textbf{V}_{k,m}$, we find the LMMSE estimator:

[TABLE]

II-B Downlink data transmission

After the first phase, the generic $m$-th AP has an estimate of the quantities $\textbf{S}_{k,m}$, for all $k=1,\ldots,K$. In order to transmit data on the downlink, a ZF precoder is considered. The precoding matrix $\mathbf{Q}_{k,m}$ is designed as follows. First of all, the matrix $\widehat{\mathcal{G}}=[\widehat{\textbf{S}}_{1,1}\dots\widehat{\textbf{S}}_{K,M}]$ is built; then, we have

[TABLE]

Finally, each precoding matrix is normalized as follows:

[TABLE]

$\forall k=1,\dots,K,\forall m=1,\dots,M$.

The previously described beamforming matrix requires a FD implementation, which presumes the use of a number of RF chains equal to the number of transmit antennas. It is well-known that at mmWave frequencies hardware complexity constraints usually prevent the use of FD architectures, and thus HY beamforming structures have been proposed, where a number of RF chains $N_{AP}^{RF}<N_{AP}$ is used. In this paper we exploit the Block Coordinate Descent algorithm [28] in order to decompose our beamformer in the cascade of a FD (baseband) one, represented by a $(N_{AP}^{RF}\times P)$-dimensional matrix and of an analog one, represented by a $(N_{AP}\times N_{AP}^{RF})$-dimensional matrix whose entries have all constant norm. While the baseband beamformer is MS-dependent, the analog beamformer at the AP is unique and is used for transmitting to all the users. In particular, at the generic $m$-th AP, the following matrix, of dimension $N_{AP}\times KP$, is formed333In the UC approach, to be detailed in the following, the matrix $\mathbf{Q}_{m}^{opt}$ is formed using the beamformers relative to the MS actually served by the $m$-th AP.

[TABLE]

and used as an input to the algorithm 1 in order to obtain the $N_{AP}\times N_{AP}^{RF}$ matrix $\mathbf{Q}_{m}^{RF}$ and the $N_{AP}^{RF}\times KP$ matrix

[TABLE]

containing the MS-specific baseband beamformers to be used at the $m$-th AP. At the generic $m$-th AP, we will thus have as many digital beamformers as the MSs to transmit to, and only one analog beamformer, that will be used to transmit jointly to all the users. Although we are here considering, for the sake of simplicity, only one single frequency, it should be noted that with a multicarrier modulation a single analog beamformer must be used for all the subcarriers, or, if complexity permits, for each properly defined subset of contiguous subcarriers [28].

II-B1 The CF approach

In this case, all the APs serve all the MSs, so the transmitted signal from the m-th AP in the n-th sample interval is:

[TABLE]

where $\textbf{x}_{k}^{DL}(n)$ is the data symbol intended for the $k$-th MS, and $\eta_{m,k}$ is a scalar coefficient taking into account the power used by the $m$-th AP to transmit towards the $k$.th MS. The $k-th$ MS receives the following ($N_{MS}\times 1$)-dimensional vector:

[TABLE]

where $\textbf{z}_{k}(n)$ is the additive thermal noise distributed as $\mathcal{CN}(0,\sigma_{z}^{2})$. A soft estimate of the $k$-th MS data symbol is thus formed as

[TABLE]

II-B2 The UC approach

In this case the APs are assumed to serve a pre-determined, fixed number of MSs, say $N$; in particular, we assume that the generic $m$-th AP serves the $N$ MSs whose channels have the largest Frobenious norms. We denote by $\mathcal{K}(m)$ the set of MSs served by the $m$-th AP. Given the sets $\mathcal{K}(m)$, for all $m=1,\dots,M$, we can define the set $\mathcal{M}(k)$ of the APs that communicate with the $k$-th user:

[TABLE]

In this case the transmitted signal from the $m$-th AP is written as:

[TABLE]

The received signal at the $k$-th MS is expressed now as:

[TABLE]

As before the $N_{MS}$-dimensional vector $\textbf{z}_{k}(n)$ represents the thermal noise at the $k$-th MS, and it is modeled as i.i.d. $\mathcal{CN}(0,\sigma_{z}^{2})$. Then it is possible to obtain a soft estimate of the data symbol $\textbf{x}_{k}^{DL}(n)$ at $k$-th MS as:

[TABLE]

Given the above equations, and assuming the use of Gaussian distributed long codewords, it is now possible to write down the achievable rate for the k-th user in the UC approach as follows:

[TABLE]

where

[TABLE]

and

[TABLE]

is the covariance matrix of the interference in the signal received at the $k$-th MS. Now, some algebraic manipulations are needed in order to express the achievable rate $\mathcal{R}_{k}$ in (15) in a form that permits solving the optimization problems considered in Sections III and IV. By letting:

[TABLE]

the covariance matrix (17) can be rewritten as

[TABLE]

thus implying that the rate for the k-th MS can be expressed as:

[TABLE]

Finally, exploiting the fact that the determinant of the product of two matrices factorizes into the product of determinants, we have:

[TABLE]

where the vector $\bm{\eta}$ is a compact notation to denote the set of all the downlink transmit powers $\eta_{m,k}$, for all the values of $k$ and $m$. All previous formulas concern the UC approach but they can be extended to the CF approach with ordinary efforts.

II-C Uplink data transmission

The third phase of the communication protocol amounts to uplink data transmission. We denote by $\textbf{x}_{k}^{UL}(n)$ the $P$-dimensional data vector to be transmitted by the $k$-th MS in the n-th sample time; the corresponding signal received at the m-th AP is expressed as:

[TABLE]

where $\widetilde{\eta}_{k}=\dfrac{P_{t,k}^{UL}}{\mbox{tr}(\textbf{L}_{k}^{H}\textbf{L}_{k})}$, and $P_{t,k}^{UL}$ is the uplink transmitted power by the $k$-th MS.

II-C1 CF approach

in this case, each AP forms the following statistic, $\forall k$:

[TABLE]

wherein now the ZF beamformer $\textbf{Q}_{k,m}$ previously defined is used as a post-coder. Then each AP sends to the CPU the vectors $\tilde{\textbf{y}}_{m,k}(n)$ via the backhaul link, and the CPU forms the following soft estimate of the data vectors transmitted by th $k$-th MS:

[TABLE]

Plugging (22) and (23) into (24) it can be checked that the signal contributions associated to the data-symbol $\textbf{x}_{k}^{UL}(n)$ are coherently summed.

II-C2 UC approach

in this case the signal transmitted by the k-th MS is decoded only by the APs belonging to the set $\mathcal{M}(k)$. Accordingly, the CPU performs the following soft estimate:

[TABLE]

Substituting (22) and (23) into (25), we get:

[TABLE]

Given the above equation, the achievable rate for the $k$-th MS can be easily shown to be expressed in the UC case as

[TABLE]

wherein the interference covariance matrix now is written as

[TABLE]

with $\mbox{card}(\cdot)$ denoting cardinality. The above expression can be used also for the CF case by letting ${\cal M}(k)=\{1,2,\ldots,M\}$, for all $k=1,\ldots,K$.

III Global Energy Efficiency Maximization

In this section we address the problem of power control for energy efficiency maximization considering both the downlink and the uplink channel of the described cellular network. The downlink case is treated in the coming Section III-A, while the uplink scenario is dealt with in Section III-B.

III-A Downlink power control

Mathematically the problem is formulated as the optimization program:

[TABLE]

where $P_{c,m}>0$ is the circuit power consumed at AP $m$, $\delta\geq 1$ is the inverse of the transmit amplifier efficiency, and $P_{max,m}$ is the maximum transmit power from AP $m$; $\bm{\eta}$ is a $KM\times 1$ vector containing all the transmit power of all AP.

Problem (29) is challenging due to its fractional objective, which has a non-concave numerator. This prevents the direct use of standard fractional programming methods such as Dinkelbach’s algorithm to solve (29) with affordable complexity. In addition, another issue is represented by the large number of optimization variables, i.e. $KM$. To counter both issues, we resort to the successive lower-bound maximization method, introduced in [29]444In [29] the method is termed successive upper-bound minimization since there the focus is on minimization problems. (it has already been used in other publications, [30] and [31]), and whose details are summarized in the Appendix. In brief, this method tackles (29) by alternatively optimizing the transmit powers of each AP, while keeping the transmit powers of the other APs fixed. However, as it will be shown, even with respect to the transmit powers of a single AP, the global energy efficiency problem remains challenging and for this reason, each subproblem is tackled by means of sequential optimization, as described next.

Consider thus problem (29) and define the variable blocks $\bm{\eta}_{m}=\left\{\eta_{m,k}\right\}_{k=1,\ldots,K}$, for the CF case, and $\bm{\eta}_{m}=\left\{\eta_{m,k}\right\}_{k\in{\cal K}(m)}$, for the UC case; the vector $\bm{\eta}_{m}$ is $K$-dimensional in the CF case, while in the UC approach it will contain $N$ entries, i.e. the powers to be used to communicate with the MSs in the set ${\cal K}(m)$. The global energy efficiency maximization problem with respect to the $m$-th variable block refers to the maximization, with respect to $\bm{\eta}_{m}$ of the quantity

[TABLE]

with $\bm{\eta}_{-m}$ representing the set of all the APs’ transmit powers except the ones in $\bm{\eta}_{m}$. It can be easily seen that the numerator of the global energy efficiency is not concave even with respect to only the variable block $\bm{\eta}_{m}$. In particular, (30) has a non-concave numerator, since $\mathcal{R}_{k}(\bm{\eta})$, which depends on all transmit powers of the base stations, is the difference of two logarithmic functions with the powers $\bm{\eta_{m}}$ of AP $m$ appearing in both logarithmic terms, and the difference of two logarithms is in general not concave. This means that (29) can not be solved with affordable complexity by either standard convex optimization nor fractional programming tools. Therefore, plain alternating maximization is not suitable, which motivates the use of the successive lower-bound maximization method. To this end, let us observe that each summand at the numerator of the global energy efficiency can be expressed as the difference between the functions $g_{1}$ and $g_{2}$ defined in (21). At this point we offer the following result.

Lemma 1**.**

Both functions $g_{1}$ and $g_{2}$ in (21) are concave in $\bm{\eta}_{m}$

Proof:

As a first step, let us show that the function $f:(x,y)\in\{\mathbb{R}_{0}^{+}\times\mathbb{R}_{0}^{+}\}\to\sqrt{xy}\in\mathbb{R}_{0}^{+}$ is jointly concave in $(x,y)$. To this end, the Hessian of $f$ can be written as

[TABLE]

which can be seen to be negative-semidefinite, since its determinant is zero and the elements on the diagonal are non-positive.

Next, we observe that both $g_{1}$ and $g_{2}$ in (21) are obtained by composing the function $f$ with the $\text{log}_{2}|(\cdot)|$ function of positive-definite argument. Then, the result follows recalling that the function $\text{log}_{2}|(\cdot)|$ is matrix-concave and matrix-increasing over the set of positive-semidefinite matrices [32, Section 3]. ∎

Based on Lemma 1, we argue that (21) is the difference between two concave function. Thus, recalling that any concave function is upper-bounded by its first-order Taylor expansion around any given point $\bm{\eta}_{m,0}$, a concave lower-bound of $\mathcal{R}_{k}$ can be obtained as:

[TABLE]

Then, a suitable approximate problem for the implementation of the sequential optimization method is

[TABLE]

It can be seen that Problem (33) fulfills the three properties P1, P2 and P3 of the sequential method detailed in the Appendix555The first-order Taylor expansion and its derivative coincide with the function and with its derivative when evaluated at the center of the expansion. Moreover, the numerator of the objective is now concave, which enables the use of fractional programming tools to globally solve (33), such as the popular Dinkelbach’s algorithm [21]. The overall power control procedure is formulated as in Algorithm 2 and, based on the properties of the successive lower-bound maximization method reviewed in the Appendix, we can state the following result.

Proposition 1**.**

Algorithm 2 monotonically improves the global energy efficiency value after each iteration and converges to a first-order optimal point of the original global energy efficiency maximization problem in (29).

III-A1 An alternative definition for the GEE on the downlink

The definition of the global energy efficiency reported in (29) considers a circuit power consumption that does not depend on the transmit power used by each base station. However, a base station that does not transmit will consume a lower circuit power, since it will switch to idle mode. We stress that such a circumstance can arise when maximizing the energy efficiency, since the marginal increase to the system achievable rate granted by their activation is overweighted by the increased network power consumption. To account for this circumstance, the terms $P_{c,m}$ can be further detailed to depend on the actual used transmit power, namely defining for all $m=1,\ldots,M$,

[TABLE]

where $P_{T}(m)=\sum_{k}\eta_{m,k}$ is the power radiated by the $m$-th AP and $\mathbbm{1}[P_{T}(m)>0]$ is the indicator function of the set $[P_{T}(m)>0]$, being $1$ when $P_{T}(m)>0$ and 0 otherwise. According to (34), we assume that the circuit power consumption, equal to $\widetilde{P}_{c,m}$ when the AP is active, is halved when it does not radiate any power.

While this more sophisticated power consumption model accounts for the lower circuit power consumption of idle access points, it leads to a non-differentiable global energy efficiency function, due to the presence of the indicator function I in (34). In this case Algorithm 2 is still guaranteed to monotonically increase the global energy efficiency value after each iteration, but no first-order optimality can be guaranteed upon convergence, owing to the non-differentiability of (34). Nevertheless, we should remark that it is also possible to approximate the indicator function in (34) by a smooth function, such as a sigmoid, thus recovering the first-order optimality property of the algorithm.

III-B Uplink power control

With regard to the uplink, the problem of global energy efficiency maximization is formulated as the optimization program:

[TABLE]

where $P_{T,max}$ denotes the maximum transmit power for the MSs666For the sake of simplicity this power is assumed to be the same for all the MSs; this assumption however can be easily relaxed.. It is easy to realize that the optimization problem (35) has the same structure as (29), thus implying that it can be solved by using the same procedure that has been developed for the downlink. It is also worth noting that while in the downlink the number of variables to be optimized is $MK$ in the CF case and $MN$ in the UC case, for the uplink only $K$ transmit powers are to be optimized, thus implying that the power optimization is much less computationally intensive on the uplink than on the downlink.

IV Sum-rate maximization

While the main focus of this work is on the maximization of the global energy efficiency function, it should be emphasized that the same approach can be adopted also for sum-rate maximization, which can be tackled as a special case of global energy efficiency maximization. Nonetheless, as detailed in the following, there are instances, when perfect CSI is available, in which the sum-rate maximization algorithm reduces to a standard convex optimization problem. Specifically, the global energy efficiency function reduces to the sum-rate function upon plugging $\delta=0$ and $P_{c}=1$ in (29a). Thus, in principle, sum-rate maximization can be performed by applying the simplified version of Algorithm 2 wherein standard convex optimization routines can be used to solve Problem (33) when $\mu=0$ and $P_{c}=1$ in each iteration777Note that when $\mu=0$ and $P_{c}=1$, (33) is no-longer a fractional program..

Nevertheless, sum-rate maximization has one peculiarity that should be highlighted with respect to the global energy efficiency maximization scenario. With reference to the case in which perfect CSI is available, it is possible to consider the ZF precoder

[TABLE]

thus removing multi-user interference. As a result, the sum-rate function simplifies to:

[TABLE]

where $\textbf{C}_{k,m,m^{\prime}}$ is:

[TABLE]

which is a concave function by virtue of Lemma 1. In this case, since the sum-rate is already concave and can be globally maximized with polynomial complexity by using standard convex programming theory.

V Performance analysis

V-A Simulation setup

We start by defining the used channel model for generating the matrices $\mathbf{H}_{k,m}$. According to the widely used clustered channel model for mmWave frequencies (see [26] and references therein), $\textbf{H}_{k,m}$ can be expressed as

[TABLE]

where $N_{cl}$ is the number of clusters, $N_{ray}$ is the number of rays that we consider for each cluster, $\gamma$ is a normalization factor defined as $\sqrt{\dfrac{N_{AP}N_{MS}}{N_{cl}N_{ray}}}$, $\textbf{H}_{LOS}$ is the line-of-sight (LOS) component, $\alpha_{i,l}$ is the complex path gain distributed as $\mathcal{CN}(0,1)$, so that its amplitude is Rayleigh-distributed, $L(r_{i,l})$ is the attenuation related to the path $(i,l)$, $\textbf{a}_{AP}$ and $\textbf{a}_{MS}$ are the array responses at the $m$-th AP and at the $k$-th MS, respectively, and they depend on the angles of arrival and departure, $\theta_{i,l,k,m}^{AP}$ and $\theta_{i,l,k,m}^{MS}$, relative to the $(i,l)$-th path of the channel between the $k$-th MS and the $m$-th AP. The path-loss is defined as [33]

[TABLE]

wherein $r$ is the distance between the transmitter and the receiver, $n$ is the path loss exponent, $X_{\sigma}$ is the shadow fading term in logarithmic units with zero mean and $\sigma^{2}$-variance and $f_{0}$ is a fixed frequency (see also table II). The $\textbf{H}_{LOS}$ in (39) is written as888For the ease of notation we omit the subscript $k,m$.

[TABLE]

In the above equation, $\eta\sim\mathcal{U}(0,2\pi)$, $I(d)$ is a 0-1 random variate indicating if a LOS link exists between the transmitter and the receiver, and $d$ is the transmitter-receiver distance, measured in meters. Denoting by $p$ the probability that $I_{LOS}(d)=1$, we have, for the UMi (Urban Microcellular) scenarios [34]:

[TABLE]

As for the number of scatterers and their positions, it should be said that, while usually for every AP-MS pair, a random and independently generated set of scatterers is considered to contribute to the channel matrix (39), in this paper, in order to model the possible channel correlation when the devices are closely spaced, we consider the same set of scatterers for the generation of all the channels. In particular, we assume that, in the considered area, there is a given number of random clusters, each one contributing with three rays999As specified in Section V, a density of 0.4 cluster/sqm. will be considered.. Given these clusters, in order to generate the generic channel $\textbf{H}_{k,m}$ between the $k$-th MS and the $m$-th AP, we consider as active only those clusters falling in an ellipse built around the position of the MS and the AP (see Fig. 2). In this way, on one hand we exclude far clusters from contributing to the channel, while, on the other hand, we are guaranteed that devices closely located will have correlated channels since they will be affected by similar sets of scatterers (see Fig. 2 for a graphical illustration).

In the following simulation we have considered a carrier frequency of $f_{0}=73\,\textrm{GHz}$, a bandwidth of $B=200\,\textrm{MHz}$ and, with regard to the above channel model, we have simulated the UMi Open Square scenario [34] of size $250\times 250\,\textrm{sqm}$. In order to generate the correlated channels, a total of 25.000 randomly deployed clusters (corresponding to a cluster density of $0.4\,\textrm{cluster/sqm}$.) has been generated. The additive white noise at the receiver has a power spectral density of $-174\,\textrm{dBm/Hz}$ and the receiver noise figure has been set to $F=6\,\textrm{dB}$. The simulated system has $M=80$ APs randomly deployed in the area to cover and equipped with an ULA with $N_{AP}=16$ antennas each. We consider both a scenario with $K=6$ users, that can be representative of a lightly-loaded network, and a scenario with $K=16$, that can be representative of a heavily-loaded scenario. We plot results for the UC networking deployment with $N=1$ and $N=3$. Each MS is equipped with and ULA with $N_{MS}=8$ antennas. Although the illustrated algorithms enable each MS to transmit and receive multiple data-streams, a multiplexing order of $P=1$ is assumed here for the sake of simplicity. The presented results show the global energy efficiency [Mbit/Joule] as defined in (29), with the circuit power consumption expressed as in (34), and the achievable rate per user [bit/s], defined as the sum-rate (i.e. the numerator of the global energy efficiency) divided by the number of MSs $K$. For benchmarking purposes, we compare the performance achieved by the proposed power control rules (maximizing the global energy efficiency and the sum-rate) with that achieved by downlink uniform power control, which assumes that $\eta_{m,k}=P_{T}/K$ in the CF case and that

[TABLE]

for the UC approach, respectively. The numerical values come from an average over 1000 independent channel scenarios as well as users and access points locations. The convergence criterion for the proposed sequential optimization algorithm is based on the computation of the relative tolerance, i.e. the norm of the difference between the current optimized vector and the optimized vector available at the previous iteration, divided by the norm of the current optimized vector. The successive convex approximation routine was directly implemented by the authors with a Matlab script, whereas the routine for maximizing the sum-rate in the case of concave rate was fmincon. All APs have been assumed to have the same maximum feasible transmit power $P_{max}$ on the downlink. The transmit amplifier efficiency of each transmitter has been assumed equal to one, i.e. $\delta=1$, while the hardware circuit power was modeled according to the model in (34) for each AP, with $\widetilde{P}_{c,m}=1\,\textrm{W}$, for all $m=1,\ldots,M$. For the uplink, we instead used $\widetilde{P}_{c,k}=0.3\,\textrm{W}$, for all $k=1,\ldots,K$. In the figures, the dashed line represent the results obtained by means of optimization procedure, instead, the solid line (it has not been reported in the legend in order to avoid a redundancy), with the same color, represent the case of uniform power allocation.

V-B Numerical results

Fig. 3 compares the global energy efficiency value versus $P_{max}$ achieved by the proposed global energy efficiency-maximizing power control scheme (labeled as OPT), considering the UC approach in the following scenarios:

• •

Perfect CSI and FD beamforming.

• •

Perfect CSI and HY beamforming, with $4$ RF chains used in the BCD-SD HY beamforming algorithm.

• •

Imperfect CSI and FD beamforming, with pilot sequences of length $\tau_{p}=64$ and uplink transmit power of $100\,\textrm{mW}$.

• •

Imperfect CSI and HY beamforming, with pilot sequences of length $\tau_{p}=64$ and uplink transmit power of $100\,\textrm{mW}$, with $4$ RF chains used in the BCD-SD HY beamforming algorithm.

The plot on the left refers to the UC rule wth $N=1$, while on the right we have the case $N=3$. Fig. 4 considers the same setting as Fig. 3 but reports the downlink average rate per user achieved by the proposed rate-maximizing power control scheme (labeled as OPT). Figs. 5 and 6 report the same results as Figs. 3 and 4, respectively, with the difference that they refer to the highly-loaded scenario, i.e. with $K=16$ users. Inspecting the figures, several considerations can be made. First of all, we see that the proposed power optimization method provides better performance than the uniform power allocation scheme. This is not always true for the case with imperfect CSI (see plot on the left in Fig. 5), as a consequence of the fact that, in the imperfect CSI case, the optimization step was performed using the estimated channels, while the plotted curve represents the true global energy efficiency, computed using the real channel coefficients. In other words, as far as the imperfect CSI scenarios are concerned, the metric that is optimized is different from the metric that is plotted. Next, as expected, we notice that FD beamforming and the availability of perfect CSI lead to better performance with respect to the practical situation that imperfect CSI and HY beamforming is to be accounted for. For instance, focusing on the plot on the left in Fig 4 and considering a maximum transmit power of $0\,\textrm{dBW}$, the average rate-per-user drops from $2.65\,\textrm{Gbit/s}$ of the ideal case to $1.65\,\textrm{Gbit/s}$, for the case in which both HY beamforming and ICSI are taken into account. Further, when comparing the lightly-loaded scenario with the highly-loaded one, we see that the system global energy efficiency does not dramatically change in the two situations, while, conversely, in a heavily-loaded system with HY beamforming and incomplete CSI the average rate-per-user at 0dBW of maximum transmit power is approximately equal to $500\,\textrm{Mbit/s}$, with about 60$\%$ loss with respect to the average rate-per-user attainable in a lightly loaded scenario. A peculiar behavior that is observed in Figs. 4 and 6 is that in the case of perfect CSI with FD beamforming with $N=1$, the uniform and optimal power allocation schemes achieve the same performance. This is due to the fact that the rate is an increasing function of the total available power $P_{max}$ and so, when only one user needs to be served, uniform power allocation coincides with the rate-maximizing power allocation strategy. Instead, this behavior is not observed in Figs. 3 and 5, where a visible gap is present between the performance with uniform power allocation and with the optimized power allocation also when $N=1$. This is expected because, unlike the rate, the global energy efficiency is not monotonically increasing with $P_{max}$.

Next, Fig. 7 shows the CDFs for the rate-per-user in case of FD and HY beamforming, respectively, considering $P_{max}=0\,\textrm{dBW}$ and the lightly loaded scenario. Again the two situations $N=1$ and $N=3$ are examined. The curves confirm the findings previously commented. Additionally, it can be seen that for the practical case of incomplete CSI and HY beamforming the scenario with $N=1$ provides much better performance than the one with $N=3$ for the vast majority of the users, i.e. it appears to be convenient to have each AP serve just one user than 3 users, except than when considering the unlucky users that happen to be situated on the left tail of the rate-per-user distribution. This behavior can be explained by noticing that in the considered setting, where $M>>K$, the majority of the MSs has in its neighborhood several APs, and letting these APs dedicate their own resources to only MS results in overall increased performance.

Finally, we turn our attention to the performance of the uplink channel. Fig. 8 reports the average uplink rate-per-user versus the uplink transmit power. No power control is performed here. Again, the scenarios $N=1$ and $N=3$ are considered. From the figure, it is seen that also for the uplink the rate is a clearly increasing function of the transmit power only in the ideal case of perfect CSI and fully digital beamforming, while, in the other situations, a saturation effect is observed. Results also show that, when focusing on the practical scheme with incomplete CSI and HY beamforming the case $N=1$ provides better results than the case $N=3$. In particular, assuming an uplink transmit power of -10dBW a data-rate of about $1.2\,\textrm{Gbit/s}$ can be achieved for the case $N=1$, versus a data-rate of about $600\,\textrm{Mbit/s}$ for the case $N=3$. Again, this behavior can be explained by noticing that in a scenario with a dense AP deployment each MS is sorrounded by several APs, and it is better to let each AP to focus its own resourced on just one MS rather then letting the AP to share its resources among multiple MSs.

VI Conclusion

This paper has considered a CF massive MIMO system operating at mmWave frequencies with a UC association between APs and MSs. Adopting a clustered channel model capable of taking into account the channel correlation for nearby devices, the paper has analyzed both the CF and UC approaches, by proposing a low-complexity power allocation rule aimed at global energy efficiency maximization. The presence of HY beamforming architectures at the APs has been considered, whereas, for the MSs, a simple channel-independent 0-1 beamforming structure has been considered. The obtained results have confirmed that the proposed resource allocation algorithms are effective at increasing the system energy efficiency and the system average rate-per-user, as well as that the use of HY beamforming architectures introduces, as expected, a considerable performance degradation. To the best of the authors’ knowledge, this is the first paper to consider the CF and UC massive MIMO architectures at mmWave frequencies. This study can be generalized along several paths. First of all, a comparison with a massive MIMO deployment with co-located antennas should be carried out. Then, given the gap observed between the HY and the FD beamforming schemes, more sophisticated HY beamformers should be considered to improve the system performance. Finally, given the intrinsic macro-diversity guaranteed by the dense AP deployment, it would be of interest to investigate on the suitability of CF architectures to circumvent the blockage effects that cause link instability: in this sense, CF massive MIMO could be envisioned as an architecture capable of making mmWave frequencies useful for ultra-reliable communications in indoor and small-sized outdoor environments.

Appendix: Successive lower-bound maximization

In this appendix we briefly review the successive lower-bound optimization technique. The method is based on the idea of merging the tools of alternating optimization [35, Section 2.7] and sequential convex programming [36]. Specifically, let us consider the generic optimization problem

[TABLE]

with $f:\mathbb{R}^{n}\to\mathbb{R}$ a differentiable function, and ${\cal X}$ a compact set. Similarly to the alternating optimization method, the successive lower-bound maximization partitions the variable space into $J$ blocks, $\bm{x}=(\bm{x}_{1},\ldots,\bm{x}_{J})$, which are cyclically optimized one at a time, while the other blocks are kept fixed. Thus, Problem (43) is decomposed into $M$ subproblems, wherein the generic subproblem is stated as

[TABLE]

with $\bm{x}_{-m}$ cotaining all variable blocks except the $m$-th. If Problem (44) is globally solved for each $m=1,\ldots,J$, then we have an instance of the alternating maximization method, which, as proved in [35, Proposition 2.7.1], monotonically improves the objective of (43), and yields a first-order optimal point if the solution of (44) is unique for any $m$, and if ${\cal X}={\cal X}_{1}\times\ldots\times{\cal X}_{J}$, with $\bm{x}_{m}\in{\cal X}_{m}$ for all $m$. However, if globally solving (44) is difficult (e.g. because (44) is not a convex problem), then implementing the alternating maximization method proves difficult. In this case, the successive lower-bound maximization method proposes to find a (possibly suboptimal) solution of (44), by means of the sequential convex programming method. Besides leading to a computationally viable algorithm, this approach is proved to preserve the optimality properties of the alternating optimization method [29], despite the fact that a possible suboptimal solution of (44) is determined.

As for sequential optimization, its basic idea is to tackle a difficult maximization problem by a sequence of easier maximization problems. To elaborate, denote by $g_{i}(\bm{x}_{m})$ the $i$-th constraint of (44), for $i=1,\ldots,C$, and consider a sequence of approximate problems $\{{\cal P}_{\ell}\}_{\ell}$ with objectives $\{f_{\ell}\}_{\ell}$ and constraint functions $\{g_{i,\ell}\}_{i=1}^{C}$, such that the following three properties are fulfilled, for all $\ell$:

1. (P1)

$f_{\ell}(\bm{x}_{m})\leq f(\bm{x}_{m})$, $g_{i,\ell}(\bm{x}_{m})\leq g_{i,\ell}(\bm{x}_{m})$, for all $i$ and $\bm{x}_{m}$; 2. (P2)

$f_{\ell}(\bm{x}_{m}^{(\ell-1)})=f(\bm{x}_{m}^{(\ell-1)})$, $g_{i,\ell}(\bm{x}_{m}^{(\ell-1)})=g_{i}(\bm{x}_{m}^{(\ell-1)})$ with $\bm{x}_{m}^{(\ell-1)}$ the maximizer of $f_{\ell-1}$; 3. (P3)

$\nabla f_{\ell}(\bm{x}_{m}^{(\ell-1)})=\nabla f(\bm{x}_{m}^{(\ell-1)})$, $\nabla g_{i,\ell}(\bm{x}_{m}^{(\ell-1)})=\nabla g_{i}(\bm{x}_{m}^{(\ell-1)})$.

In [36] (see also [37, 29]) it is shown that, subject to constraint qualifications, the sequence $\{f(\bm{x}_{m}^{(\ell)})\}_{\ell}$ of the solutions of the $\ell$-th Problem ${\cal P}_{\ell}$, is monotonically increasing and converges. Moreover, every convergent sequence $\{\bm{x}_{m}^{(\ell)}\}_{\ell}$ attains a first-order optimal point of the original Problem (44). Thus, sequential optimization guarantees at the same time to monotonically improve the objective function, and to fulfill the Karush Kuhn Tucker (KKT) first-order optimality conditions of the original problem. Nevertheless, in order to be able to use the method, it is necessary to find lower bounds of the original objective function, which fulfill all three properties P1, P2, P3, while at the same time leading to an approximate problem that can be solved with affordable complexity.