Sum Rate Optimization of IRS-Aided Uplink Muliantenna NOMA with Practical Reflection

Recently, intelligent reflecting surfaces (IRSs) have drawn huge attention as a promising solution for 6G networks to enhance diverse performance metrics in a cost-effective way. For massive connectivity toward a higher spectral efficiency, we address an intelligent reflecting surface (IRS) to an uplink nonorthogonal multiple access (NOMA) network supported by a multiantenna receiver. We maximize the sum rate of the IRS-aided NOMA network by optimizing the IRS reflection pattern under unit modulus and practical reflection. For a moderate-sized IRS, we obtain an upper bound on the optimal sum rate by solving a determinant maximization (max-det) problem after rank relaxation, which also leads to a feasible solution through Gaussian randomization. For a large number of IRS elements, we apply the iterative algorithms relying on the gradient, such as Broyden–Fletcher–Goldfarb–Shanno (BFGS) and limited-memory BFGS algorithms for which the gradient of the sum rate is derived in a computationally efficient form. The results show that the max-det approach provides a near-optimal performance under unit modulus reflection, while the gradient-based iterative algorithms exhibit merits in performance and complexity for a large-sized IRS with practical reflection.


Introduction
Intelligent reflecting surfaces (IRSs) have drawn enormous attention from the academy and industry as a cost-effective building block for 6G wireless communication networks demanding high spectral and energy efficiency [1][2][3]. An IRS constructed with a large number of passive reflection elements can reconfigure a wireless propagation channel to be favorable for information or energy transfer by controlling their reflecting patterns. In doing so, the IRS avoids a large power consumption with passive elements and achieves a full-duplex gain without complicated signal processing such as interference cancellation and demodulation. In this aspect, various IRS-assisted wireless communications have been explored for their own purposes, from basic multiuser or/and multiantenna communication systems [4][5][6][7][8][9] to more complicated system configurations, as surveyed in Ref. [10].
For multiuser communications, nonorthogonal multiple access (NOMA) craving for higher spectral efficiency and massive connectivity has been considered as a promising candidate for 6G networks [11][12][13][14][15][16]. Hence, many recent IRS studies have been devoted to IRS-assisted NOMA for further improvement in spectral efficiency, energy efficiency, and reliability [17][18][19][20][21][22][23][24][25][26][27]. One of the major concerns of IRS-aided NOMA networks resides in how to reflect a superimposed NOMA signal appropriately to meet the design goal of a given network. In the downlink, joint optimization of IRS reflection and power allocation was studied with a single-antenna base station (BS) to minimize the transmit power of a two-user NOMA signal [18] and maximize the sum rate of a multiuser NOMA signal [19]. For a multiantenna BS in the downlink, transmit beamforming and IRS reflection were optimized without or with power allocation to maximize the sum rate [20], minimize the transmit power [21][22][23], and maximize the minimum rate [24]. On the other hand, in the uplink, IRS-aided NOMA was studied to maximize the sum rate achieved with a singleantenna receiver by optimizing the IRS reflection vector [27]. The sum rate maximization problem was also extended to IRS-aided NOMA accompanied by wireless power transfer with a single-antenna receiver [28,29] and a multiantenna receiver [30].
The recent studies on IRS-aided NOMA networks have assumed unit modulus reflection for IRS since the phase control is more readily implementable than the amplitude control in implementing passive IRS elements. It should also be noted that IRS optimization under unit modulus reflection has resorted to semidefinite relaxation (SDR) solving a semidefinite program (SDP) [31] for the downlink with a given transmit beamforming and for the uplink with a single-antenna BS, where the SDP deals with linear functions of semidefinite matrices [20,23,24,[27][28][29]. However, the SDR approach is less favorable for an IRS requiring a large number of passive elements to compensate for double fading due to its complexity increasing polynomially with the number of elements. To reduce the complexity of the SDR, an iterative algorithm based on a second-order surrogate function similar to the gradient descent was proposed in Ref. [5] for the simple objective function given by the trace of a matrix.
For the uplink NOMA with a multiantenna BS [30], the sum rate is given in a more complicated log-determinant expression similar to a multiple-input multiple-output (MIMO) capacity [7]. In this case, the SDP, having linear functions of semidefinite matrices in the objective and constraints, is no longer applicable. A few studies have dealt with such a complicated objective function in IRS reflection optimization [7,30], which was optimized through a suboptimal sequential optimization method, optimizing one IRS element while fixing the other elements without knowing its performance gap to the optimal one. In addition, practical IRS reflection models with phase-dependent amplitudes observed in practical circuits [32] have not been studied for IRS-aided NOMA networks yet. In this regard, we consider a sum rate maximization problem of an IRS-aided uplink NOMA network with a multiantenna BS and tackle the problem under a generalized IRS reflection model encompassing unit modulus and practical reflection [32]. The main contributions are summarized as follows: • We formulate the sum rate optimization problem by incorporating the optimal receive beamforming for a given IRS phase vector in the objective function, which is characterized by the log-determinant of a matrix as in the MIMO capacity expression [7,30]. As a result, only the IRS phase vector needs to be optimized in the uplink, while the transmit beamforming and IRS phase vectors have been optimized alternately in the downlink [5]; • For a moderate-sized IRS with unit modulus reflection, we propose an extended SDR approach that converts the sum rate maximization problem into a determinant maximization (max-det) problem [33]. The max-det solution not only provides an upper bound on the optimal sum rate of the IRS-aided uplink multiantenna NOMA, but also leads to a rank-one feasible solution resulting in a near-optimal performance. This approach can also be employed to obtain an upper bound on the IRS-aided MIMO capacity under unit modulus reflection; • For a large-sized IRS under generalized reflection, we reformulate the problem as an unconstrained nonlinear optimization problem that can be solved by using gradient-based iterative algorithms [34]. In particular, we use the more sophisticated Broyden-Fletcher-Goldfarb-Shanno (BFGS) and limited-memory BFGS (L-BFGS) algorithms [34,35], while the gradient-descent approach is used in [5]. For an efficient implementation of such iterative algorithms, we derive the gradient of the complicated objective function in a computationally efficient form under the generalized reflection; • We analyze the computational complexity of the iterative algorithms when the derived gradient is used to update the search point. The results show that the iterative algorithms reduce the complexity of the extended SDR with max-det optimization significantly. In addition, the iterative algorithms provide a performance gain over the conventional sequential optimization method [7] at a reduced computational time.
Notation: The sets of n × m complex-valued and real-valued matrices are denoted by C n×m and R n×m , respectively, with C n = C n×1 and R n = R n×1 , while the set of n × n positive semidefinite Hermitian matrices is denoted by S n + . The transpose, Hermitian, and trace are denoted by (·) T , (·) H , and tr(·), respectively. We use diag(a) for the diagonal matrix with a diagonal vector a, [a] n for the nth entry of a vector a, and, [A] n,m for the (n, m)th entry of a matrix A, and CN (µ, Σ) for complex Gaussian distribution with mean vector µ and covariance matrix Σ.

System Model
We consider the uplink of a single-cell network described in Figure 1. The network consists of a BS equipped with M antennas, K devices equipped with a single antenna, and an IRS comprising N reflection elements. The channels from device k to the BS and to the IRS are denoted by v k ∈ C M and f k ∈ C N , respectively, for k ∈ K {1, 2, · · · , K}. The channel from the IRS to the BS is denoted by G ∈ C M×N . The IRS reflection vector is denoted by θ = [θ 1 , θ 2 , · · · , θ N ] T ∈ C N , where θ n = e jφ n for n ∈ N {1, 2, · · · , N} under unit modulus reflection. To address the amplitude distortion of practical IRS control circuits depending on the phase, we express the IRS reflection in a generalized form as [32] where with α ≥ 0, β min ≥ 0, and φ 0 ≥ 0. The values for parameters α, β min , and φ in (2) are determined by the specific circuit implementation, where (2) with α = 0 represents the unit modulus reflection with β n (φ n ) = 1. For the uplink transmission, we allow K devices to transmit their symbols simultaneously, where the number K of devices is larger than the number M of receiving antennas for NOMA. However, the following results are also applicable to space division multiple access with K ≤ M. The signal received at the BS is then written as where s k and p k are the symbol and transmit power of device k ∈ K, respectively, and z ∼ CN (0, σ 2 I M ) is the noise vector added at the BS. We can express the received signal (3) in a concise form as where H k is the equivalent channel from device k to the BS, which is given by Without loss of generality, we assume that the devices are indexed in the successive interference cancellation (SIC) order. The BS detects the device symbols from s 1 to s K sequentially by applying receive beamforming w k ∈ C M to the received signal after SIC in detecting s k , specifically by applying the receive beamforming w k to the received signal after {s 1 , s 2 , · · · , s k−1 } being detected and canceled, which is given bŷ whereŷ 1 = y. We obtaiñ from which s k can be detected. Thus, the signal-to-interference-and-noise ratio (SINR) in detecting s k fromỹ k is given by The optimal receive beamforming that maximizes the SINR is given by the minimum mean square error (MMSE) beamforming, expressed as [36] for k ∈ K and B K+1 = σ 2 I M . The maximum SINR of device k achieved with the optimal beamforming w o k is given by which leads to the achievable rate as follows: using the matrix determinant lemma, det(B + uu H ) = det(B) det(1 + u H B −1 u) for an invertible matrix B [15,36], and , which leads to (12). From (12), the sum rate of all devices is given by Refs. [15,36] where to the SIC order. Note that device rates in (8) depend on the SIC order since the SIC order affects B k for 2 ≤ k ≤ K. However, the sum rate (13) determined by B 1 and B K+1 does not depend on the SIC order. Finally, the sum rate is expressed as a function of transmit power p = [p 1 , p 2 , · · · , p K ] T of the devices and IRS phase shifts whereθ is a function of φ.

Problem Formulation
This paper aims to maximize the sum rate of the IRS-aided uplink NOMA by optimizing the transmit power p = [p 1 , p 2 , · · · , p K ] T of the devices and the phase shifts Since the sum rate (14) is a non-decreasing function of p k irrespective of φ, the optimal power of problem (15) is given by that results in the sum rate as with and Finally, the sum rate optimization problem (15) becomes It should be noted that problem (19) under unit modulus reflection is equivalent to a subproblem of the IRS-aided MIMO capacity optimization problem [7], which was solved by the customized sequential optimization method. In the following, we will provide alternative methods providing either a better performance under unit modulus reflection or a faster computation under generalized reflection than the conventional method [7].

Determinant Maximization for a Moderate-Sized IRS
This subsection tackles the problem (19) for a moderate N by extending the SDR approach. For this purpose, let us rewrite the signal matrix (18) as where Ξ = diag([ξ 1 , ξ 2 , · · · , ξ K ] T ), H = [H 1 , H 2 , · · · , H K ] ∈ C M×K(N+1) , and ⊗ denotes the Kronecker product. The signal matrix is a Hermitian semidefinite matrix and is linear withX whenX =θθ H . We defineX =θθ H ∈ S N+1 + , of which the diagonal entries satisfỹ X n,n = |θ n | 2 ≤ 1. In this case, we can transform (19) s.t.X n,n ≤ 1, n = 1, 2, · · · , N + 1, LetX o andR o =R(X o ) denote the optimal solution of (21) and the corresponding optimal sum rate, respectively. Since findingX o is intractable due to the rank constraint, we resort to an approximate solution by finding a rank-relaxed solution first and then estimating a rank-one solution from the rank-relaxed one as follows.
In the first step, by relaxing the rank constraint (21c), we can approximate (21) to the max-det problem defined in [33], which is expressed in the standard from as follows: Since (22) is known to be a convex optimization problem with linear matrix inequalities in (22b) and (22c), it can be solved with an existing convex optimization solver, which leads to the optimal solutionX andR =R(X ). Here, the optimal valueR of the rank-relaxed problem (22) can serve as an upper bound on the optimal valueR o of the rank-constrained problem (21), since the constraint set of (22) includes that of (21).
In the second step, from the solutionX of (22), we obtain a rank-one solution close toX o =θ oθ H o though the Gaussian randomization procedure [24,31]. Since any rank-one solutionX of (21) is decomposed asX =xx H , the Gaussian randomization procedure generates L zero-mean complex Gaussian samples as {x l ∼ CN (0,X )} L l=1 for a rank-one solution so thatx lx We then obtain a feasible candidateX l subject to [X l ] n,n ≤ 1 for the solution of (21), which is equivalent To generateθ l , we obtain the phases as [φ l ] n = ∠[x l ] n for n ∈ N. The phase vectors {φ l } L l=1 generated by the phases of complex numbers are feasible for problem (19) and {X l =θ lθ H l } L l=1 are also feasible for problem (21) since [X l ] n,n ≤ 1 and rank(X l ) = 1. Among all the feasible candidates {φ l } L l=1 (or equivalently {X l } L l=1 ), the Gaussian randomization procedure finds the best candidate as l = arg max 1≤l≤L R(φ l ) = arg max 1≤l≤LR (X l ) for its output. The sum rateR(X l ) with the best candidateX l is a non-decreasing function of the number L of random samples so thatR(X l ) is likely to move closer to the optimal sum rateR(X o ) as L increases. Later, we will empirically demonstrate that Gaussian randomization with a sufficient L provides a good approximate solution close to the optimal valueR(X o ) of our problem, as in the other SDR applications [31]. (22) consists of X with l = (N + 1)(N + 2)/2 complex variables, whilst Y ∈ S M + and F ∈ S N+1 + are in the constraints. Thus, the problem can be solved by an interior-point algorithm with O((M 2 + (N + 1) 2 )l 2 ) operations per search point and a worstcase complexity of O( √ N + 1) iterations [33]. In short, the complexity solving (22) is given by O(N 6.5 ), which becomes unacceptably large as N increases.

Gradient-Based Iterative Algorithms for a Large-Sized IRS
This subsection provides a suboptimal approach to solving problems (19) in a faster way for a large N by transforming the problem into an equivalent unconstrained nonlinear optimization problem. For this purpose, we remove the phase constraints in (19) as which has the same optimal value with (19) due to the periodicity of R(φ) for all the entries of φ with period 2π. For ease of exposition, we convert (23) to a minimization problem as min φ∈R N f (φ) (24) by defining f (φ) = −R(φ). Due to the nonconvexity of R(φ) and its complicated form, it is almost impossible to obtain the optimal solution of (24) even for unit modulus reflection. Instead of the sequential optimization optimizing one IRS element at a time [7,30], we solve the problem through an iterative algorithm minimizing a local approximation (or a surrogate function) at each iteration to update its IRS phases φ t for t = 0, 1, · · · simultaneously [34]. To accommodate a large N, we adopt the algorithms based on second-order Taylor series approximations but relying on the gradient ∇ , · · · , ∂ f ∂φ N ] T in their implementation. The algorithms include the gradient descent (GD) [5,34] and quasi Newton methods such as BFGS and L-BFGS [34,35], which are briefly summarized in the following.
The algorithms are based on second-order Taylor approximations that can be expressed in a generic form as where ∆φ t = φ − φ t , g t = ∇ f (φ)| φ=φ t , and A t ∈ C N×N is chosen by an algorithm. The solution is updated by minimizing f A (φ, φ t ) as The GD with the update rule for δ t > 0 at complexity O(N) is obtained with a choice of A t = 1 δ t I M , where the step size δ t is determined by the Armijo rule [5]. The BFGS and L-BFGS update the search point as where Q t is an estimate of the inverse Hessian A −1 t with A t = ∇ 2 f (φ)| φ=φ t to reduce the complexity of the Newton method computing the Hessian and its inverse. The BFGS estimates Q t+1 with Q t , u t = φ t+1 − φ t , and r t = g t+1 − g t at complexity O(N 2 ) [34] whilst the L-BFGS having m B memories estimates Q t+1 with {u i , r i } t i=t−m B +1 at complexity O(m B N) for [35].

Simulation Results
The performance of the IRS-aided uplink NOMA is evaluated when the maximum transmit power is set to P max k = 23 dBm for k ∈ K and the noise power is set to σ 2 = −100 dBm. Practical IRS reflection is modeled with α = 1.6, β min = 0.2, and φ 0 = 0.43π as in Ref. [32]. We set the tolerance of the algorithms to 10 −5 up to maximum 500 iterations. The simulation setup is illustrated in Figure 2, where the (x, y, z) coordinates are given in meter. The BS and IRS are located at (0, 0, 10) and (50, 50, 10), respectively, whilst the devices are uniformly distributed in the shaded rectangular region bounded by (100, −20, 0), (100, 20, 0), (250, −20, 0), and (250, 20, 0). The channels are modelled as where the subscripts LoS and NLoS represent the line-of-sight (LoS) and non-LoS (NLoS) components, respectively, and ω x,y and κ x,y denote the path loss and Rician factor between nodes x and y for x, y ∈ K ∪ {R, B} with R for the IRS and B for the BS. The path loss is given by ω x,y = 10 −3 d −ν x,y x,y at distance d x,y with path loss exponent ν x,y , where ν R,B = 2.2, ν k,R = 2.8, and ν k,B = 4. The Rician factor is set to κ R,B = κ k,R = 2 and κ k,B = 0. We model the LoS components with a uniform linear array for the BS and an N v × N h uniform planar array constructed by N v = 8 and N h = N/8 for the IRS as [7] where a B (ϕ) ∈ C M and a R (ϕ, ϑ) ∈ C N are the array response at the BS and IRS, respectively, defined in [7] with the azimuth (elevation) angle-of-arrival ϕ A x,y (ϑ A x,y ) and angleof-departure ϕ D x,y (ϑ D x,y ) from node x to node y. Specifically, for the antenna and IRS with half-wavelength element spacing, we have for ϕ ∈ [0, 2π) and ϑ ∈ [−π/2, π/2). The NLoS components are modeled to be uncorrelated complex Gaussian.
The average sum rate of the network is shown as the number N of IRS elements increases in Figure 3 when M = 2 and K = 4. We provide the results with unit modulus reflection and practical reflection drawn with solid and dashed lines, respectively. Bound represents the upper boundR =R(X ) obtained with the solutionX of (22), and Maxdet denotes the performance of a feasible solution φ l derived fromX through Gaussian randomization with L = 50. Bound and Max-det are shown up to N = 128 due to the formidable computational time in solving (22) for a large N. BFGS, L-BFGS with m B = 10, and GD implemented with the derived gradient are compared with ConvSeq denoting the sequential optimization in Ref. [7]. The results with random IRS phases denoted by Random are also added to serve as a lower bound. Clearly, the average sum rate increases steeply with N by optimizing the IRS reflection. In cases of unit modulus reflection, Max-det provides the best performance close to the optimal one estimated by Bound for a moderate N. The iterative algorithms provide almost the same performance under unit modulus reflection, but BFGS, L-BFGS, and GD provide a slight gain over ConvSeq increasing with N under practical reflection. However, the gradient-based iterative algorithms are observed to reduce the computational time, as shown in Figure 4, which provides the average evaluation time per sample in obtaining the results of unit modulus reflection in Figure 3 with Intel(R) Xeon(R) Gold 6226R CPU @ 2.90Hz. Clearly, the iterative algorithms exhibit a significant reduction in computational time over Max-det at the cost of a performance loss. The computational time of GD is comparable to that of ConvSeq since both GD and ConvSeq require a large number of iterations until their convergence. L-BFGS provides the best computational time, with less complexity in updating the solution than BFGS and with a smaller number of iterations than GD by finding a better search point through inverse Hessian estimation. Thus, L-BFGS with the derived gradient would be a choice of practical merit for a large N.   Figure 5a and N = 256 in Figure 5b. The average sum rate increases as the number K of devices increases. In addition, the gradient-based methods and ConvSeq provide a similar performance except for practical reflection with N = 256. When N = 64 in Figure 5a, the performance of the gradient-based algorithms becomes close to the optimal one as K increases for unit modulus reflection. When N = 256 in Figure 5b, the gradient-based algorithms outperform ConvSeq by about 0.3 dB for practical reflection. Again, L-BFGS provides a good performance among the gradient-based algorithms. Hence, we compare the performance of L-BFGS and ConvSeq for different numbers M of antennas with N = 256 in Figure 6; we set M = 2 in Figure 6a and M = 4 in Figure 6b. The sum rate becomes almost doubled by doubling the number of antennas. Again, L-BFGS provides a similar or slightly improved performance compared with ConvSeq, which is obtained at a computational time about 11% and 18% of ConvSeq with M = 2 and 4, respectively, for most K values in the figures.
Discussion : It is noteworthy that the nonlinear optimization problem in (24) with respect to the IRS phase vector φ is a non-convex optimization problem that exhibits multiple local minima in general. Hence, the iterative algorithms considered herein do not guarantee a convergence to the global optimal point, but to one of the local minima. The gradient-based algorithms updating N variables simultaneously for the next search point tend to find a similar local minimum at a different convergence rate. However, ConvSeq, updating one variable at a time, tends to find a worse point in particular for practical reflection since it resorts to limited information for the next search point. To improve the performance, we may run an iterative algorithm with different initial points, resulting in different local minima so that a better solution can be found. However, it is observed that the gain is trivial for this problem. From this, devising a new algorithm filling the gap to the optimal performance with a complexity between those of Max-det and gradient-based algorithms would be an interesting topic for further study.

Concluding Remarks
We have considered a sum rate maximization problem for the IRS-aided uplink multiantenna NOMA under a generalized reflection model including unit and phase-dependent amplitudes. We have solved the problem through extended SDR to obtain an upper bound on the sum rate and a near-optimal solution for a moderate-sized IRS. We have applied the gradient-based iterative algorithms for a large-sized IRS by providing the gradient in an explicit form under generalized reflection. The results show that, among the gradient-based algorithms, L-BFGS implemented with the derived gradient provides a more competitive solution than the conventional method in both computational time and performance.