{'fraction_non_alphanumeric': 0.03728611949858175, 'fraction_numerical': 0.012260957086650197, 'mean_word_length': 4.463942007249094, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 28, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Energy Efficient Coordinated Beamforming for Multicell System: Duality-Based Algorithm Design and Massive MIMO Transition

In this paper, we investigate joint beamforming and power allocation in multicell multiple-input single-output (MISO) downlink networks. Our goal is to maximize the utility function defined as the ratio between the system weighted sum rate and the total power consumption subject to the users' quality of service requirements and per-base-station (BS) power constraints. The considered problem is nonconvex and its objective is in a fractional form. To circumvent this problem, we first resort to an virtual uplink formulations of the the primal problem by introducing an auxiliary variable and applying the uplink-downlink duality theory. By exploiting the analytic structure of the optimal beamformers in the dual uplink problem, an efficient algorithm is then developed to solve the considered problem. Furthermore, to reduce further the exchange overhead between coordinated BSs in a large-scale antenna system, an effective coordinated power allocation solution only based on statistical channel state information is reached by deriving the asymptotic optimization problem, which is used to obtain the power allocation in a long-term timescale. Numerical results validate the effectiveness of our proposed schemes and show that both the spectral efficiency and the energy efficiency can be simultaneously improved over traditional downlink coordinated schemes, especially in the middle-high transmit power region.


I. INTRODUCTION
M ULTIPLE-INPUT multiple-output (MIMO) technologies have attracted much research interest in the past decade due to the potential of significantly enhancing the link reliability and improving the spectral efficiency (SE) without requiring additional transmit power and bandwidth resources [1]. Among those, multiuser MIMO and coordinated multipoint transmission and reception (CoMP) have been intensively studied for advanced cellular systems [2]. The cooperation between the base stations (BSs) offers great advantages in suppressing the inter-cell interference and enabling efficient resource allocation to significantly increase the SE by exploiting the spatial diversity and multiplexing gain [3]- [5].
Recently, besides the objective of improving the SE, energy efficiency (EE) has become a main focus in mobile communications, due to the explosive growth of energy consumption resulting from the exponential increase in mobile multimedia data traffic and mobile terminals [6], [7]. Motivated by these, advanced technologies aiming at much improved EE have emerged as an important topic both in industry and academia [8].
EE is usually defined as the ratio of the system weighted sum rate (WSR) to the system total power consumption. Recently energy-efficient resource allocation has been widely studied for orthogonal frequency division multiple access (OFDMA) downlink networks with single antenna transceiver or fixed transmit beamformers [9]- [13]. On the other hand, there have been tremendous research efforts to improve the system EE by using multiple antenna technologies [14]- [19]. In [14] the authors studied the power allocation problem in the uplink of multiuser MIMO systems and showed that the EE is maximized by adaptively switching off some of the user antennas. The EE of multiuser MIMO downlink systems was also studied in [15]- [19]. In particular, in [15], the energy-efficient multiuser beamforming optimization problem was solved by approximating the objective with a convex and tight upper bound function. In [16] a zero-gradient-based energy-efficient iterative approach was developed which was guaranteed to converge to a local maximum for MIMO interference channel. The authors in [17] addressed the EE optimization problem by solving three subproblems and choosing the best one. Focusing on the downlink of multicell multiuser systems, we addressed in [18], [19] the energy-efficient beamforming problem with per-BS power constraints by using jointly fractional programming [20]- [22] and the relationship between the user rate and the minimum mean square error (MMSE). The research results in the above literature have shown that in many cases, the system EE is improved at the cost of decreasing system SE. This means that it is crucial to achieve a reasonable tradeoff between EE and SE in the design of a wireless communications system.
Recently, massive MIMO technology has emerged as a promising solution to improve both SE and EE and is well analyzed especially for the uncorrelated channel model [23]- [25]. Particularly its achievable performance in metrics such as the SINR and sum rate (SR) has been derived with closedform expressions for some special cases. For instance, the deterministic quantity of the SINR for zero-forcing (ZF) and regularized ZF (RZF) precoding were derived under channel uncertainty for the multiple-input single-output (MISO) broadcast system [26], showing that the performance of massive MIMO systems mainly depends on the channel statistics and the channel small fading has little effect. By exploiting the asymptotic behavior of the achievable performance of massive MIMO system, recently a number of coordinated beamforming algorithms have been developed for multicell systems, which usually have reduced inter-cell communication burden. In [27], [28], the signal to interference plus noise ratio (SINR) balancing problem for the multicell massive MIMO system was solved by deriving the asymptotic expression of the achievable SINR based on large-dimension random matrix theory (RMT) [29]- [32]. Furthermore, the asymptotic properties of the SINR were also used to address the power minimization problem subject to some predefined SINR constraints [33]- [35]. More recently, energy-efficient transmission design for massive MIMO systems has also drawn much attention. In [36], the energy-efficient resource allocation problem was studied in OFDMA downlink network with a large number of transmit antennas and fixed maximum ratio transmission (MRT) beamforming. The tradeoff between SE and EE was investigated in a massive MIMO uplink system with fixed receiver, such as maximum-ratio combining (MRC) or zero-forcing (ZF) [37]. More recently, the authors in [38]- [40] investigated the energy efficiency of massive MIMO systems with non-ideal hardware. The research results illustrated that the maximal EE is achieved at some finite number of transceiver antennas when there is a non-zero circuit power per antenna. Similarly, the studies in [41] also demonstrated that it is optimal from the perspective of EE to have a few transceiver antennas when the transceiver design is power inefficient.
In this paper, we investigate joint beamforming optimization and power allocation (JBOPA) for coordinated multicell downlink systems. First, the problem of interest is formulated as maximizing the network utility function defined as the ratio of system weighted sum rate to total power consumption subject to predefined user rate requirements and per-BS power constraints. Due to the non-convex nature of the user rate and the fractional form in the objective, the problem is hard to solve directly. To circumvent this problem, we first resort to an virtual uplink formulations of the the primal problem by introducing an auxiliary variable and the uplink-downlink duality theory. By exploiting the property that the dual problem has an analytic beamforming solution, we propose to further solve the power allocation using a convex approximation. By this means, the dual problem can be efficiently solved and its solution is converted into the downlink using the conventional uplink-downlink theory. We show that the developed algorithm has provable convergence. Then, in order to further reduce the exchange of signalling overhead between coordinated BSs, especially in a large-scale system, a new algorithm used to update the downlink power based on the statistical channel information is developed by exploiting the property of a massive MIMO system. Numerical results validate the effectiveness of our proposed schemes and show that both the spectral efficiency and the energy efficiency can be simultaneously improved traditional downlink coordinated schemes, especially in the middle-high transmit power region.
The remainder of this paper is organized as follows. The system model is described in Section II. Section III reveals the dual uplink form of the primal optimization problem via some necessary optimization transformations. An effective optimization algorithm based on instantaneous channel information is proposed in Section IV by exploiting the newly derived uplinkdownlink duality. An extended effective optimization algorithm based on statistical channel information is further proposed in Section V to reduce the overhead of information exchange between coordinated BSs. In Section VI, numerical evaluations of these algorithms are carried out and conclusions are finally drawn in Section VII.
The following notations are used throughout this paper. Bold lowercase and uppercase letters represent column vectors and matrices, respectively. The superscripts (·) H , and (·) −1 represent the conjugate transpose operator, and the matrix inverse, respectively. A denotes the Euclidean norm for vectors and the Frobenius norm of matrix A. The probability density function (pdf) of a circular complex Gaussian random vector with mean μ and covariance matrix is denoted as CN (μ, ). The function ceil (x) rounds the elements of x to the nearest integers greater than or equal to x. The function f loor (x) rounds the elements of x to the nearest integers less than x. (x, y) + = max (x, y). 1 l N denotes N dimensions column vector with all elements equal one and a.s. − − → denotes almost sure convergence.

II. SYSTEM MODEL
Consider a K -cell MU-MISO downlink system where the BS in cell j ( j = 1, . . . , K ) is equipped with M j transmit antennas and serves N j single-antenna user. We assume M j ≥ N j , ∀ j, thus, user scheduling is not taken into account. Let the k-th user in cell j be User-( j, k) and the BS in cell m be BSm, respectively. Then, the received signal of the User-( j, k) is denoted as where h m, j,k ∈ C M m ×1 = σ m, j,k h m, j,k denotes the flat fading channel coefficient from BS-m to User-( j, k) 1 , in which σ m, j,k represents the large-scale channel information from BS-m to User-( j, k), and h m, j,k has i.i.d. complex entries of zero mean and variance 1 M m [26], [28], w j,k ∈ C M m ×1 denotes the unitnorm beamforming vector for User-( j, k), p j,k and x j,k denote respectively the transmit power and the information signal intended for User-( j, k), and z j,k is a zero-mean circularly symmetric complex Gaussian random noise with power spectrum density ω 2 j,k . Without loss of generality, in the rest of our paper, we assume that each BS serves the same number of users in each cell, i.e., N j = N , ∀ j 2 . For as the normalized channel flat fading where W represents the occupied bandwidth, · and · denote respectively the function ceil (·) and the function f loor (·), n, m = 1, . . . , N , N = K N. Then, the instantaneous rate of where and − → γ m denotes the SINR of User-m N , m − N m N and is given as In this paper, the EE of interest is defined as the ratio of the weighted sum rate to the total power consumption, given by where w denotes the collection of all beamforming vectors, p = p 1 , . . . , p m , . . . , p N T , ξ ≥ 1 is a constant which accounts for the inefficiency of the power amplifier, P c is the constant circuit power consumption per antenna which includes the power dissipations in the transmit filter, the mixer, the frequency synthesizer, and the digital-to-analog converter, P 0 is the basic power consumed at the BS which is independent of the number of transmit antennas 3 [36], α m = α m N ,m−N m N denotes the 2 The developed algorithms can be implemented in the case where different BSs serve different number of users in each cell, i.e., N i = N j , i = j. 3 Note that in communication system, a higher R m implies a higher code rate or higher order modulation being used, thereby may increase the power consumption of the encoder and decoder circuits. This leads to that the power consumption in the baseband is a function of R m , which is usually nonlinear and makes the EE optimization problem extremely complex. For simplicity, similar to the existing literature [7]- [17], here this fact is ignored and we assume that the power consumption at the baseband is a constant and is also included in P 0 . This is reasonable due to the fact that the variation of power consumption for different levels of code rates or modulations is in general marginal in contrast to the basic power consumption.
weighting factor that is used to represent the priority of Userm N , m − N m N in the system, which is typically included to achieve a certain fairness index among the users. For example, to obtain the proportional fairness, we can set where R m is the average data rate of User-m N , m − N m N in the previous time slots [42], [43]. In order to balance the SE and the EE, our goal is to maximize the EE of the K -cell MU-MISO downlink transmission system subject to some predefined user rate requirements and per-BS transmit power constrains, i.e., where K = {1, . . . , K }, U = 1, . . . , N , P j is the transmit power constraint of BS-j, and γ m represents the target SINR of User-m N , m − N m N . Note that, different from the EE optimization problem considered in the existing literature [14]- [19], the per-user minimum rate requirements are considered in our problem formulation making more challenging. In what follows, we attempt to find its solution from the perspective of uplink-downlink duality and then design an effective algorithm. Also, it is worth mentioning that though there is a feasibility issue in the above optimization problem (5), it is similar to the power minimization problem subject to a given minimum target rate demand set as in [34]. Throughout this paper, we assume that the set of per-user target rate requirements is feasible and only focus on finding the solution to the problem.

III. PROBLEM TRANSFORMATIONS
It is easy to understand that problem (5) is non-convex in the optimization variables and therefore it is difficult to find the globally optimal solution, since the user rate is non-convex in general for the MISO interference channel. Furthermore, the objective function in a fractional form brings more obstacles to solve the problem directly. In what follows, we focus on finding an efficient optimization method to solve problem (5) by extending the conventional uplink-downlink duality with reformulating it into the following form based on the per-BS constraints.
where β is a newly introduced auxiliary variable, β * denotes the ratio of the minimum total transmit power that needs to satisfy the target SINR constraints to the total power constraints 4 , and The equivalent relation between problem (5) and problem (6) is given by the following theorem. Theorem 1: Problems (5) and (6) have the same optimal objective value and the same optimal solutions.
Proof: Let w opt and p opt be the optimal solution to problem (5) and let w opt , p opt and β opt be the optimal solution to problem (6), respectively. One can easily prove by contradiction that N m=1 p opt m = β opt K j=1 P j . Otherwise, β opt can be scaled down to achieve a better objective value. On the one hand, it is easy to know that w opt and p opt are also a feasible solution to problem (5) due to the fact that the constraint sets of problem (5) is a subset of the constraint sets of problem (6). On the other hand, letβ = 1 K j=1 P j N m=1 p opt m , recalling the per-BS transmit power constraints, we have 0 ≤β ≤ 1 and easily know that w opt and p opt are also a feasible solution to problem (6) in conjunction with the rate constraints. With these two observations, we can prove by contradiction that − → f w opt , p opt = − → g w opt , p opt , β opt must hold and the same objective is achieved with the same beamforming and power solution. Otherwise, the objective with a lower value can be improved by replacing the beamforming matrix and power vector with the solution of the other problem, which contradicts the assumption that these two objective values are both optimal ones.
It is easy to verify that the achievable objective value of (6) monotonically increases with the auxiliary variable β in the range β * ≤ β ≤ β opt , and monotonically decreases with the auxiliary variable β in the range β opt ≤ β ≤ 1. Thus the optimal value of β can be efficiently obtained through a one dimension search method such as the simulated annealing algorithm, or the hill-climbing method in the range β * ≤ β ≤ 1 [44], [45]. Based on these analysis, in what follows, the original EE optimization problem is first transformed into a sum-rate maximization problem subject to per-BS transmit power constraints, per-user target rate requirements and additional sum power consumption, i.e., (8), where the value of the auxiliary variable β is fixed. Then, we solve iteratively (8) by updating β via a one dimension search method to obtain the solution to the energy efficiency optimization problem (6).
In contrast to the conventional WSR maximization (WSRMax) problem, problem (8) has additional total power constraint and minimum SINR requirements, which is named as extended WSRMax problem, making the problem more difficult to solve due to the coupling between the optimization variables. In order to design an effective algorithm, in the sequel, the uplink-downlink duality is applied to solve problem (8) which can be rewritten as the following problem.
where λ = {λ j ≥ 0, K +1 j=1 λ j = 0} are auxiliary variables, and Note that a feasible solution to problem (8) is also feasible for problem (9). It also means that the optimal objective value of problem (9) for any given λ ≥ 0, is an upper bound on the optimal objective value of problem (8). Furthermore, it has been shown in [46], [47] that the solution to problem (8) can be equivalently achieved by solving the problem min λ ψ (λ). One can easily verify that the function ψ (λ) is a convex function but is not necessarily differentiable 5 , so we can solve the problem min λ ψ (λ) via the sub-gradient algorithm or ellipsoid algorithm. Furthermore, it is shown in [48] that such an algorithm is guaranteed to converge to its optimal solution. The sub-gradient of the function ψ (λ) is given by the following proposition, that can be proved with a similar method as used in [46].
Proposition 1: The sub-gradient g of the function ψ (λ) with fixed parameters λ where p ( * ) is the optimal solution of problem (9) with fixed parameters λ = λ and β.
The conventional uplink-downlink duality theory in [49], [50] shows that problem (9) is dual to the following virtual uplink problem: in which q ∈ R N ×1 + denotes the virtual uplink transmit power vector with q m = q m where and is calculated as According to the uplink-downlink duality theory [46], [47], the uplink and the downlink have the same achievable SINR region with the same beamforming vectors and the same total transmit power constraint. This means that once the uplink problem is solved, the beamforming solution to the downlink problem can be obtained at the same time, while the downlink power solution remains to be calculated. Using the fact that the uplink and the downlink achieve the same SINR, we have By rearranging the terms in (14), we have The above equalities can be used to calculate the downlink power solution from the corresponding uplink solution.

IV. JBOPA ALGORITHM WITH INSTANTANEOUS CHANNEL INFORMATION
In this section, we focus on finding the solution to the rate constrained EE maximization problem (5). Section III has revealed that solving (5) is equivalent to iteratively solving the extended WSRMax problem (8) for fixed β and searching for the optimal value of β, which can be further recast to solve the virtual uplink problem (11).

A. Algorithm Design
It is well known that the optimal receiver beamforming vector that maximizes the uplink SINR expressed in (13) is the MMSE filter, which is in general calculated as 6 [50] In what follows, we focus on the maximization over q of problem (11) for fixed w and β that can be formulated as follows: 6 Note that though the constants P 0 , P c and ζ which model the power consumption in transceiver circuits and the power amplifier efficiency affect the energy efficiency through the denominator of the energy efficiency defined in (4), P 0 and P c are independent of the expression of the normalized beamformer w, as shown in (17). In other words, they are only related with the search over β and thus has an impact on the value of energy efficiency. It also means that a possible mismatch between the true value of P c , P 0 and the value that is used to design the proposed beamformer only impacts the value of energy efficiency but does not impact the expression of the beamformer. Furthermore, it is worth mentioning that the considered problem (5) is more complex than the energy efficient design in [15]- [19] where the user rate requirements is not taking into account. Different from them, we solve the complex energy efficiency optimization problem from the perspective of uplink-downlink duality by iteratively updating an auxiliary variable. In addition, we would like to point out that the solution of the beamforming vector w given by (17) is the optimal receiver to the dual problem (11) of the extended WSRMax problem with a fixed β, as a result ζ is not included in the structure of the normalized beamformer w.
where γ denotes the collection of the auxiliary variables { γ m }. The optimization problem (18) is a signomial optimization problem and is non-convex in general. Fortunately, note that, all constraints are monomial function in problem (18). If therefore the objective function can reformulated as a monomial function, problem (18) become a GP problem in standard form. Motivated by this observation, a successive convex approximation approach is first presented in the following proposition, for the proof please see [52].
To overcome the non-convex nature of problem (18), we start with reformulating problem (18) by applying the local approximation of Lemma 1 in its objective function, yielding where δ is a small constant used to control the desired approximation accuracy, typically set to be δ = 0.1. Notice that problem (20) is lower boundary problem of the original problem (18) that is difficult to solve directly. It has been proven in [53] that problem (20) can be reformulated as a convex GP problem in standard form by using jointly the logarithmic change of the variable and a logarithmic transformation of the objective function and constraints, respectively. In other words, problem (20) can be easily solved by using the powerful GP optimization tool packets [54]. Furthermore, the results in [52] have shown that by starting from an initial point, we can search for a close local optimum by solving a sequence of GPs that locally approximate the original problem (18). Based on the above analysis, an effective JBOPA optimization algorithm is developed to solve the virtual uplink problem (8) with fixed β and is summarized as Algorithm 1 7 . According to the uplink-downlink duality, it is known that the optimal uplink beamforming vectors w also serve as the optimal downlink ones, while the optimal downlink power allocation p can be calculated from the uplink solution. The details are discussed in the following subsection.
Remark 1: Note that when the value of β equals to one, Algorithm 1 is equivalent to solving the conventional WSRMax problem. In order to solve the original EE maximization problem, we need to iteratively run Algorithm 1 to obtain the solution of β by using a one dimension search method that usually only requires a limited number of iterations. Let C be the computational complexity of Algorithm 1, the computational complexity of our proposed EE maximization algorithm is O N β C where N β is the required number of iterations in searching for an optimal β [44], [45].

B. Downlink Transmit Power Computation
It is necessary to point out that the sub-gradient given in (10) is calculated based on a given downlink transmit power p. Therefore, it is important to compute the dowlink transmit power from the obtained uplink transmit power q and transmit beamforming vector w. Combining (14) with (16) where P max = N m=1 q m . After some basic manipulations, we where matrices G and D are given by, respectively Multiplying both sides of (22a) by λ T = λ 1 , . . . , λ m , . . . , λ N T and using (22), we have Defining an extended power vectorp = p 1 and an extended coupling matrix Thus, we have the following eigensystem by combining (22a) and (26):p According to the conclusions in [51], we can easily obtain the optimal power vectorp as the first components of the dominant eigenvector of Q w, λ, ← − γ , P max , which can be scaled so that its last component equals one.

C. Algorithm Convergence and Complexity Analysis
Notice that at each iteration the achieved virtual uplink energy efficiency is monotonically non-decreasing due to the fact that steps 3 and 5 in Algorithm 1 all aim to maximize the objective function in problem (7), i.e., ← − g (0) ≤ ← − g (1) ≤ . . . ≤ ← − g (n) ≤ . . .. Since the achievable rate region under the sum power constraint is bounded, the convergence of Algorithm 1 is guaranteed by the monotonic convergence theorem [56]. Meanwhile, the studies in [46] have revealed that the convergence of the proposed update method based on the sub-gradient ellipsoid method for the auxiliary variables λ is guaranteed to minimize ψ (λ). Thus, the convergence of the whole algorithm is guaranteed and the convergence speed of our proposed algorithm is similar with that of the proposed algorithm developed in [46]. It is easy to see that the major computational load involes steps 3 and 5 in Algorithm 1. The GP operation at step 3 fortunately has very fast convergence speed and low computational complexity [53]. The matrix inversion at step 5 can be computed efficiently with In this section, the joint beamforming optimization and power allocation are studied based on the statistical channel information to reduce the frequency of power updates and the overhead of signal exchange between coordinated BSs 8 . The results in [26]- [28] have shown that the SINR of each user converges to different deterministic values when BS is equipped with a large number of transmit antennas. In what follows, these observations can be used to obtain an approximated power allocation problem related only to the statistical channel information which is used to optimize the transmit power. By plugging (17) in the expression of the virtual uplink SINR of User-m N , m − N m N , the asymptotic approximation for ← − γ m is given in the following Theorem 2.
Theorem 2: The instantaneous random variable ← − γ m can be approximated by a deterministic quantity where m = Proof: Please refer to Appendix A. In order to reduce the overhead of signal exchange between coordinated BSs and apply the conclusion obtained in Theorem 2, in the sequel, we resort to achieve the power allocation solution by solving an approximated new power allocation problem directly defined as (30) which only relies on the large scale fading coefficients in place of solving problem (11) which is related with the instantaneous CSI 9 .
It is easy to verify that problem (31) is convex and thus can be tackled using standard convex approaches [45]. Particularly, the power constraint meets with equality in the optimal solution. The corresponding Lagrange function is given by where ξ and ς are respectively the Lagrange multipliers associated with the SINR constraints and the power constraint. The Karush-Kuhn-Tucker (KKT) conditions are then given by Note that ξ ≥ 0 acts as a slack variable in the last equation and therefore can be removed. The above equations are rewritten as Using a similar approach as in [45, pp. 245-246], we have or, simply expressed as q m = α m ln (2) It is seen that the power solution in the energy-efficient largescale system is water-filling-like power allocation algorithm. As noted in Proposition 1 and Algorithm 1 the update of auxiliary variables λ needs to know the downlink transmit power p obtained withp = Q w, λ, ← − γ , P max p. A similar procedure for analyzing the instantaneous random variable ← − γ m can be applied to the matrices G and D such that the asymptotic approximation that only depends on the statistical CSI of the matrix Q w, λ, ← − γ , P max can be obtained. The asymptotic approximation of G and D are described as the following Theorem. Proof: Please refer to Appendix A. In what follows, in order to realize the ideas of power allocation based on statistical channel information, the extended coupling matrix Q w, λ, ← − γ , P max is approximated by the following matrix.
[D] m,n = Based on the above analysis, the power allocation problem (30) can be solved via the alternating optimization approach. To be more specific, we solve the power allocation problem (30) by sequentially fixing two of the three variables β, λ, q and updating the third. The details of the alternating optimization algorithm based on duality theorem that is used to solve the power allocation problem (30) with fixed β based on statistical channel information is summarized as Algorithm 2, where ρ and denote the large-scale system EE and the objective value of problem (31), respectively. Remark 2: In Algorithm 2, the power allocation is optimized based on statistical channel information and thereby is carried out in a long-term timescale, in which only the large-scale channel information and a few real numbers are exchanged between the BSs. While the beamforming vector w m can be updated based on local instantaneous channel state information, i.e., h m,n , ∀n and thereby on a short-term basis. Compared with Algorithm 1, the computational complexity can be significantly reduce due the long-term timescale power allocation update. In addition, Algorithm 2 only needs one time beamforming vector computation for each change of the instantaneous CSI and only needs the local CSI to compute the beamforming vector. Table I lists the achieved SR in terms of 10 7 bps of Algorithm 1 and Algorithm 2 for several antenna configurations over 1000 random channel realizations. It is easy to see that the difference between the achieved SR of Algorithm 1 and that of Algorithm 2 decreases with an increasing number of transmit antennas indicated by the percentage between the achieved SR of Algorithm 2 and that of Algorithm 1.

VI. NUMERICAL RESULTS
In this section, we investigate the performance of the proposed JBOPA algorithms via numerical simulations. We consider a cooperative cluster of K = 3 hexagonal adjacent cells each consisting of one BS and multiple users. In cell j, BS-j is equipped with M transmit antennas and serves N single antenna users. The cell radius is set to 500 m and each user randomly locates in the cooperation region and has a least 400 m distance from its serving BS as shown in Fig. 1. The channel vector h m, j,k from BS-m to User-( j, k) is generated based on the formulation h m, j,k ∈ C M×1 = σ m, j,k h w m, j,k = β m, j,k h w m, j,k , where h w m, j,k denotes the small scale fading part and is assumed to be Gaussian distributed with zero mean and variance 1 M , and β m, j,k denotes the large scale fading factor which in decibels is given as 10 log 10 (β m, j,k ) = −38 log 10 (d m, j,k ) − 34.5 + η j,k , where η j,k represents the lognormal shadow fading with zero mean and standard deviation 8 dB, d m, j,k denotes the distance between the BS-m and User-( j, k) [57]. The circuit power per antenna is P c = 30 dBm, and the basic power consumed at the BS is P 0 = 40 dBm [36]. The noise power spectrum density ω 2 j,k = −174 dBm/Hz, ∀ j, k and the occupied bandwidth W = 10 MHz. For simply, we assume that each BS has the same transmit power budget P, i.e., P j = P, ∀ j. The noise figure is 9 dB [4]. The weighting factor α and the inefficiency factor of power amplifier ξ are all set to unit.
For comparison, several related existing algorithms are also evaluated in our numerical simulations, including the power minimization algorithm that aims to minimize the transmit power subject to the same constraints, the conventional WSRMax algorithm that aims to maximize the SE with the same constraints, the performance of an extension of the algorithm developed in [15] that is described in Appendix D, and the normalized MRT beamforming transmission with equal power allocation at each BS and subject to no rate constraints. In our simulations, the target rate of each user is set respectively to be 70% of the user rate achieved by the normalized MRT beamforming with equal power allocation. Fig. 2 illustrates the average EE of the above coordinated beamforming algorithms for finite-size system over 1000 random channel realizations. It can be seen that our proposed Algorithm 1 achieves the best EE performance among all algorithms, with a considerable gain at the high transmit power region. It is interesting to see that the Algorithm 1 and the SRMax algorithm achieve the same EE performance in the lower transmit power region such as 26-38 dBm. This implies that in this region it is possible to simultaneously achieve the optimal EE and SE, while in other regions a reasonable tradeoff between SE and EE is desired, which is the goal of our algorithm design. The results also show that our proposed Algorithm 1 outperforms the power minimization algorithm in all regions. This suggests that the power minimization subject to rate constraints in general cannot achieve the optimal EE, in which all the achieved user rates exactly meet with the constraints. In other words, it is very possible that achieving rates greater than the constraints also bring improved SE. Fig. 3 illustrates the average WSR performances of some coordinated beamforming algorithms for finite-size system over 1000 random channel realizations. It can be seen that Algorithm 1 and the SRMax algorithm achieve the same SE performance in the lower transmit power region. In other words, the optimal EE algorithm also achieves the optimal system SE in the lower transmit power region, and vice versa. However, in terms of the SR, in the middle-high transmit region, the SRMax algorithm achieves the best WSR performance at the cost of a certain EE loss. Fig. 4 illustrates the average EE of the proposed algorithm and the WSRMax algorithm achieved with imperfect CSI,  under configuration {M, N } = {4, 2}. In our numerical simulation, only an imperfect estimate h m, j,k of the true channel h m, j,k is available at the transmitter which is modeled as [26] where h e m, j,k denotes the estimation error, h m, j,k and h e m, j,k have i.i.d entries of zero mean and variance 1 M and independent of each other and independent of z j,k , β m, j,k denotes the large scale fading factor which in decibels is given as 10 log 10 (β m, j,k ) = −38 log 10 (d m, j,k ) − 34.5 + η j,k , where η j,k represents the log-normal shadow fading with zero mean and standard deviation 8 dB, d m, j,k denotes the distance between the BS-m and User-( j, k) [57]. The parameter δ ∈ [0, 1] reflects the accuracy or quality of the channel estimate h m, j,k , i.e., δ = 0 corresponds to perfect channel state information, whereas for δ = 1 the imperfect estimate h m, j,k is completely uncorrelated to the true channel. Numerical results show that the EE performance achieved by the related algorithms become worse with an increasing value of δ. It is see that under different levels of imperfect CSI, the EE performance achieved by Algorithm 1 is always better than that of the WSRMax Algorithm at the high transmit power region. Fig. 5 illustrates the performance comparison of the two developed coordinated beamforming algorithms versus the number of transmit antennas over 1000 random channel realizations, where the the ratio between the number of transmit antennas and the number of served users is set to be 16. Numerical results show that the performance achieved by the proposed Algorithm 2 is very closed to that of the proposed Algorithm 1 with lower computational complexity and limited scale information between the coordinated BSs. This suggests that the proposed Algorithm 2 based on statistical state channel information is feasible for using the determination approximate expression of SINR to optimize the transmit power while the beamformers are calculated based on the instantaneous channel information.

VII. CONCLUSIONS
In this paper we have studied joint beamforming and power allocation in multicell downlink systems to maximize the system energy efficiency subject to per-BS transmit power constraints and per-user target SINR requirements. The considered objective function in fractional form is non-convex and it is hard to find the globally optimal solution. In order to obtain a tractable form, the primal problem was investigated from the viewpoint of the uplink-downlink duality by introducing some auxiliary variables. Furthermore, the uplink power allocation problem with fixed beamformers can be efficiently solved by using standard GP solvers. The convergence of the proposed JBOPA Algorithm 1 has been proven with the monotonic boundary theorem and the property of sub-gradient of convex problem. To reduce the overhead of the CSI collection, a new power allocation algorithm based on statistical channel information has been further developed. Numerical results have validated the effectiveness of our proposed schemes and shown that both SE and EE performance can be simultaneously improved over traditional downlink coordinated schemes, especially in the middle-high transmit power region.

A. Related Proof
Useful Results from Random Matrix Theory: We reproduce the following Lemmas [29]- [32].
Lemma 2: Let x ∈ C N , independent and identically distributed (i.i.d.), with zero mean and variance 1 N , A ∈ C N ×N Hermitian with bounded spectral norm whose elements are independent of x, then , j=1 whose the entries are given by: , with the following assumptions hold. S1: The complex random variables S2: There exists a real number σ max < ∞ such that: sup whose elements are the unique solutions of the deterministic system of N + M equations where 1 N tr ( (z)) is the Stieltjes transform [31] of the probability measure.
The values of i (z) and j (z), ∀i, j can be obtained by initializing them to known values and iterating over equations (43) and (44) until their values converge. Furthermore, the differential of the Stieltjes transform of the matrix AA H can be calculated with the method used in [33].
where M m (z) denotes the Stieltjes transform of the matrix m evaluated at the point z which can be computed using Lemma 3.
[D] m,n = Plugging (48) in (51) yields the deterministic uplink-dowmlink power transformation matrix G in Theorem 3. This completes the proof of Theorem 3.

B. JBOPA Initialization Algorithm With Instantaneous Channel Information
Before proceeding to propose an effective method to solve problem (11), we would like to discuss the feasiblity of the target SINRs in the considered systems [51]. In general, the PMP subject to SINR constraints can be adopted to validate the feasiblity of some given SINR subject to transmit power constraint in wireless communications. Herein, we propose to solve the PMP that is formulated as (17) for fixed feasible target SINRs, given by By using the similar derivation method used in [27], [50], we can obtain the dual problem of problem (14) that is given by the following optimization problem: where γ m denotes the SINR of User-m N , m − N m N and is calculated as where q m represents the dual uplink power of Userm N , m − N m N , ν m N denotes the virtual uplink noise variance of BS-m N . Similarly, it is easily known that the maximization over w in (53) can be solved explicitly for fixed power q. In particular, the optimal receiver beamforming vector w k that maximizes the SINR is the MMSE filter. According to the Sherman-Morrison identity, the MMSE receiver filter can be equivalent to the following form: Bases on this, problem (53) can be equivalent to the following optimization problem by substituting the expression of beamformers w given by (55) into the target SINR constraints: It is easy to see that the inner optimization of problem (53) Similarly, we have the following eigensystem where P max = N m=1 q m , ν is a N dimension column vector and is defined as ν m = ν m N , ∀m. The details of the optimal uplink power allocation algorithm for finite-size system are summarized as Algorithm 3 where g denotes the subgradient of variables ν given by g = g 1 , . . . , g j , . . . , g K with g j = j N m=( j−1)N +1 p m and S ν = ν : ν ≥ 0, K j=1 1 − ν j P j = 0 [50]. Furthermore, the convergence of Algorithm 3 can be guaranteed by the fixed function theorem [49]. k + θ g k with p ( * ) , and obtain ν ( * * ) , let q (0) = q ( * ) , and ν ( * ) = ν ( * * ) , then go to step 2.

C. JBOPA Initialization Algorithm With Statistical Channel Information
Similar to the JBOPA Algorithm 1, we also need to solve the initialization problem of the virtual uplink transmit power before applying the developed power allocation optimization Algorithm 2 to solve effectively problem (6). Similar procedure for obtaining the deterministic quantity ← − ϒ m of ← − γ m can be applied also to γ m in order to examine the asymptotically optimal solution of the PMP (53) based on statistical CSI.  For fixed ν, the optimal solution of problem (60) with respect to variables q can be easily obtained by using fixed equation theory. In addition, the optimization of problem (60) with respect to variables ν can also be solved by using the subgradient method. The details of the algorithm that is used to solve problem (60) in large-scale system is summarized as Algorithm 4. k + θ g k , then obtain ν ( * * ) , let q ( * ) = q (0) , q (0) = q ( * * ) , ν ( * ) = ν ( * * ) and go to step 2.

D. An Extended JBOPA Optimization Algorithm
It is interesting to find that the authors in [15] also adopted the convex approximation method to investigate the energy efficiency optimization problem for single cell multiuser downlink broadcast channel. It is easy to see that the developed algorithm in [15] can be easily extended to solve the considered energy efficiency optimization problem. To optimize the transmit power, the following lemma is applied to approximate the objective function of problem (5).
where μ = ξ N j=1 p m + K j=1 M j P c + P 0 . Similar to problem (23), problem (62) can also be solved by using the powerful this, the convergence of Algorithm 5 is also guaranteed by the monotonic boundary theorem [56].
In addition, problem (60) has M = 2N K j=1 M j real optimization variables, N + 1 second-order-cone (SOC) constraints where each of them consists of M real dimensions. It is known from [48] that the computational complexity of the power minimization problem in terms of number of iterations is upper bounded by O N + 1 , and the complexity of each iteration is within the order of O 2 M 3 N . Thus, the total worst-case computational complexity of PMP is given by O 2 N + 1 M 3 N . Also, since both the Algorithm 1 and the extended Algorithm 5 include a GP step, we ignore this GP step in the complexity analysis.