Joint Precoding and RRH selection for User-centric Green MIMO C-RAN

TL;DR

This paper presents a two-stage optimization framework for joint precoding and RRH selection in user-centric C-RANs, reducing power consumption while satisfying user rate and power constraints.

Contribution

It introduces a low-complexity user selection and precoding algorithm using re-weighted l1-norm minimization and WMMSE, achieving near-optimal power efficiency.

Findings

Rapid convergence of the proposed algorithms.

Multiple antennas at users improve performance.

Near-optimal network power consumption achieved.

Abstract

This paper jointly optimizes the precoding matrices and the set of active remote radio heads (RRHs) to minimize the network power consumption (NPC) for a user-centric cloud radio access network (C-RAN), where both the RRHs and users have multiple antennas and each user is served by its nearby RRHs. Both users' rate requirements and per-RRH power constraints are considered. Due to these conflicting constraints, this optimization problem may be infeasible. In this paper, we propose to solve this problem in two stages. In Stage I, a low-complexity user selection algorithm is proposed to find the largest subset of feasible users. In Stage II, a low-complexity algorithm is proposed to solve the optimization problem with the users selected from Stage I. Specifically, the re-weighted $l_1$-norm minimization method is used to transform the original problem with non-smooth objective function into a series of weighted power minimization (WPM) problems, each of which can be solved by the weighted minimum mean square error (WMMSE) method. The solution obtained by the WMMSE method is proved to satisfy the Karush-Kuhn-Tucker (KKT) conditions of the WPM problem. Moreover, a low-complexity algorithm based on Newton's method and the gradient descent method is developed to update the precoder matrices in each iteration of the WMMSE method. Simulation results demonstrate the rapid convergence of the proposed algorithms and the benefits of equipping multiple antennas at the user side. Moreover, the proposed algorithm is shown to achieve near-optimal performance in terms of NPC.

Equations

$$\begin{array}[]{l}{\bf{y}}_{k}={{\bf{\bar{H}}}_{k,k}}{{\bf{\bar{V}}}_{k}}{{\bf{s}}_{k}}+\sum\nolimits_{j\in{\cal U},j\neq k}{{\bf{\bar{H}}}_{j,k}}{{\bf{\bar{V}}}_{j}}{{\bf{s}}_{j}}+{\bf{n}}_{k},\end{array}$$

$$\begin{array}[]{l}{R_{k}}({\bf{V}})={\rm{lo}}{{\rm{g}}}\;\left|{{\bf{I}}+{{\bf{\bar{H}}}_{k,k}}{{\bf{\bar{V}}}_{k}}{\bf{\bar{V}}}_{k}^{H}{\bf{\bar{H}}}_{k,k}^{H}{\bf{J}}_{k}^{-1}}\right|,\end{array}$$

$$C1:R_k(\mathbf{V})\ge R_{k,\min },\forall k\in \mathcal{U}.$$

$$P_{i}^{{\rm{rrh}}}({\bf{V}})=\left\{\begin{array}[]{l}{{{\eta_{i}}}}P_{i}^{{\rm{tr}}}({\bf{V}})+MP_{i}^{{\rm{a,rrh}}},\ {\rm{if}}\ P_{i}^{{\rm{tr}}}({\bf{V}})>0\\ MP_{i}^{{\rm{s,rrh}}},\qquad\qquad\quad\ {\rm{if}}\ P_{i}^{{\rm{tr}}}({\bf{V}})=0\end{array}\right.$$

$$C2:P_i^{\mathrm{tr}}(\mathbf{V})\le P_{i,\max },\forall i\in \mathcal{I},$$

$$P_{i}^{{\rm{fr}}}({\bf{V}})=\left\{\begin{array}[]{l}{\rho_{i}}\sum\nolimits_{k\in{{\cal U}_{i}}}{R_{k}}({\bf{V}})+P_{i}^{a,{\rm{fr}}},\ {\rm{if}}\ P_{i}^{{\rm{tr}}}({\bf{V}})>0,\\ P_{i}^{s,{\rm{fr}}},\qquad\qquad\qquad\quad\quad\ {\rm{if}}\ P_{i}^{{\rm{tr}}}({\bf{V}})=0{\rm{.}}\end{array}\right.$$

$$\hat{P}(\mathcal{A},\mathbf{V})$$

$$=$$

$$\begin{array}[]{l}\mathop{\max}\limits_{{\bf{V}},{\cal U}\subseteq{\overline{\cal U}}}\quad\left|\cal U\right|\\ \quad{\rm{s}}.{\rm{t}}.\quad{\kern 1.0pt}{\rm{C1}},{\rm{C2}}.\end{array}$$

$$\underset{\mathcal{A},\mathbf{V}}{\min }$$

$$C1,{\sum }_{k\in {\mathcal{U}}_i^{\star }}{\Vert {\mathbf{V}}_{i,k}\Vert }_F^2\le P_{i,\max },i\in \mathcal{A},$$

$${\sum }_{k\in {\mathcal{U}}_i^{\star }}{\Vert {\mathbf{V}}_{i,k}\Vert }_F^2=0,i\in \mathcal{I}\backslash \mathcal{A},$$

$$\begin{array}[]{*{20}{l}}{\mathop{{\rm{min}}}\limits_{{{\left\{{{\alpha_{k}}}\right\}}_{k\in{\cal U}}},{\bf{V}}}\;\sum\nolimits_{k\in{\cal U}}{{{\left({{\alpha_{k}}-1}\right)}^{2}}}}\\ {\;{\rm{s}}.{\rm{t}}.\;\;\;{\rm{C2}},{R_{k}}({\bf{V}})\geq\alpha_{k}^{2}{R_{k,{\rm{min}}}},\forall k\in{\cal U},}\\ \end{array}$$

$${\hat{\mathbf{s}}}_k={\mathbf{U}}_k^H{\mathbf{y}}_k,\forall k\in \mathcal{U}.$$

$${\mathbf{E}}_k$$

$$h_k\left(\mathbf{V},{\mathbf{U}}_k,{\mathbf{W}}_k\right)=\log \left|{\mathbf{W}}_k\right|-\mathrm{Tr}\left({\mathbf{W}}_k{\mathbf{E}}_k\right)+d,\forall k.$$

$${\mathbf{U}}_k^{\star }={\left({\sum }_{j\in \mathcal{U}}{\bar{\mathbf{H}}}_{j,k}{\bar{\mathbf{V}}}_j{\bar{\mathbf{V}}}_j^H{\bar{\mathbf{H}}}_{j,k}^H+{\sigma }_k^2\mathbf{I}\right)}^{-1}{\bar{\mathbf{H}}}_{k,k}{\bar{\mathbf{V}}}_k,{\mathbf{W}}_k^{\star }=\mathbf{E}{}_k^{\star -1},\forall k,\forall k$$

$${\mathbf{E}}_k^{\star }={\mathbf{I}}_d-{\bar{\mathbf{V}}}_k^H{\bar{\mathbf{H}}}_{k,k}^H{\left({\sum }_{j\in \mathcal{U}}{\bar{\mathbf{H}}}_{j,k}{\bar{\mathbf{V}}}_j{\bar{\mathbf{V}}}_j^H{\bar{\mathbf{H}}}_{j,k}^H+{\sigma }_k^2\mathbf{I}\right)}^{-1}{\bar{\mathbf{H}}}_{k,k}{\bar{\mathbf{V}}}_k.$$

$$\begin{array}[]{*{20}{l}}{\mathop{{\rm{min}}}\limits_{{{\left\{{{\alpha_{k}}}\right\}}_{k\in{\cal U}}},{\bf{V}},{\bf{W}},{\bf{U}}}\;\sum\nolimits_{k\in{\cal U}}{{{\left({{\alpha_{k}}-1}\right)}^{2}}}}\\ {\;{\rm{s}}{\rm{.t}}{\rm{.}}\;\;\;{\rm{C2}},{h_{k}}\left({{\bf{V}},{{\bf{U}}_{k}},{{\bf{W}}_{k}}}\right)\geq\alpha_{k}^{2}{R_{k,{\rm{min}}}},\forall k\in{\cal U},}\end{array}$$

$$\begin{array}[]{*{20}{l}}{\mathop{\min}\limits_{\{\alpha_{k}\}_{k\in\cal U},{\bf{V}}}\;\sum\nolimits_{k\in{\cal U}}{{{\left({{\alpha_{k}}-1}\right)}^{2}}}}\\ \begin{array}[]{l}{\rm{s}}.{\rm{t}}.\;\ {\rm{C2}},{\rm{C5}}:{\rm{Tr}}\left({{{\left({{\bf{U}}_{k}^{H}{{\bf{\bar{H}}}_{k,k}}{{\bf{\bar{V}}}_{k}}-{{\bf{I}}_{k}}}\right)}^{H}}{{\bf{W}}_{k}}\left({{\bf{U}}_{k}^{H}{{\bf{\bar{H}}}_{k,k}}{{\bf{\bar{V}}}_{k}}-{{\bf{I}}_{k}}}\right)}\right)\\ \qquad\quad\quad+\sum\nolimits_{j\in{\cal U},j\neq k}{{\rm{Tr}}}\left({{\bf{\bar{V}}}_{j}^{H}{\bf{\bar{H}}}_{j,k}^{H}{{\bf{U}}_{k}}{{\bf{W}}_{k}}{\bf{U}}_{k}^{H}{{\bf{\bar{H}}}_{j,k}}{{\bf{\bar{V}}}_{j}}}\right)+{\alpha_{k}^{2}}{R_{k,\min}}\leq{t_{k}},\forall k\in\cal U,\end{array}\end{array}$$

$$\begin{array}[]{l}\mathop{\min}\limits_{\{\alpha_{k}\}_{k\in\cal U},{\bf{V}}}\;\sum\nolimits_{k\in{\cal U}}{{{\left({{\alpha_{k}}-1}\right)}^{2}}}\\ {\rm{s}}{\rm{.t}}{\rm{.}}\ {\left\|{{{\bf{x}}_{k}}}\right\|_{2}}\leq\sqrt{{t_{k}}},\forall k\in\cal U,\\ \quad\ \ {\left\|{{{\bf{y}}_{i}}}\right\|_{2}}\leq\sqrt{{P_{i.\max}}},\forall i\in\cal I,\end{array}$$

$$\begin{array}[]{l}{{\bf{x}}_{k}}=\left[{{\rm{vec}}{{\left({{\bf{\bar{V}}}_{1}^{H}{\bf{\bar{H}}}_{1,k}^{H}{{\bf{U}}_{k}}{\bf{W}}_{k}^{1/2}}\right)}^{H}},\cdots,}\right.{\rm{vec}}{\left({\left({{{\bf{\bar{V}}}_{k}^{H}}{{\bf{\bar{H}}}_{k,k}^{H}}{\bf{U}}_{k}-{{\bf{I}}_{k}}}\right){\bf{W}}_{k}^{1/2}}\right)^{H}},\\ \qquad\qquad\qquad\qquad\qquad\qquad{\left.{{\rm{}}\cdots,{\rm{vec}}{{\left({{\bf{\bar{V}}}_{K}^{H}{\bf{\bar{H}}}_{K,k}^{H}{{\bf{U}}_{k}}{\bf{W}}_{k}^{1/2}}\right)}^{H}},\alpha_{k}\sqrt{{R_{k,\min}}}}\right]^{H}}\end{array}$$

$${\mathbf{y}}_i={\left[\mathrm{vec}{\left({\mathbf{V}}_{i,q_1^i}\right)}^H,\cdots ,\mathrm{vec}{\left({\mathbf{V}}_{i,q_{|{\mathcal{U}}_i|}^i}\right)}^H\right]}^H.$$

$${\Vert {\sum }_{k\in {\mathcal{U}}_i}{\Vert {\mathbf{V}}_{i,k}\Vert }_F^2\Vert }_0\approx a_i^{(n)}{\sum }_{k\in {\mathcal{U}}_i}{\Vert {\mathbf{V}}_{i,k}\Vert }_F^2,\text{vspace}-0.2cm$$

$$a_i^{(n)}=\frac{1}{{\sum }_{k\in {\mathcal{U}}_i}{\Vert {\mathbf{V}}_{i,k}^{(n)}\Vert }_F^2+\delta },\forall i,$$

$$\begin{array}[]{l}\mathop{\min}\limits_{\bf{V}}\ \sum\nolimits_{i\in\cal I}{\omega_{i}^{(n-1)}\sum\nolimits_{k\in{{\cal U}_{i}}}{\left\|{{{\bf{V}}_{i,k}}}\right\|_{F}^{2}}}\\ {\rm{s}}{\rm{.t.}}\quad{\rm{C1,C2}},\end{array}$$

$$\begin{array}[]{l}\mathop{{\rm{min}}}\limits_{{\bf{V,W,U}}}\sum\nolimits_{k\in{\cal U}}{{\rm{Tr}}\left({{\bf{\bar{V}}}_{k}^{H}}{{\bf{G}}_{k}}{{{\bf{\bar{V}}}}_{k}}\right)}\\ \;\text{s.t.}\;\;\;{h_{k}}\left({{\bf{V}},{\bf{U}}_{k},{\bf{W}}_{k}}\right)\geq{R_{k,{\text{min}}}},\forall k\in\cal U,\\ \quad\quad\;\sum\nolimits_{k\in{{\cal U}_{i}}}{\left\|{{{\bf{B}}_{i,k}}{{{\bf{\bar{V}}}}_{k}}}\right\|}_{F}^{2}\leq{P_{i,\max}},\forall i\in{\cal I},\end{array}\vspace{-0.2cm}$$

$${\mathbf{G}}_k=\mathrm{blkdiag}\left\{{\omega }_{s_1^k}{\mathbf{I}}_{M\times M},\cdots ,{\omega }_{s_{|{\mathcal{I}}_k|}^k}{\mathbf{I}}_{M\times M}\right\}\text{vspace}-0.2cm$$

$${\mathbf{B}}_{i,k}=\mathrm{diag}\left\{\stackrel{s_1^k}{\overbrace{0_{1\times M}}},\cdots ,\stackrel{s_j^k}{\overbrace{1_{1\times M}}},\stackrel{s_{j+1}^k}{\overbrace{0_{1\times M}}},\cdots ,\stackrel{s_{|{\mathcal{I}}_k|}^k}{\overbrace{0_{1\times M}}}\right\},\mathrm{if}{s}_j^k=i,\forall i\in \mathcal{I},k\in \mathcal{U}.\text{vspace}-0.2cm$$

$$\begin{array}[]{l}\mathop{{\rm{min}}}\limits_{{\bf{V}}}\ \sum\nolimits_{k\in{\cal U}}{{\rm{Tr}}\left({{\bf{\bar{V}}}_{k}^{H}}{{\bf{G}}_{k}}{{{\bf{\bar{V}}}}_{k}}\right)}\\ {\rm{s}}{\rm{.t}}{\rm{.}}\ \sum\nolimits_{j\in{\cal U}}{{\rm{Tr}}\left({{\bf{\bar{V}}}_{j}^{H}{\bf{\bar{H}}}_{j,k}^{H}{{\bf{U}}_{k}}{{\bf{W}}_{k}}{\bf{U}}_{k}^{H}{{{\bf{\bar{H}}}}_{j,k}}{{{\bf{\bar{V}}}}_{j}}}\right)}\!\!-\!\!{\rm{Tr}}\left({{{\bf{W}}_{k}}{\bf{U}}_{k}^{H}{{{\bf{\bar{H}}}}_{k,k}}{{{\bf{\bar{V}}}}_{k}}}\right)\!\!-\!\!{\rm{Tr}}\left({{\bf{\bar{V}}}_{k}^{H}{\bf{\bar{H}}}_{k,k}^{H}{{\bf{U}}_{k}}{{\bf{W}}_{k}}}\right)\!\leq\!{c_{k}},\forall k\\ \quad\ \ \sum\nolimits_{k\in{{\cal U}_{i}}}{{\rm{Tr}}}\left({{\bf{\bar{V}}}_{k}^{H}{{\bf{B}}_{i,k}}{{{\bf{\bar{V}}}}_{k}}}\right)\leq{P_{i,\max}},\forall i,\end{array}$$

$$T_{\mathrm{SOCP}}=O\left(\sqrt{I+K}{\left(2lMKd\right)}^2\left(2K^2{d}^2+2dMKl\right)\right).$$

$$\mathcal{L}\left(\mathbf{V},\mathbf{\lambda },\mathbf{\mu }\right)$$

Full text

I Introduction

Mobile communications has been developing very rapidly [2, 3, 4]. In recent years, C-RAN has been proposed as a promising solution to support the exponential growth of mobile data traffic [5, 6]. In C-RAN, all the baseband processing is performed at the baseband unit (BBU) pool with powerful computation capacity, while the remote radio heads (RRHs) perform the basic functionalities of signal processing [7, 8]. The RRHs are geographically distributed away from each other, but connected to the BBU pool through optical fiber transport links. Under the C-RAN architecture, centralized signal processing technologies can be realized. Hence, significant performance gains can be achieved. In addition, the RRHs can be densely deployed in the network with low operation cost due to their simple functionalities. This will significantly reduce the average access distance for the users, and thus lowers the transmission power.

On the other hand, it was reported that the total energy consumption of wireless communications contributes more than 3 percent of the worldwide electrical energy consumption [9], and this portion is expected to grow in the near future due to the explosive growth of high-data-rate applications and mobile devices. Hence, energy efficiency has attracted extensive interest and becomes one of the main performance metrics in the future fifth generation (5G) systems [10, 11, 12, 13, 14]. When a large number of RRHs are deployed in the network, the network power consumption (NPC) of C-RAN will become considerable due to the increasing circuit power consumption of the RRHs. Fortunately, it was reported in [15] that the traffic load varies substantially over both time and space due to user mobility and varying channel state. Hence, the NPC can be significantly reduced by putting some RRHs with light load into sleep mode while maintaining the quality of service (QoS) requirements of the users, which is the focus of this paper.

Recently, the NPC minimization problem for C-RAN has been extensively studied in [16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26]. These papers formulated the joint RRH selection and beamforming vector optimization problem as a mixed-integer non-linear programming (MINLP) problem, which has a nonconvex discontinuous $l_{0}$-norm in the objective function or constraints. We summarize the existing approaches to solve the MINLP problem as follows. The first approach was proposed in [16], which first reformulated the problem as an extended mixed integer second-order cone programming (SOCP) and then applied the branch-and-cut method to obtain the optimal solution. In the second approach in [17, 18], the MINLP was first decomposed into a master problem and a beamforming subproblem. Then, an iterative algorithm based on the Benders decomposition was derived to find the optimal solution. Although these two approaches yield the optimal solution, they have an exponential complexity. The third approach is the smooth function method, where the $l_{0}$-norm was approximated as Gaussian-like function in [19], the exponential function in [20], and arctangent function in [21]. However, the smooth function cannot produce sparse solutions in general. The last approach was inspired by the compression sensing, named re-weighted $l_{1}$-norm minimization method [27]. This method has been widely adopted in the literature [28, 22, 23, 24, 25, 26] due to its low computational complexity and sparsity guarantee, which will also be applied in this paper.

All of the above papers only considered the single-antenna user (SAU) case. With the increasing development in antenna technology [29, 30], it is possible to equip the wireless devices with multiple antennas. When both the transmitter and the receiver are equipped with multiple antennas, multiple streams can be transmitted simultaneously, rather than only one stream in the SAU case. Simulation results show that with the increasing number of receive antennas, more users can be admitted. Therefore, in this paper, we consider the multiple-antenna user (MAU) case and jointly optimize the precoding matrices and the set of active RRHs to minimize the NPC subject to users’ rate requirements and per-RRH power constraints.

Unfortunately, the techniques in [16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26] dealing with the SAU case cannot be extended directly to the MAU case. The reasons are as follows. Firstly, since the rate constraints and power constraints are conflicting with each other, this problem may be infeasible. In the SAU networks, the rate requirements can be equivalently represented as signal-to-interference-plus-noise ratio (SINR) constraints, which can be transformed into an SOCP problem. Hence, the feasibility of the original problem can be easily checked by solving the SOCP feasibility problem. However, the rate constraints in the MAU case is non-convex and much more complex due to the complicated rate expression, which cannot be transformed into the SOCP formulation as in the SAU case. Hence, new techniques need to be developed to check the feasibility of the original problem. Secondly, even though the original problem is checked to be feasible, how to solve it is still difficult, since it cannot be transformed into an SOCP problem as in the SAU case. [31] proposed the weighted minimum mean square error (WMMSE) method to solve the rate maximization problem for MIMO interfering broadcast channels, where the rate expression is in the objective function. Recently, there have been some work in applying the WMMSE method to solve the energy efficiency (measured in bit/s/Joule) optimization problems under rate constraints [32, 33]. However, these researches have not addressed the feasibility problem due to the incorporated rate constraints. Only in [33], a heuristic method was proposed to check the feasibility based on the interference alignment technique, under the assumption that the transmit power is approaching infinity, which in not practical. Since the problem considered in this paper imposes power constraints at each RRH, the heuristic method developed in [33] is not applicable. More importantly, they have not revealed the hidden property of applying WMMSE method to the optimization problem with rate constraints, such as the convergence property and the optimality of the solutions.

To the best of our knowledge, this paper is the first attempt to solve the joint RRH and precoding optimization problem to minimize the NPC in the MAU based user-centric C-RAN, where each user can be served by an arbitrary subset of RRHs. Due to the conflicting constraints, this problem may be infeasible. Some users should be removed or rescheduled for the next transmission to guarantee the rate requirements of other users. We provide a comprehensive analysis for this problem by considering two stages: user selection in Stage I and algorithm design in Stage II. The main contributions of this paper are summarized as follows:

In Stage I, a low-complexity user selection approach is proposed to maximize the number of admitted users that can have their QoS requirements satisfied. Specifically, in each step we solve an alternative problem by introducing a series of auxiliary variables. This alternative problem is always feasible. By replacing the rate expression in the constraints with its lower-bound, an iterative algorithm is proposed to solve this problem along with the complexity and convergence analysis of the algorithm. The alternative problem should be solved at most $K$ times, where $K$ is the total number of users. Its complexity is much lower than the optimal exhaustive user selection method that has an exponential complexity. Simulation results show that both algorithms achieve similar performance. 2. 2.

In Stage II, a low-complexity algorithm is proposed to solve the NPC minimization problem with the users selected from Stage I. Specifically, the re-weighted $l_{1}$-norm minimization method is adopted to convert the non-smooth optimization problem into a series of smooth weighted power minimization (WPM) problems. We again replace the rate expression with its lower-bound and adapt the WMMSE algorithm originally designed for a rate maximization problem to solve the WPM problem. In addition, we strictly prove that when the WMMSE algorithm is initialized with a feasible solution, the sequences of precoder matrices generated in the iterative procedure will finally converge to the Karush-Kuhn-Tucker (KKT) point of the WPM problem. 3. 3.

In each iteration of the WMMSE algorithm, there is a subproblem for the precoder matrices being updated with some other fixed variables. Most existing papers [32, 33, 34, 25, 35] directly transform it into an SOCP problem and apply the interior point method [36] to solve it, which may incur high computational complexity. In this paper, we go one step further and develop a low-complexity algorithm to solve this subproblem by exploiting its special structure. Specifically, we equivalently solve its dual problem because the subproblem is a convex problem. Fortunately, the objective function of the dual problem is differentiable, and the block coordinate descent (BCD) method is adopted to solve the dual problem. In each iteration of the BCD method, Newton’s method and the gradient descent method are applied to update the Lagrangian multipliers. It is strictly proved that the BCD method can obtain the globally optimal solution of the subproblem. Complexity analysis in conjunction with the simulation results show that the BCD method has a much lower computational complexity than the interior point method.

This paper is organized as follows. In Section II, we introduce the system model and formulate the optimization problems. In Section III, a new approach is introduced to select the maximum number of admitted users. An iterative algorithm with low complexity is provided in Section IV. Simulation results are presented in Section V. Conclusions are drawn in Section VI.

Notations: Uppercase and lowercase boldface denote matrices and vectors, respectively. For a matrix ${\bf{A}}$, ${\left\|{\bf{A}}\right\|_{F}}$ denotes the Frobenius norm of ${\bf{A}}$ and ${{\bf{A}}^{H}}$ represents the Hermitian transpose of ${\bf{A}}$. ${\bf{I}}_{m}$ denotes a $m\times m$ identity matrix. For a vector ${\bf{a}}$, ${\rm{diag}}(\bf{a})$ denotes the diagonal matrix with diagonal elements given by $\bf{a}$. ${\rm{blkdiag}}(\bf{\cdot})$ represent the block-diagonal matrices. $\mathbb{E}(\cdot)$, and ${\rm{Tr}}(\cdot)$ represent expectation, trace operators, respectively. ${\bf{A}}\succeq{\bf{B}}$ means ${\bf{A}}-{\bf{B}}$ is a positive semidefinite matrix. For vector ${\bf{a}}\in{\mathbb{C}^{n\times 1}}$, ${\left\|{\bf{a}}\right\|_{2}}$ is the Euclidean norm. ${\cal C}{\cal N}\left({{\bf{0}},\sigma^{2}{{\bf{I}}}}\right)$ represents the complex circularly symmetric Gaussian distribution with zero mean vector and covariance matrix $\sigma^{2}{{\bf{I}}}$. For a vector ${\bf{x}}$, ${\left\|{\bf{x}}\right\|_{0}}$ is ${l_{0}}$-norm, means the number of nonzero entries in a vector.

II System Model and Problem Formulation

II-A System model

Consider a downlink C-RAN consisting of $I$ RRHs and $K$ users 111In dense networks, the number of RRHs may be larger than the number of users so that the average distance between serving RRHs and users can be significantly reduced, leading to improved performance. In some extreme cases, each user may be served by its dedicated RRHs as in [37, 38], where each RRH serves only one user., where each RRH is equipped with $M$ transmit antennas and each user has $N$ receive antennas, as shown in Fig. 1. Denote the set of RRHs and users as ${\cal I}=\left\{{1,\cdots,I}\right\}$ and $\bar{\cal U}=\left\{{1,\cdots,K}\right\}$, respectively. It is assumed that each RRH is connected to the BBU pool via fronthaul link and the BBU pool has access to all users’ CSI and data information.

Let ${{\cal U}}\subseteq\bar{\cal U}$ be the set of users that can be admitted to this networks. To reduce the computational complexity of the dense network, the user-centric clustering method is adopted, where each user $k\in{\cal U}$ is assumed to be served by its nearby RRHs since the distant RRHs contribute less to user’s signal quality due to the large path loss. The unselected RRHs are turned into idle mode, such as RRH 3 and RRH 5 in Fig. 1. Let ${{\cal I}_{k}}\subseteq{\cal I}$ and ${{{\cal U}}_{i}}\subseteq{{\cal U}}$ be the candidate set of RRHs for serving user $k$ and candidate set of users served by RRH $i$, respectively. Note that the set of RRHs serving the users may overlap with each other. For example, in Fig. 1, RRH 12 jointly serves user 1 and user 6.

Denote ${\bf{V}}_{i,k}\in\mathbb{C}^{M\times d}$ as the precoding matrix used by the $i$th RRH to transmit data vector ${\bf{s}}_{k}\in\mathbb{C}^{d\times 1}$ to the $k$th user, where $d$ is the number of data streams for each user, and ${\bf{s}}_{k}$ satisfies $\mathbb{E}\left[{{{\bf{s}}_{k}}{\bf{s}}_{k}^{\rm{H}}}\right]={\bf{I}}_{d}$ and $\mathbb{E}\left[{{{\bf{s}}_{k}}{\bf{s}}_{l}^{\rm{H}}}\right]={\bf{0}},{\rm{for}}\ l\neq k$. Let ${{\bf{\bar{V}}}_{k}}={\left[{{\bf{V}}_{i,k}^{\rm{H}},\forall i\in{{\cal I}_{k}}}\right]^{\rm{H}}}\in{\mathbb{C}^{\left|{{{\cal I}_{k}}}\right|M\times d}}$ be the big precoding matrix for user $k$ from all RRHs in ${\cal I}_{k}$. In addition, define a set of new channel matrices ${\bf{\bar{H}}}_{j,k}=[{{\bf{H}}_{i,k},\forall i\in{{\cal I}_{j}}}]\in{{\mathbb{C}}^{N\times\left|{{{\cal I}_{l}}}\right|M}}$, representing the overall CSI from RRHs in ${{\cal I}_{j}}$ to user $k$, where ${{{\bf{H}}_{i,k}}}\in{\mathbb{C}^{N\times M}}$ denotes the channel matrix from the $i$th RRH to the $k$th user. Then, the received signal vector at the $k$th user, denoted as ${\bf{y}}_{k}\in\mathbb{C}^{N\times 1}$, is given by

[TABLE]

where ${\bf{n}}_{k}$ is the noise vector at the $k$th user, which satisfies ${\cal C}{\cal N}\left({{\bf{0}},\sigma_{k}^{2}{{\bf{I}}_{N}}}\right)$. Then, the achievable rate (nat/s/Hz) of the $k$th user is given by[39]

[TABLE]

where ${\rm{log}}(\cdot)$ is the base of natural logarithm, ${\bf{J}}_{k}=\sum\nolimits_{j\in{\cal U},j\neq k}{{\bf{\bar{H}}}_{j,k}}{{\bf{\bar{V}}}_{j}}{{\bf{\bar{V}}}_{j}^{H}}{{\bf{\bar{H}}}_{j,k}^{H}}+{\sigma^{2}_{k}}\bf{I}$ is the interference-plus-noise covariance matrix, and ${\bf{V}}$ is the collection of all precoding matrices. Each user’s data rate should be larger than the minimum requirement:

[TABLE]

With densely deployed RRHs, the power consumption on the RRHs and the corresponding fronthaul links may be significant. Switching off some RRHs and the corresponding fronthual links may be a good option to reduce the NPC. To this end, it is critical to model the NPC.

II-B NPC model

The realistic NPC model should consist of three parts: power consumption at the RRHs, that at the fronthaul links and that at the BBU pool.

As in [22], the power consumption at RRH $i$ can be modeled as follows:

[TABLE]

where $\eta_{i}>1$ accounts for the inefficiency of the power amplifier of RRH $i$, $P_{i}^{{\rm{tr}}}({\bf{V}})$ is the total transmit power of RRH $i$ given by $P_{i}^{{\rm{tr}}}({\bf{V}}){\rm{=}}\sum\nolimits_{k\in{{\cal U}_{i}}}{\left\|{{{\bf{V}}_{i,k}}}\right\|_{F}^{2}}$ that satisfies the power constraint:

[TABLE]

$P_{i}^{{\rm{a,rrh}}}$ and $P_{i}^{{\rm{s,rrh}}}$ represent the power consumption for each antenna (or each RF chain) when RRH $i$ is in active mode and sleep mode, respectively. In practical systems, $P_{i}^{{\rm{active}}}$ is much higher than $P_{i}^{{\rm{sleep}}}$, which motivates us to switch off some RRHs.

In general, more power consumption will be consumed on the fronthaul links when they support high data rates. In [26], this power was modeled to be proportional to the total fronthaul data rate. We modify the model in [26] to account for the power when the fronthaul links are in the sleep mode as follows:

[TABLE]

where $\rho_{i}$ is the proportional factor for fronthaul link $i$. The power consumed in the BBU pool mainly depends on the computational complexity for signal processing. However, how to accurately model this kind of power consumption is still not fully understood. As in most papers [16, 22, 23, 26], the BBU power consumption is modeled as a constant ${P_{{\rm{BBU}}}}$ for simplicity. Let ${\cal A}$ denote the active RRH set. Then, the NPC can be modeled as

[TABLE]

where $P_{i}^{c}$ and $P_{i}^{s}$ are two constants, given by $P_{i}^{c}=M(P_{i}^{a,{\rm{rrh}}}-P_{i}^{s,{\rm{rrh}}})+P_{i}^{a,{\rm{fr}}}-P_{i}^{s,{\rm{fr}}}$ and $P_{i}^{s}=MP_{i}^{s,{\rm{rrh}}}+P_{i}^{s,{\rm{fr}}}$.

II-C Problem Formulation

Due to the power constraints C2, the rate requirements C1 may not be satisfied for all users. Some users should be removed to make the optimization problem feasible. Hence, we formulate a two-stage optimization problem. In Stage I, one should maximize the number of admitted users that can be supported by the system; in Stage II, one should jointly select some RRHs and optimize the precoding matrices to minimize the NPC with the selected users from Stage I.

Specifically, the optimization problem in Stage I can be formulated as

[TABLE]

Then in Stage II, we aim to jointly select the RRHs and optimize the precoding matrices to minimize the NPC with the users selected from Stage I, which can be formulated as222In general, the number of transmit antennas should be optimized to additionally reduce the NPC as seen in the RRH power consumption model in (4). However, the resulting problem will be much more difficult to solve, and will be left for future work.

[TABLE]

where ${\cal U}_{i}^{\star}$ is the solution from Stage I. Note that when the system parameters are given, the last two terms in (8) are constants, and are omitted in the objective function.

Both the optimization problems in the two stages are MINLP problems and are difficult to solve. The intuitive approach to solve this kind of problems is through the exhaustive search. For example, to solve the NPC minimization problem in Stage II, one must solve the precoding matrices that minimizes the NPC with each given $\cal A$ and obtain the corresponding objective value. Finally, the $\cal A$ that achieves the minimum NPC together with the corresponding precoding matrices is the optimal solution of Problem (10). However, the exhaustive search has exponentially prohibitive complexity with respect to the number of RRHs, which is hard to be implemented in practice in dense C-RANs. The same issue holds for the user selection problem in Stage I, where the exhaustive search method has an exponential complexity of the number of users. Hence, this motivates us to develop low-complexity algorithms to solve these two Problems.

III Stage I: Low-complexity User Selection Algorithm

In this section, we provide a low-complexity user selection algorithm to guarantee the rate requirements of other users. Specifically, for an arbitrary given subset of users $\cal U$, we construct an alternative problem by introducing a series of auxiliary variables $\{\alpha_{k}\}_{k\in\cal U}$:

[TABLE]

Obviously, Problem (11) is always feasible and the optimal $\alpha_{k}$ for each user $k$ should be no larger than one. This can be easily proved by contradiction. Moreover, user $k$ can be admitted if and only if the optimal $\alpha_{k}$ is equal to one. Hence, maximizing the number of admitted users is equal to finding the largest subset of users $\cal U$, in which all $\{\alpha_{k}\}_{k\in\cal U}$ are equal to one.

Based on the above analysis, we provide a low-complexity user selection (USC) algorithm to solve Problem (9) in Stage I. The main idea is to remove each user with the least $\alpha_{k}<1$ in each iteration. It is intuitive since the user with the least $\alpha_{k}$ has the largest gap to its rate target.

III-A Algorithm to solve Problem (11)

In step 2 of Algorithm 1, Problem (11) needs to be solved. Due to constraints C3 in (11), Problem (11) is a non-convex problem, which is difficult to solve. To handle this difficulty, we apply the relationships between WMMSE and the rate expression.

We consider the linear receiver filter so that the estimated signal vector is given by

[TABLE]

where ${\bf{U}}_{k}\in\mathbb{C}^{N\times d}$ is the receiver filter of the $k$th user. Since the signal vectors ${\bf{s}}_{k}$’s and noise ${\bf{n}}_{k}$’s are mutually independent, the mean square error (MSE) matrix at the $k$th user is given by

[TABLE]

By introducing a set of auxiliary matrices $\{{{\bf{W}}_{k}}\succeq{\bf{0}}\}$, we define the following functions

[TABLE]

where ${\bf{E}}_{k}$ is the MSE matrix of user $k$ given in (13). The following lemma establishes the relationships between the rate expression and function ${h_{k}}\left({\bf{V}},{{\bf{U}}_{k},{\bf{W}}_{k}}\right)$.

Lemma 1 [31] : ${h_{k}}\left({\bf{V}},{{\bf{U}}_{k},{\bf{W}}_{k}}\right)$ is a concave function for each set of the matrices ${\bf{V}}$, ${\bf{U}}_{k}$ and ${\bf{W}}_{k}$ when the other two are given. Given ${\bf{V}}$, ${h_{k}}\left({\bf{V}},{{\bf{U}}_{k},{\bf{W}}_{k}}\right)$ is the lower-bound of the data rate ${R_{k}}({\bf{V}})$ in (2). The optimal ${{\bf{U}}_{k},{\bf{W}}_{k}}$ for ${h_{k}}\left({\bf{V}},{{\bf{U}}_{k},{\bf{W}}_{k}}\right)$ to achieve the data rate is given by

[TABLE]

where ${\bf{E}}_{k}^{\star}$ is obtained by plugging the expression of ${\bf{U}}_{k}^{\star}$ into the $k$th user’s MSE in (13)

[TABLE]

$\Box$

By replacing the first set of constraints in (11) with its lower-bound ${h_{k}}\left({\bf{V}},{{\bf{U}}_{k},{\bf{W}}_{k}}\right)$, we have the following optimization problem

[TABLE]

where ${\bf{U}}$ and ${\bf{W}}$ are the collection of matrices ${{\bf{U}}_{k},\forall k}$ and ${{\bf{W}}_{k},\forall k}$, respectively.

To solve Problem (17), we apply the block coordinate descent method: given ${\bf{V}}$, update ${\bf{U}}$ and ${\bf{W}}$ by using (15); update $\{\alpha_{k}\}_{k\in\cal U}$ and ${\bf{V}}$ with given ${\bf{U}}$ and ${\bf{W}}$. We only need to solve the latter one. Putting the MSE expression in (13) into constraints C4 in Problem (17) yields

[TABLE]

where ${t_{k}}=\log\left|{{{\bf{W}}_{k}}}\right|+d-\sigma_{k}^{2}{\rm{Tr}}\left({{\bf{U}}_{k}^{H}{{\bf{U}}_{k}}{{\bf{W}}_{k}}}\right)$.

Without loss of generality, we assume ${\cal U}={\bar{\cal U}}=\{1,\cdots,K\}$ and define the indices of ${{{\cal U}_{i}}}$ as ${{\cal U}_{i}}=\{q_{1}^{i},\cdots,q_{|{{\cal U}_{i}}|}^{i}\}$. Problem (18) can be equivalently transformed into the following problem

[TABLE]

where ${{\bf{x}}_{k}}$ is given by

[TABLE]

and ${\bf{y}}_{i}$ is given by

[TABLE]

Problem (19) is an SOCP problem for which a globally optimal solution can be obtained by existing techniques such as interior point method [36].

Based on the above analysis, the iterative algorithm for solving Problem (11) is formally described in Algorithm 2.

**Theorem 1: **Algorithm 2 will converge during the iterative procedure.

**Proof: **Please see Appendix A. $\Box$

III-B Overall complexity to solve Problem (9) in Stage I

We first analyze the complexity of Algorithm 2 to solve Problem (11). For simplicity, we assume that candidate size for each user is equal, i.e., $\left|{{{\cal I}_{k}}}\right|=l$, and ${\cal U}=\bar{\cal U}$. In each iteration of Algorithm 2, the main complexity lies in step 2, where the SOCP Problem (19) is solved. This problem has $2MKld+K$ real variables, $K$ SOC constraints where each has $2Kd^{2}+1$ real dimensions, and $I$ SOC constraints where each has $2Md\left|{{{\cal U}_{i}}}\right|$ real dimensions. From [page 196, [40]], the complexity is $O\left({{{\left({2MKld+K}\right)}^{2}}\left({2{K^{2}}{d^{2}}+K+2Md\sum\nolimits_{i\in{\cal I}}{\left|{{{\cal U}_{i}}}\right|}}\right)}\right)$, and the total number of iterations required is $O\left({\sqrt{I+K}}\right)$. Note that $\sum\nolimits_{i\in{\cal I}}{\left|{{{\cal U}_{i}}}\right|}{\rm{=}}\sum\nolimits_{k\in{\cal U}}{\left|{{{\cal I}_{k}}}\right|}=Kl$, the total complexity to solve the SOCP Problem (19) is given by $O\left({\sqrt{I+K}{{\left({2MKld+K}\right)}^{2}}\left({2{K^{2}}{d^{2}}+K+2MdKl}\right)}\right).$ Finally, Algorithm 2 should be run at most $K$ times, then the overall complexity to solve Problem (9) in Stage I is at most ${T_{{\rm{StageI}}}}=O\left({K\sqrt{I+K}{{\left({2MKld+K}\right)}^{2}}\left({2{K^{2}}{d^{2}}+K+2MdKl}\right)}\right).$

IV Stage II: A Low-complexity Algorithm to Solve Problem (10)

In this section, we provide a low-complexity algorithm to solve Problem (10) with the selected users from Stage I. First, we adopt the re-weighted $l_{1}$-norm method [27] to transform the original non-smooth optimization problem into a series of WPM problems. Then, the WPM problem is solved by the WMMSE algorithm. In each iteration of the WMMSE algorithm, there is a subproblem that the precoder matrices should be optimized. We exploit the special structure of the subproblem and develop a low-complexity algorithm to solve it.

IV-A Reweighted $l_{1}$-norm minimization

For simplicity, the subscript in ${\cal U^{\star}}$ is omitted. It is easy to see that the minimum rate constraints in C1 of Problem (10) hold with equality at the optimal point, i.e., ${R_{k}}({\bf{V}})={R_{k,\min}},\forall k$. Then, defining $\tilde{P}_{i}^{c}\buildrel\Delta\over{=}{\rho_{i}}\sum\nolimits_{k\in{\cal U}_{i}^{\star}}{{R_{k,\min}}}+P_{i}^{c}$ and using the ${l_{0}}$-norm, the objective function of Problem (10) is equivalent to $\sum\nolimits_{i\in{\cal I}}{\left({{\eta_{i}}\sum\nolimits_{k\in{{\cal U}_{i}}}{\left\|{{{\bf{V}}_{i,k}}}\right\|_{F}^{2}+{{\left\|{\sum\nolimits_{k\in{{\cal U}_{i}}}{\left\|{{{\bf{V}}_{i,k}}}\right\|_{F}^{2}}}\right\|}_{0}}\tilde{P}_{i}^{c}}}\right)}.$ This rewritten expression enables us to apply the compressive sensing techniques [41], where the non-smooth ${l_{0}}$-norm objective can often be approximated by a re-weighted ${l_{1}}$-norm, i.e.,

[TABLE]

where $a_{i}^{(n)}$ is a weight factor of the $i$th RRH at the $n$th iteration that is iteratively updated as

[TABLE]

where $\delta$ is a small constant regularization parameter and ${{\bf{V}}_{i,k}^{(n)}}$ is the solution in the $n$th iteration. The above updating rule shows that those RRHs with lower transmit power in the previous iteration will have larger weights, which force them to be shut off in the future iterative procedure.

By using the approximation in (21), we have the following optimization problem that should be solved in the $n$-th iteration

[TABLE]

where $\omega_{i}^{(n-1)}={\eta_{i}}+{a_{i}^{(n-1)}}{\tilde{P}_{i}^{c}}.$

Based on the above analysis, the re-weighted $l_{1}$-norm based (RLN) algorithm to solve Problem (10) is given in Algorithm 3. The convergence of the RLN algorithm is proved in [28]. In addition, [28] showed that the RLN algorithm is guaranteed to achieve sparse solutions, while the other smooth approximations cannot produce sparse solutions in general.

IV-B Algorithm to Solve Problem (23)

For simplicity, the subscript of $\omega_{i}^{(n-1)}$ in Problem (23) is omitted. The main difficulty in solving Problem (23) lies in the rate requirement, which is non-convex. To handle this difficulty, we again apply the relationship between WMMSE and the rate expression. Based on Lemma 1, we replace the rate constraints in (23) with its lower bound ${h_{k}}\left({\bf{V}},{{\bf{U}}_{k},{\bf{W}}_{k}}\right)$. Define the indices of ${{\cal I}_{k}}$ as ${{\cal I}_{k}}=\left\{{s_{1}^{k},\cdots,s_{\left|{\cal I}_{k}\right|}^{k}}\right\}$, we have the following optimization problem

[TABLE]

where ${{\bf{G}}_{k}}$ and ${{\bf{B}}_{i,k}}$ are both diagonal matrices, given by

[TABLE]

and

[TABLE]

By solving Problem (24), we can find a solution that satisfies the KKT conditions of Problem (23). To solve it, we again apply the block coordinate descent method. Matrices ${\bf{U}}$ and ${\bf{W}}$ can be updated with (15). The remaining task is to update ${\bf{V}}$ with given ${\bf{U}}$ and ${\bf{W}}$. Plugging the MSE expression in (13) into the first set of Problem (24) yields

[TABLE]

where ${c_{k}}=\log\left|{{{\bf{W}}_{k}}}\right|+d-{R_{k,\min}}-{\rm{Tr}}\left({{{\bf{W}}_{k}}}\right)-\sigma_{k}^{2}{\rm{Tr}}\left({{\bf{U}}_{k}^{H}{{\bf{U}}_{k}}{{\bf{W}}_{k}}}\right)$.

Then, an WMMSE algorithm is proposed to solve Problem (23) in Algorithm 4. In the following theorem, we show the property of the WMMSE algorithm.

**Theorem 2: **The sequence of ${\bf{V}}$ generated by the WMMSE algorithm will converge to the KKT point of Problem (23).

**Proof: **Please see Appendix B. $\Box$

IV-C Low-complexity Algorithm to Solve Problem (27)

Since $\omega_{i}>0,\forall i$, matrices $\{{\bf{G}}_{k},\forall k\in\cal U\}$ are positive definite matrices. Then, Problem (27) can be similarly transformed an SOCP problem as in (19). Using the same method in Subsection III-B, the total complexity to solve this problem by using the interior point method is

[TABLE]

In the following, we go one step further to design an algorithm with lower complexity. Obviously, Problem (27) is a convex problem, and it satisfies the Slater’s condition [36]. Hence, the duality gap between Problem (27) and its dual problem is zero [36]. Then we can solve its dual problem instead of directly solving it.

With some simple manipulations, the Lagrangian function of Problem (27) is given by

[TABLE]

where ${\bm{\lambda}}=[{{\lambda_{k}},\forall k\in\cal U}]^{H}$ and ${\bm{\mu}}=[{{\mu_{i}},\forall i\in\cal I}]^{H}$ are the Lagrangian multipliers associated with the first and second sets of constrains of Problem (27), respectively, and ${\bf{\bar{G}}}_{k}$ is given by

[TABLE]

The dual function is given by

[TABLE]

Note that matrices $\{{\bf{G}}_{k},\forall k\in\cal U\}$ are positive definite matrices. Problem (29) is a convex problem, and the optimal solution can be obtained from its first-order derivative condition as:

[TABLE]

By inserting this solution into (29), the dual function becomes

[TABLE]

Hence, the dual problem of Problem (27) is given by

[TABLE]

where $f({\bm{\lambda}},{\bm{\mu}})=-g({\bm{\lambda}},{\bm{\mu}})$.

Fortunately, the objective function of the dual problem in (32) is differentiable and the dual problem is convex [36], the descent methods such as the gradient descent method and Newton’s method [36, 42] can be applied to solve it. In the following, we also utilize the block coordinate descent method to solve the dual problem (32): Optimize $\left\{{{\lambda_{k}}},\forall k\right\}$ with $\left\{{{\mu_{i}}},\forall i\right\}$, and vice versa.

Given $\left\{{{\mu_{i}}},\forall i\right\}$, Newton’s method is applied to find the optimal $\left\{{{\lambda_{k}}},\forall k\right\}$ of the dual problem, which is summarized in Algorithm 5. 333Since $\left\{{{\mu_{i}}},\forall i\right\}$ are given, $f({\bm{\lambda}})$ is short for $f({\bm{\lambda}},{\bm{\mu}})$ and the same for $f({\bm{\mu}})$ later.

In step 4 of Algorithm 5, the backtracking line search method is used to find the step size, where $\xi$ is typically chosen as a very small value and $\varphi$ is chosen between [math] and $1$. The step $\kappa^{(t)}$ starts with one and then reduces by a factor of $\varphi$ until the stop condition (33) is satisfied. Note that in each iteration of Algorithm 5, the step value $\kappa^{(t)}$ may be different. The constant $\xi$ can be regarded as the acceptable fraction of the decrease in the objective value of $f$ that is predicted by the line search method.

However, to make this algorithm work, there are still problems to be solved: how to calculate the gradient and how to compute the Hessian matrix. To derive the expressions of the gradient and the Hessian matrix, we first introduce some useful results in the matrix differential calculus. Given a matrix function $\bm{\Gamma}(x)$, one has [43, 44]

[TABLE]

[TABLE]

In addition, to simplify the expressions of the gradient and the Hessian matrix, one defines some matrices:

[TABLE]

Based on the above results and definitions, the gradient can be derived as follows:

[TABLE]

with

[TABLE]

The Hessian matrix of $f(\bm{\lambda})$ can be calculated as:

[TABLE]

Next, given $\left\{{{\lambda_{k}}},\forall k\in\cal U\right\}$, we solve the dual problem (32) to update $\left\{{{\mu_{i}}},\forall i\in\cal I\right\}$. Here, the gradient descent method [36] is applied. Although Newton’s method converges faster than the gradient descent method, simulation results show that the gradient method also converges within five iterations but it has much lower computational complexity than Newton’s method since it does not require the calculations of the Hessian matrix and the inverse of the Hessian matrix. The gradient descent method to update $\left\{{{\mu_{i}}},\forall i\in\cal I\right\}$ is given in Algorithm 6.

In Algorithm 6, the gradient $\nabla f(\bm{\mu})$ is required. Define ${{\bf{D}}_{k}}={{\bf{C}}_{k}}{\bf{C}}_{k}^{\rm{H}}$. Then, by using the results in (34) and (35), the gradient $\nabla f(\bm{\mu})$ can be calculated as

[TABLE]

with

[TABLE]

Finally, based on the above analysis, the method to solve the dual problem (32) is given in Algorithm 7, which is named as Block Coordinate Descent (BCD) method.

**Theorem 3: **The sequences of $\bm{\mu}$ and $\bm{\lambda}$ generated by Algorithm 7 will converge to the globally optimal solution of the dual problem (32).

**Proof: **Please see Appendix C. $\Box$

When the optimal $\bm{\mu}$ and $\bm{\lambda}$ are obtained by using Algorithm 7, the optimal solution to Problem (27) is given by (30). As there is zero duality gap between the primal problem (27) and dual problem (32), which means that this solution is the globally optimal solution of Problem (27).

IV-D Overall Complexity to Solve Problem (27) in Stage II

In this subsection, we analyze the overall complexity to solve Problem (27). It mainly includes three layers of iterations: the first layer is the RLN algorithm to deal with the non-smooth $l_{0}$ norm, the second layer is the WMMSE algorithm to deal with the non-convex rate constraints, and the third layer is the BCD algorithm to solve Problem (27).

We first analyze the complexity of the third layer for BCD algorithm. Note that ${{{\bf{\tilde{H}}}}_{j,k}}$, {{{\bf{\mathord{\buildrel\lower 3.0pt\hbox{\scriptscriptstyle\smile}\over{H}}}}}_{j,k}}, and ${{{\bf{\hat{H}}}}_{j,k}}$ can be calculated before the iterations of the BCD Algorithm. The main complexity of the BCD Algorithm lies in step 2 and step 3, where Newton’s method and gradient descent method are used to update ${\bm{\lambda}}$ and $\bm{\mu}$, respectively.

We first analyze the computational complexity of Newton’s method under the same assumption in Subsection III-B. The main complexity in each iteration of Newton’s method lies in step 2 and step 4 of Algorithm 5. We first analyze step 2 of Algorithm 5. According to [45], the complexity of calculating $\{{{\bf{\tilde{G}}}}_{k},\forall k\in\cal U\}$ is on the order of $O\left({{{K\left({Ml}\right)}^{2.376}}}\right)$. For any two matrices ${\bf{X}}\in{{\mathbb{C}}^{m\times n}},{\bf{Y}}\in{{\mathbb{C}}^{n\times p}}$, the complexity of computing ${\bf{XY}}$ is on the order of $O\left({mnp}\right)$ [36]. In general, $d\ll MI$. Then, the total complexity of computing $\left\{{{{\bf{C}}_{k}},{{\bf{F}}_{k}},\forall k\in\cal U}\right\}$ is on the order of $O\left({K{M^{2}}{l^{2}}d}\right)$. Similarly, the total complexity of computing $\left\{{{{\bf{Y}}_{j,k}},{{{\bf{\tilde{Y}}}}_{j,k}},{{\bf{Z}}_{j,k}}}\forall j,k\in\cal U\right\}$ is on the order of $O\left({{K^{2}}{M^{2}}{l^{2}}d}\right)$. Hence, the total complexity of computing $\left\{{{{\bf{C}}_{k}},{{\bf{F}}_{k}}},\forall k\right\}$ and $\left\{{{{\bf{Y}}_{j,k}},{{{\bf{\tilde{Y}}}}_{j,k}},{{\bf{Z}}_{j,k}},\forall j,k\in\cal U}\right\}$ is on the order of $O\left({{K^{2}}{M^{2}}{l^{2}}d}\right)$. With a similar analysis, the total complexity of computing (39) is on the order of $O\left({{K^{3}}Ml{d^{2}}}\right)$. In addition, the complexity of computing the inverse of $\nabla^{2}f(\bm{\lambda}^{(t)})$ is on the order of $O\left({{K^{2.376}}}\right)$ [45]. Hence, the total complexity of step 2 of Algorithm 5 is on the order of $O\left({\max\left\{{{K^{3}}Ml{d^{2}},{{\left(K{Ml}\right)}^{2.376}},{K^{2}}{{\left({Ml}\right)}^{2}}d}\right\}}\right)$. In the $t$th iteration of step 4 of Algorithm 5, $f({{\bm{\lambda}}^{(t+1)}})$ is required to calculate $m^{(t)}$ times. The complexity in each time is on the order of $O\left({\max\left\{{{{\left({Ml}\right)}^{2.376}},K{{\left({Ml}\right)}^{2}}d}\right\}}\right)$. Thus, the total complexity of step 4 of Algorithm 5 is on the order of $O\left({m^{(t)}}{\max\left\{{{{\left({Ml}\right)}^{2.376}},K{{\left({Ml}\right)}^{2}}d}\right\}}\right)$. Simulation results show that in general $m^{(t)}$ is always equal to one, which means that $f({{\bm{\lambda}}^{(t+1)}})$ only needs to be computed for once. Hence, the complexity of step 4 of Algorithm 5 can be approximately by $O\left({\max\left\{{{{\left({Ml}\right)}^{2.376}},K{{\left({Ml}\right)}^{2}}d}\right\}}\right)$. As a result, the total complexity of Newton’s method is

[TABLE]

Simulation results show that Newton’s method converges very rapidly and in general five iterations are enough for the algorithm to converge.

By using the similarly analytical technique to Newton’s method, the total complexity of the gradient descent method is given by

[TABLE]

The simulation results in the next section show that the gradient descent method usually converges within five iterations. Hence, in each iteration of the BCD Algorithm, the complexity of Newton’s method dominates the complexity of the gradient descent method.

Based on the above analysis, the overall complexity to solve Problem (27) in Stage II is

[TABLE]

where $t_{{\rm{RLN}}}$, $t_{{\rm{WMMSE}}}$ and $t_{{\rm{BCD}}}$ represent the average number of iterations required by the RLN, WMMSE, and BCD algorithms, respectively. Simulation results show that these three algorithms converge very fast and generally five iterations are enough to achieve large portion of the final performance.

V Simulation Results

In this section, we present simulation results to evaluate the performance of the proposed algorithms. To be more realistic, we consider a wrap-around system model shown in Fig. 3 as in [46], where the C-RAN network is deployed in the central square with $[-1000\ {\rm{1000}}]\times[-1000\ {\rm{1000}}]$ meters, surrounded by eight uncoordinated square macrocells. It is assumed that all the users and RRHs are uniformly and independently distributed in the C-RAN region. We adopt the channel model that consists of four parts: 1) the long term evolution (LTE) standard path loss model: $PL_{i,k}=148.1+37.6{\log_{10}}d_{i,k}\ ({\rm{dB}})$, where $d_{i,k}$ (in km) is the distance from the $i$th RRH to the $k$th user; 2) Log-normal shadowing with zero mean and 8 dB standard derivation; 3) Rayleigh fading with zero mean and unit variance; 4) transmit antenna power gain of 9 dBi. Each user is assumed to have the same rate requirement, i.e., $R_{\rm{min}}=R_{k,\rm{min}},\forall k$, and each RRH has the same power constraint, i.e., $P_{\rm{max}}=P_{i,\rm{max}}=4{\rm{W}},\forall i\in\cal I$. It is assumed that each user is potentially served by its nearest $X$ RRHs, i.e., $|{{\cal I}_{k}}|=X,\forall i$. Unless stated otherwise, the system parameters are set as follows: error tolerance is $\varepsilon=10^{-3}$, thermal noise power is $\sigma^{2}=-104\ {\rm{dBm}}$, $I=12$, $K=8$, $X=3$, $M=2$, $N=2$, $d={\rm{min}}\{M,N\}$, $\eta_{i}=4$ [47], $\rho_{i}=0.5$ [26], $P_{i}^{a,{\rm{rrh}}}=3.4{\rm{W}}$, $P_{i}^{s,{\rm{rrh}}}=2.15{\rm{W}}$, $P_{i}^{a,{\rm{fr}}}=3.85{\rm{W}}$, $P_{i}^{s,{\rm{fr}}}=0.75{\rm{W}}$, ${P_{{\rm{BBU}}}}=20W$ [48, 22]. Moreover, let $\cal L$ be the set of uncoordinated base stations (BSs) in C-RAN’s nearby eight macrocells. The noise power at user $k$ can be modeled as $\sigma_{k}^{2}={\sigma^{2}}+\sum\nolimits_{m\in{\cal L}}{{P_{\max}}P{L_{m,k}}{S_{m,k}}{G_{m}}}$ [46], where $P{L_{m,k}}$ and $S_{m,k}$ are the large-scale fading and shadowing respectively from the BS in macrocell $m$ to user $k$, ${G_{m}}$ represents the antenna gain.

V-A Properties of the Proposed Algorithms

V-A1 Convergence behavior of Algorithm 2

Fig. 3 shows the convergence behaviour of Algorithm 2 for different numbers of receive antennas. The results are obtained by averaging over 100 channel realizations. Due to the non-convexity of Problem (11), different initial points for Algorithm 2 may yield different solutions. To investigate this effect, we consider two initialization schemes: 1) SVD-initial, in which the beam directions for each user are chosen as the unitary matrices obtained by the singular value decomposition (SVD) of channel matrices and the total power at each RRH is equally allocated to the users potentially served by each RRH; 2) Rand-initial, in which both the beam directions and power allocations are randomly generated. It can be seen from Fig. 3 that the objective value of Problem (11) monotonically decreases during the iterative procedure for two initialization schemes. In addition, the algorithm converges very fast and in general six iterations are sufficient to achieve a large proportion of the converged value for different numbers of receive antennas and different initialization schemes. It is interesting to find that the algorithm under two different initialization schemes will converge to almost the same value. As expected, the converged objective value decreases with the number of receive antennas since more degrees of freedom are available.

V-A2 User selection performance of USC algorithm

Fig. 5 compares the performance of the USC algorithm with two algorithms: greedy search method and exhaustive search method. For the greedy search method, in each time we compute the objective value of Problem (11) when excluding one user, then the user yielding the smallest objective value will be removed. This procedure continues until all remaining users are feasible. Note that this algorithm increases quadratically with $K$. The exhaustive search method checks all feasible sets of users and chooses the largest one. Its complexity increases exponentially with $K$. As expected, the number of admitted users decreases with the rate requirements for all algorithms. The greedy search method achieves almost the same performance as the exhaustive one, and the performance gap between the exhaustive search algorithm and the proposed USC algorithm can be negligible. However, the complexity of our proposed USC algorithm only increases linearly with $K$. The impact of initial points is also studied and we find that both initialization schemes (SVD-initial and rand-initial) have similar performance, which is not shown here for clarity.

V-A3 Convergence behaviour of the RLN algorithm

The convergence behaviours of the RLN algorithm are shown in Figs. 5 (a) and (b) for the NPC and the number of the remaining RRHs in each iteration, respectively. Three different values of $\delta$ are tested, i.e., $\delta={10^{-4}},{10^{-5}}$ and ${10^{-6}}$. One randomly generated channel is used to obtain the convergence behaviour, where the USC algorithm is first executed to find the largest feasible set of users. In this example, User 8 is removed to guarantee the feasibility of the other users as seen in Fig. 7. It can be seen from the figures that for all values of $\delta$, both the number of active RRHs and the NPC decrease rapidly and there is no additional decrease after the fifth iteration. At the converged state, only six RRHs are active. Compared to the full cooperation strategy where all RRHs are active, we can save large amount of power as seen from Fig. 5 (a). Fig. 7 illustrates the converged state of the system. It can be seen that RRH 2 is switched off since it is far from the users and User 8 is not selected as it is far from the RRHs. We also study the impact of initialization schemes on the performance of the RLN algorithm. The initial precoders for the RLN algorithm are the outputs of the USC algorithm which is initialized with the SVD-initial and rand-initial schemes. The simulation results show they achieve almost the same performance, which is not shown here for clarity.

V-A4 Convergence behaviour of the WMMSE algorithm

In step 2 of each iteration of the RLN algorithm, we need to solve Problem (23) by using the WMMSE algorithm. Fig. 7 shows the convergence performance of the WMMSE algorithm for the first three iterations of the RLN algorithm. It is observed that the WMMSE algorithm converges within ten iterations for the first iteration of RLN algorithm. However, the objective values stay almost fixed for the second and third iterations of the RLN algorithm. This means that only in the first iteration of RLN algorithm, some iterations are required for the WMMSE algorithm.

V-A5 Convergence behaviour of the BCD algorithm

In step 2 of each iteration of the WMMSE algorithm, Problem (27) should be solved to update the precoding matrices by using the BCD algorithm. Fig. 9 shows the convergence behaviour of the BCD algorithm for the first iteration of the WMMSE algorithm. It is seen that the algorithm converges very fast and one iteration is sufficient to achieve a large portion of the converged value (99.2% in this example).

V-A6 Convergence behaviour of Newton’s method and the gradient descent method

In each iteration of the BCD algorithm, Newton’s method is required to update $\left\{{\lambda_{k}},\forall k\right\}$ and the gradient descent method is applied to update $\left\{{\mu_{i}},\forall i\right\}$. The convergence behaviours of these two algorithms for the first three iterations of the BCD algorithm are shown in Figs. 9 (a) and (b), respectively. Newton’s method requires several iterations to converge only in the first iteration of the BCD algorithm, while stays almost constant for the second and third iterations of the BCD algorithm. Interestingly, the gradient descent method only requires one iteration to converge in the first iteration of the BCD algorithm and keeps fixed during the rest of the iterations of the BCD algorithm. By combining the complexity analysis in (42), (43) and the above convergence behaviours, we can conclude that the BCD algorithm has a much lower computational complexity than directly solving the SOCP problem.

V-A7 Impacts of the number of data streams

In Fig. 11, the impact of the number of data streams on the number of admitted users is studied. As expected, the number of admitted users decreases with the rate requirements and larger number of data streams can support more users. We find significant performance gains can be achieved when the number of data streams increases from 1 to 2, especially for the high data rate requirements. However, only marginal performance gains are achieved by the case of $d=4$ over the case of $d=2$, which comes at the higher cost of computational complexity. This reveals that the performance saturates with the increase of data streams $d$. In Fig. 11, the impacts of data streams on the NPC and on the number of active RRHs are studied with the same setup in Fig. 11. For fair comparison, we only consider the set of users that can be supported under the case of $d=1$, so that all cases can support the selected users. Fig. 11 (a) shows that the NPC first increases with the rate requirements when $R_{\rm{min}}\leq 3\ {\rm{nats/s/Hz}}$ and then decreases significantly when $R_{\rm{min}}>3\ {\rm{nats/s/Hz}}$. The reason can be explained as follows. When $R_{\rm{min}}$ increases from $1$ to $3\ {\rm{nats/s/Hz}}$, the number of admitted users almost keeps stable as shown in Fig. 11, while the fronthaul power increases when the rate requirement increases and the number of active RRHs increases to support the higher rate requirements as seen in Fig. 11 (b), which in turn consumes more power consumption. On the other hand, when $R_{\rm{min}}$ increase from $3$ to $6\ {\rm{nats/s/Hz}}$, the number of admitted users decreases dramatically as shown in Fig. 11, which leads to reduced transmit power and a reduced number of active RRHs as shown in Fig. 11 (b). Again, it is observed from Fig. 11 (a) that a greater number of data streams requires lower NPC, but the performance gain shrinks with the number of data streams.

V-A8 Impacts of the number of transmit antennas

In Fig. 13, the impact of the number of transmit antennas on the number of admitted users is studied. As expected, the number of admitted users increases with the number of transmit antennas due to more degrees of freedom. Significant performance gains can be achieved by the case of $M=2$ over the case of $M=1$, especially in the high rate regime. However, the performance gain shrinks for the case of $M=4$ over the case of $M=2$. In Fig. 13, the impacts of the number of antennas on the NPC and the number of active RRHs are investigated. For fair comparison, it is also assumed that the set of users selected from the case of $M=1$ are the input of Stage II for all cases of different values of $M$ so that the selected users are the same and feasible for all cases. It is interesting to find that when $M$ increases, the NPC increases while the number of active RRHs decreases. This is mainly due to the fact that the RRH power consumption model in (4) increases linearly with $M$, and this increased power consumption dominates the reduced power consumption resulting from the reduced number of active RRHs. It should be emphasized that in some other cases with different values of system parameters, the NPC may not increase with $M$ and the counter part happens, such as the case of the low circuit power consumption for each antenna and high power consumption associated with the fronthaul power consumption.

V-A9 Impacts of the candidate size

The impact of candidate size on the number of admitted users is illustrated in Fig. 15 for a dense network with 20 RRHs and 12 users. As expected, larger candidate sizes can support more users due to the increased degrees of freedom. However, the performance gains decreases with the candidate sizes, which implies that there is no need to consider distant RRHs for each user since they contribute less to their signal strength. In general, the candidate size should be no larger than 4 to achieve a good tradeoff between performance and complexity. Similarly to the trend observed in Fig. 11 (a), it is seen from Fig. 15 (a) that the NPC increases in the low rate regime, while decreasing significantly in the high rate regime. For the former part, the reason is that the increased fronthaul power dominates the reduced circuit power for the reduced active RRHs. While for the latter part, the reason is the opposite. Also, it is observed that the NPC performance gain for larger candidate size is more obvious in the low rate regime, while the performance is almost the same in the high rate regime. This is mainly due to the fact that in the high rate regime, only a small number of users can be admitted, and these users are separated far away. As a result, the multiuser interference is not so significant and each user’s nearest RRH is able to serve it with the rate requirement.

V-B Performance comparison

We compare the performance of the RLN algorithm with the following RRH selection methods:

• •

Exhaustive search (Exhau-search) method: For each given active RRH set ${\cal A}$, this method first checks its feasibility. If feasible, the method will use the WMMSE algorithm to solve the corresponding transmit power minimization problem. The complexity of this method increase exponentially with $I$, which is served as the performance benchmark for our proposed algorithm.

• •

Successive RRH selection (Succesive-sel) method: This method first lets all the RRHs be active and check its feasibility. If feasible, the method applies the WMMSE algorithm to solve the transmit power minimization problem. Then, the method gradually removes the RRHs according to their transmit power from the lowest to the highest until the problem becomes infeasible. The complexity of this scheme increases linearly with $I$.

• •

Greedy search method: In each step, we exclude each RRH and calculate the NPC when the remaining RRHs are active. Then, we remove the RRH so that the remaining RRHs yield the least NPC. This procedure terminates until the problem becomes infeasible. The complexity of this scheme increases quadratically with $I$.

• •

Full cooperative (Full-coop) method: In this method, all the selected RRHs in cluster-formation stage are active and the WMMSE algorithm is used to solve the transmit power minimization problem.

For fair comparison, we assume in the following simulation results, only the channel realizations that are feasible for all users are considered.

V-B1 Impact of the rate requirements

Figs. 17 (a) and (b) illustrate the average NPC and the corresponding number of active RRHs versus the rate requirements, respectively. Fig. 17 (a) shows that the RLN algorithm outperforms the ‘Succesive-sel’ method and ‘Full-coop’ method for all rate regimes. However, the performance of the ‘Greedy search’ method is slightly better than the RLN algorithm when $R_{\rm{min}}\leq 3{\rm{nats/s/Hz}}$, while the RLN algorithm outperforms the ‘Greedy search’ method in the high rate regime and the performance gain increases with the rate requirements. Fig. 17 (b) shows a similar trend in terms of the number of active RRHs. Compared with the optimal ‘Exhau-search’ method, the performance loss in power consumption is at most $8\%$ when $R_{\rm{min}}=1$ nats/s/Hz, and this gap gradually diminishes with the increase of rate requirements. In particular, the performance gain provided by the ‘Exhau-search’ method over the RLN algorithm is negligible when $R_{\rm{min}}=5$ nats/s/Hz. As expected, the ‘Full-coop’ method consumes the highest power since all selected RRHs are active.

V-B2 Impact of the number of RRHs

Figs. 17 (a) and (b) illustrate the average NPC and the corresponding number of active RRHs versus the total number of RRHs, respectively. It is seen that the NPC achieved by all schemes decreases with $I$ due to the fact that when there are more RRHs, the average access distance between users and RRHs decreases significantly and thus leads to more reduced transmit power. It is again observed that the performance of the RLN algorithm is superior to that of the ‘Succesive-sel’ method. This implies that selecting the RRHs only based on the transmit power is not enough, and may incur significant performance loss. Note that the ‘Greedy search’ method requires higher power consumption than the RLN algorithm for all numbers of RRHs, especially when $I=6$. Also, the performance of ‘Exhau-search’ method is slightly better than the RLN algorithm. Note that although the number of active RRHs increases slightly with the total number of RRHs as seen in Figs. 17 (b), the NPC decreases. This may due to the fact that the overall transmit power reduction overwhelms the increase of circuit power.

VI Conclusion

In this paper, a joint selection of active RRHs and optimization of the precoding matrices which minimizes the NPC for the MIMO C-RAN, while guaranteeing users’ rate requirements and per-RRH power constraints, has been studied. A low-complexity user selection was proposed to guarantee the feasibility of the other users. Then a low-complexity iterative algorithm, based on the reweighted $l_{1}$-norm minimization method, WMMSE algorithm, Newton’s method, and gradient descent method, was proposed to solve the network power minimization problem for the selected users. Simulation results show that the proposed algorithms converge fast, which is attractive for practical implementation. Also, more antennas at the user side can admit more users. The proposed user selection algorithm was shown to achieve the similar performance as the optimal exhaustive search method. Moreover, our proposed algorithm was shown to achieve much greater power savings than the full cooperation method, and the performance loss compared with the optimal approach is insignificant.

Appendix A Proof of Theorem 1

In step 2 of the $n$th iteration, we solve Problem (19) to obtain the optimal $\{\alpha_{k}^{(n)}\}_{k\in\cal U}$ and ${\bf{V}}^{(n)}$ with given ${\bf{U}}^{(n-1)}$ and ${\bf{W}}^{(n-1)}$. Hence, we have ${h_{k}}\left({{\bf{V}}^{(n)},{{\bf{U}}_{k}^{(n-1)}},{{\bf{W}}_{k}^{(n-1)}}}\right)\geq\left({\alpha_{k}^{(n)}}\right)^{2}{R_{k,{\rm{min}}}},\forall k$. In step 3 of the $n$th iteration, we update ${\bf{U}}^{(n)}$ and ${\bf{W}}^{(n)}$ as in (15) with ${\bf{V}}^{(n)}$. According to Lemma 1, we have $R_{k}\left({{{\bf{V}}^{(n)}}}\right)={h_{k}}\left({{{\bf{V}}^{(n)}},{\bf{U}}_{k}^{(n)},{\bf{W}}_{k}^{(n)}}\right)\geq{h_{k}}\left({{{\bf{V}}^{(n)}},{\bf{U}}_{k}^{(n-1)},{\bf{W}}_{k}^{(n-1)}}\right)$. Hence, we have

[TABLE]

In step 2 of the $(n+1)$th iteration, we obtain $\{\alpha_{k}^{(n+1)}\}_{k\in\cal U}$ and ${\bf{V}}^{(n+1)}$ with given ${\bf{U}}^{(n)}$ and ${\bf{W}}^{(n)}$ by solving Problem (19). Then we have $\sum\nolimits_{k\in{\cal U}}{{{\left({\alpha_{k}^{(n+1)}-1}\right)}^{2}}}\leq\sum\nolimits_{k\in{\cal U}}{{{\left({\alpha_{k}^{(n)}-1}\right)}^{2}}}$. The reason is that from (A.1), $\{\alpha_{k}^{(n)}\}_{k\in\cal U}$ and ${\bf{V}}^{(n)}$ is just a feasible solution for Problem (19) with given ${\bf{U}}^{(n)}$ and ${\bf{W}}^{(n)}$. Hence, the objective value of Problem (11) is monotonically decreasing. Obviously, the objective value is lower bounded by zero. Hence, Algorithm 2 will converge.

Appendix B Proof of Theorem 2

We first prove that the sequence of ${\bf{V}}$ generated by the WMMSE algorithm (i.e. Algorithm 4) always satisfies the rate requirements of Problem (23). In step 2, we obtain ${\bf{V}}^{(l)}$ with ${\bf{U}}^{(l-1)}$ and ${\bf{W}}^{(l-1)}$. Hence, ${h_{k}}\left({{\bf{V}}^{(l)},{\bf{U}}_{k}^{(l-1)},{\bf{W}}_{k}^{(l-1)}}\right)\geq{R_{k,{\text{min}}}},\forall k$ hold since ${\bf{V}}^{(l)}$ is feasible for Problem (24). According to Lemma 1, ${h_{k}}\left({{\bf{V}}^{(l)},{\bf{U}}_{k}^{(l-1)},{\bf{W}}_{k}^{(l-1)}}\right)$ is a lower-bound of ${R_{k}}({\bf{V}}^{(l)})$, i.e., ${R_{k}}({\bf{V}}^{(l)})\geq{h_{k}}\left({{\bf{V}}^{(l)},{\bf{U}}_{k}^{(l-1)},{\bf{W}}_{k}^{(l-1)}}\right)$. Hence, ${R_{k}}({\bf{V}}^{(l)})\geq R_{k,{\text{min}}}$ holds. Thus, the sequence of ${\bf{V}}$ generated by the WMMSE algorithm satisfies the rate requirements of Problem (23).

Next, we show that the value of the objective function of Problem (23) monotonically decreases during the iterative process of the WMMSE algorithm. Denote ${\rm{Obj}}({\bf{V}}^{(l)})$ as the objective value of Problem (23) when ${\bf{V}}={\bf{V}}^{(l)}$. Step 2 of the WMMSE algorithm updates ${\bf{V}}^{(l)}$ by solving Problem (27) with ${\bf{U}}^{(l-1)}$ and ${\bf{W}}^{(l-1)}$. The objective value of this step, ${\rm{Obj}}({\bf{V}}^{(l)})$, will be no larger than ${\rm{Obj}}({\bf{V}}^{(l-1)})$, i.e., ${\rm{Obj}}({\bf{V}}^{(l)})\leq{\rm{Obj}}({\bf{V}}^{(l-1)})$. The reason is that ${\bf{V}}^{(l-1)}$ is a feasible solution for Problem (27) with ${\bf{U}}_{k}^{(l-1)}$ and ${\bf{W}}_{k}^{(l-1)}$ since ${h_{k}}\left({{\bf{V}}^{(l-1)},{\bf{U}}_{k}^{(l-1)},{\bf{W}}_{k}^{(l-1)}}\right)={R_{k}}({\bf{V}}^{(l-1)})\geq{R_{k,{\text{min}}}}$ holds as proved above. In step 3 of the WMMSE algorithm, we update ${\bf{U}}^{(l)}$ and ${\bf{W}}^{(l)}$ by using (15) with ${\bf{V}}^{(l)}$. This step increases the value of ${h_{k}}\left({{\bf{V}},{\bf{U}}_{k},{\bf{W}}_{k}}\right)$ while maintaining the same objective value of Problem (23). Therefore, this step provides “room” for the next iteration to decrease the objective value. In addition, the objective value is lower bounded by zero. Hence, the WMMSE algorithm converges.

Then, we prove that given the initial set of precoders, the WMMSE algorithm converges to a unique solution. Obviously, when ${\bf{V}}$ is given, ${\bf{U}}$ and ${\bf{W}}$ can be uniquely determined by (15). The remaining task is to prove that given ${\bf{U}}$ and ${\bf{W}}$, the BCD algorithm can obtain the unique globally optimal solution ${\bf{V}}$. Since $\{{\bf{G}}_{k},\forall k\}$ are positive definite matrices, the objective function in Problem (27) is a strictly convex function with respect to (w.r.t.) $\bf{V}$. Obviously, the constraints in Problem (27) are convex w.r.t. $\bf{V}$ [36]. Hence, Problem (27) is a strictly convex problem [36]. According to [Page 137 in [36]], the globally optimal solution of Problem (27) is unique. On the other hand, Theorem 3 proves that the BCD algorithm can obtain the globally optimal solution to the dual problem (32). As Problem (27) is a convex problem and it satisfies the Slater’s condition [36], the duality gap between Problem (27) and its dual problem (32) is zero [36]. As a result, the BCD algorithm can obtain the unique globally optimal solution ${\bf{V}}$. Finally, by alternatively updating step 2 and step 3, the WMMSE algorithm will converge to a unique solution. It should be emphasized that as Problem (24) is non-convex, it may have many locally optimal solutions, and the unique solution of the WMMSE algorithm depends on the initial point. However, given the initial points of precoders, the WMMSE algorithm will converge to a unique solution.

Finally, we prove that the unique solution satisfies the KKT conditions of Problem (23). Denote the converged solution of the WMMSE algorithm as ${\bf{V}}^{\star}$, ${\bf{U}}^{\star}$ and ${\bf{W}}^{\star}$. With given ${\bf{U}}^{\star}$ and ${\bf{W}}^{\star}$, the Lagrange function of Problem (24) can be written as

[TABLE]

where ${{\bm{\lambda}}}=\left\{{\lambda_{k},\forall k\in\cal U}\right\}$ and ${\bm{\mu}}=\left\{{\mu_{i},\forall i\in\cal I}\right\}$ are the corresponding Lagrange multipliers.

According to Theorem 3, the BCD algorithm can obtain the globally optimal solution of Problem (27) (also Problem (24)) with given ${\bf{U}}^{\star}$ and ${\bf{W}}^{\star}$, there must exist ${{\bm{\lambda}}^{\star}}$ and ${\bm{\mu}^{\star}}$ such that $\{{\bf{V}}^{\star},{\bm{\lambda}}^{\star},\bm{\mu}^{\star}\}$ satisfy the following KKT conditions

[TABLE]

Since ${\bf{U}}^{\star}$ and ${\bf{W}}^{\star}$ are updated by using (15), we have ${h_{k}}\left({{\bf{V}}^{\star},{\bf{U}}_{k}^{\star},{\bf{W}}_{k}^{\star}}\right)={R_{k}}({\bf{V}}^{\star})$ according to Lemma 1. By substituting it into the equations (B.6), (B.6) and (B.6), we find that the set of equations (B.6)-(B.6) are just the KKT conditions of Problem (23).

Appendix C Proof of Theorem 3

According to [36], the dual problem of any optimization problem is a convex problem. Thus, the dual problem (32) is jointly convex with respect to ${\bm{\lambda}}$ and $\bm{\mu}$. Assuming that the constraint of this problem satisfies the Slater’s condition, the KKT condition of this problem is sufficient and necessary for optimality. For given $\bm{\mu}$, the dual problem (32) is a convex problem w.r.t. ${\bm{\lambda}}$. According to [36], Newton’s method can obtain the globally optimal solution of dual problem (32) for given $\bm{\mu}$. In addition, for given ${\bm{\lambda}}$, the dual problem (32) is convex w.r.t. $\bm{\mu}$, and the gradient descent method can be applied to obtain the globally optimal solution. Then by adopting the same idea as in the proof of Theorem 1 in [49], we can prove that the converged solution also satisfies the KKT condition of Problem (32). Since Problem (32) is a convex optimization problem, Algorithm 6 can attain the globally optimal solution of Problem (32).