Joint Clustering and Resource Allocation Optimization in Ultra-Dense Networks with Multiple Drones as Small Cells Using Game Theory

In this study, we consider the combination of clustering and resource allocation based on game theory in ultra-dense networks that consist of multiple macrocells using massive multiple-input multiple-output and a vast number of randomly distributed drones serving as small-cell base stations. In particular, to mitigate the intercell interference, we propose a coalition game for clustering small cells, with the utility function being the ratio of signal to interference. Then, the optimization problem of resource allocation is divided into two subproblems: subchannel allocation and power allocation. We use the Hungarian method, which is efficient for solving binary optimization problems, to assign the subchannels to users in each cluster of small cells. Additionally, a centralized algorithm with low computational complexity and a distributed algorithm based on the Stackelberg game are provided to maximize the network energy efficiency (EE). The numerical results demonstrate that the game-based method outperforms the centralized method in terms of execution time in small cells and is better than traditional clustering in terms of EE.


Introduction
The global mobile data traffic is witnessing an exponential increase, and the number of global mobile connected devices is forecast to increase rapidly from 8.8 billion in 2018 to 13.1 billion in 2023 [1]. The deployment of ultra-dense networks (UDNs) as a type of cellular network is an extremely promising technology for supporting the huge number of connections, especially in the applications of massive Internet of Things (IoT). In comparison with present network technologies, UDNs help to considerably improve the network capacity, spectral efficiency, as well as coverage with a high quality-of-service (QoS) by reducing the distance between the base stations (BSs) and user equipment (UE) [2,3]. However, the highly dense deployment of small cells in UDNs causes high interference, high power consumption, and high computational complexity, which must be addressed to satisfy the stringent constraints of the beyond fifth-generation networks (B5G) [4].
In UDNs, small cells are smaller than ever, and therefore the distance of neighboring small cells decreases. As a result, the intercell interference is extremely high. In addition, network elements. Energy efficiency (EE), which is defined as the ratio of total data rate to total power consumption, is a useful metric for modeling the trade-off between throughput and energy. However, EE maximization problems are nonconvex optimization problems and much harder than throughput maximization problems or power minimization problems. In [24,25], two optimization problems (i.e., SCA and power allocation (PA)) were considered with PA problems modeled by the Stackelberg game (SG) in order to maximize the EE. Although there were some considerable improvements in terms of EE, there was a limitation in the applicability of the used model of UDNs, which had only one macro base station (MBS) and one BS in a network. In the previous works, either cooperation or noncooperation between network elements was considered, but these elements have the characteristics of both collaboration and competition,. In detail, small-cell base stations (SBSs) prefer to be in joint clusters for cooperatively transmitting signals to UEs. This leads to a decrease in intracluster interference. On the other hand, the SBSs in different clusters tend to conflict regarding power usage. This means that increasing the power of an SBS in a cluster also increases the interference for SUEs in adjacent clusters. Therefore, there is a lack of a game-based paradigm that has the characteristics of both cooperation and competition of network elements in order to concurrently mitigate interference, maximize the network EE in drone-aided UDNs with multiple MBSs, and handle a massive number of randomly distributed drones.
Motivated by the aforementioned discussion, this paper proposes a mathematical model to jointly design the clustering and resource allocation for maximizing the total network EE in UDNs consisting of multiple MBSs and many small-cell drones (SCDs). Cells served by MBSs are macrocells, and cells served by drones are small cells. We propose a coalition game for clustering SCDs to restrict the intercell interference. The SCDs in each cluster cooperatively transmit the signals to UEs. The SCA methods for two tiers are provided to choose the subchannel with the highest channel gain for each UE. To reduce the computational complexity, centralized and distributed PA methods are proposed in terms of maximizing the EE. The main contributions of our paper are summarized as follows: • Adjacent SCDs tend to join together and cooperatively serve UEs because they highly interfere with each other. Thus, we form this cooperation with a coalition game, with the utility function being the signal-to-interference ratio (SIR). • To fully eliminate the intracluster interference, we model the SCA for small cells as a binary optimization problem. Then, the Hungarian method is proposed to obtain the solution. Meanwhile, the SCA for macrocells chooses the subchannel with the highest channel gain without any constraint. • To reduce the complexity of nonconvex PA optimization problems, we propose an iterative centralized PA algorithm to obtain at least a locally optimal solution in terms of maximizing the EE. Inequalities are used for relaxing the complex objective functions of EE maximization into convex functions that are easy to solve with programming tools. • The distributed PA algorithm based on SG consists of two subcooperative games with MBSs as leaders and SCDs as followers. This algorithm decouples the very complex PA optimization problems into multiple low-complextiy and convex ones that are extremely useful for solving optimization problems with a huge number of variables.
Compared with the conference version [26], a macrocell layer consisting of multiple MBSs is added to the networks. As such, the interference environment is more complex due to intercell interference. Additionally, the objective function is to minimize the network EE, which is more practical to use energy more efficiently compared with the objective of maximizing the sum-rate in UDNs but is more complex to solve because of the fractionaltype function. Moreover, we propose a centralized approach to solve PA as a benchmark and the distributed one using the Stackelberg game, and multiple simulations with different aspects such as EE, sum rate, and execution time were performed to prove the efficiency of the proposed methods. The rest of this paper is organized as follows: Section 2 presents the system model and optimization problems used in this paper. A small-cell clustering method based on the coalition game is described in Section 3. Then, the SCA for both macrocells and small cells as well as the beamforming method for MBSs are discussed in Section 4. In addition, Section 5 describes both centralized and distributed optimization methods for solving PA problems. Numerical simulations proving the efficiency of our proposed methods are outlined in Section 6. Section 7 provides the conclusions of this study.

System Model
In this paper, we consider the downlink transmission in a two-tier UDN that consists of multiple macrocells and a massive number of UAVs serving as small-cell base stations with highly dense deployments. MBSs and SCDs share the same N orthogonal subchannels with bandwidth BW per subchannel to serve their UEs. The advantage of using multiple orthogonal subchannels is that multiple UEs can be simultaneously served in these subchannels without interference. One typical macrocell with two subchannel is illustrated in Figure 1. In the macrocell tier, each F MBS is equipped with T array antennas to transmit the signal to M MUEs by using massive multiple-input multiple-output (mMIMO) techniques. In the small-cell tier, to serve U UEs, many B SCDs equipped with one omni-directional antenna each are randomly deployed according to the independent homogeneous Poisson point processes (PPP) Φ with density λ in the coverage of MBSs (B > U) [24]. When UEs are located in the coverage of one or some SCDs, they are served by SCDs, i.e., SUEs. On the other hand, if the UEs are not in the coverage of any SCD, then the power of SCDs is not high enough to serve these UEs or causes extremely high interference with surrounding SCDs, especially in UDNs. Therefore, they are served by MBSs, i.e., MUEs. We used the main mathematical notation, which is shown in Table 1, to design tbe optimization problems. At each MBS, because of the mMIMO, the number of antennas of MBS f is much greater than the number of its authorized MUEs (T M f ). Each signal s (n) m, f is multiplied by a beamforming vector w m, f ∈ C T×1 , and then the aggregation of processed signals of MUEs that are served by MBS f is transmitted.

Notation Definition
The interference power from other clusters except cluster C i The aggregation x (n) f at MBS f in subchannel n is the sum of M f processed signals: where the binary indicator of SCA is defined as The received signal y (n) m, f at MUE m served by MBS f in subchannel n is expressed by interference from other MBSs interference from clusters +n m , (2) where the aggregation of signals at MBS f is x The received signal y (n) u,C i of SUE u served by cluster C i in subchannel n is given by interference from other clusters u,b is a complex Gaussian random variable with zero mean and unit variance that represents the small-scale fading of the channel from SCD b to SUE u, and β (n) u,b is the large-scale coefficient of this channel. Therefore, the SINRs at SUE u and MUE m are, respectively, written as

Problem Formulation
As the deployment of small cells in UDNs becomes extremely dense, the problem of effectively using resources becomes more important than ever. In this section, we formulate the optimization problems for the resource allocation of both spectrum and power according to the processing operations in Figure 2. The overall architectures of the proposed optimization solutions are shown in Figure 2 with the names of the 6 problems in blocks. The detailed description and algorithms are presented in Sections 3-5. The objective functions at different stages are different to obtain an optimal solution for multiple proposed objectives, such as interference mitigation in clustering, a guarantee of high-gain channels in subchannel allocation, and efficient energy usage in power allocation.  Due to the high efficiency of spectral usage, co-channel assignment is usually used in practical systems [27]. However, this causes high interference between the two tiers in UDNs and is hard to control, especially with the random deployment of small cells. Thus, the optimization problems with the objective to choose the highest channel gain guarantee that the BSs serve their UEs in good subchannels.

SCA for macrocells:
We consider MBS f serving M f MUEs. By using mMIMO technology at each MBS, the number of MUEs in each subchannel causes a minor interference effect. The target is to find a subchannel with the highest gain for each MUE. The optimization problem can be expressed as SCA for small cells: In Figure 2, after clustering the small cells, the SCDs in each cluster cooperatively transmit the same signal to their SUEs. Therefore, the power of the desired signal at each SUE is improved compared with that using only one serving SCD. Considering cluster C i , to find the appropriate subchannels for the SUEs in set U C i , an optimization problem is formulated as where VC u,n = ∑ b∈C i |h (n) u,b | 2 is the sum of channel gains from SCDs in cluster C i to SUE u. (7b), (7c) indicates that each subchannel n is authorized by at most one SUE u ∈ U C i and vice versa. (7d) is the constraint of the binary values of SCA. If a (n) u,C i = 1, then SUE u is served by the SCDs in cluster C i through subchannel n, and vice versa if a Optimization problem (7) is denoted by (P2) in Figure 2. Problem (7) is integer linear programming because (7a)-(7c) are affine functions, and (7d) is the constraint of two discrete values of variables. The solving methods for (P2) and (P3) are proposed in Section 4.1.

Power Allocation (PA)
In UDNs, the distance between the SCDs and SUEs as well as the distance between small cells are very small: from a few dozen to a few hundred meters. Therefore, increasing the power at an SCD to improve the signal to its UEs causes extremely high interference to other small cells compared with that of conventional networks. In addition, different tiers have different objectives that need to be optimized. For that reason, after SCA, we formulate two optimization problems for PA in macrocells and small cells ( Figure 2).
PA for MBSs (EE maximization): Each MBS, equipped with massive antennas, is connected directly to a core network. MBSs, using high power with wide coverage areas, guarantee the serving of all UEs that are not connected to any small cell. Therefore, an important problem for MBSs is maximizing the data rate of the MUEs to satisfy the standards of the new generations of future networks (B5G) while efficiently using energy. EE, which is the ratio of data rate to power consumption, is widely used for evaluating the efficiency of energy usage. The optimization problem of PA for MBSs is as follows where the variable p (n) m, f represents the transmit power of MBS f to MUE m in subchannel n. P total m, f is the sum of power fo ther hardware circuits and power for transmitting signal in all subchannels of MBS f , P max MBS is the maximum power of each MBS, and T pm min is the minimum throughput of each MUE. (8b) and (8c) are the constraints of the transmit power of MBSs, which range from 0 to P max MBS . (8d) guarantees the QoS of all UEs.
PA for SCDs (EE maximization): Due to the short distance between the SCDs and SUEs, the channel links experience the line of sight (LoS). However, with a large number of small cells, the power consumption rises considerably. Thus, the target of SCDs is effectively using limited energy. The optimization problem for maximizing the EE of small cells is formulated as where p u,b is the sum of power for hardware circuits and power for transmitting signal in all subchannels of SCD b ∈ C i , and P max SCD is the maximum transmit power of each SCD, T ps min is the minimum throughput of each MUE. (9b) and (9c) represent the range of transmit power. To guarantee the QoS, (9d) is used for all SUEs.

SCDs Clustering by Coalition Game (P1)
In this section, we introduce a coalition game for clustering small cells. In UDNs, an effective clustering method is essential because of three main reasons, as follows: • Mitigating interference: Intercell interference is a big challenge in UDNs due to the ultra-dense density of small cells, especially for SUEs in the overlapping areas of small cells. Formulating clusters that combine small cells seriously interfering with each other is an appropriate solution. The SCDs in the same cluster cooperatively transmit the signal to UEs so as there is no intracluster interference. • Enhancing the quality of the signal: In UDNs, the SUEs in the overlapping areas of two or more small cells are usually far from their primal SCDs, while they receive high interference from the others. In many cases, using some of these SCDs along with primary SCDs to serve UEs helps improve the power of the desired signal. Furthermore, multiple SCDs transmitting the same signal to a UE also improves the diversity gain. • Decreasing handover processes: The smaller the size of small cells, the larger the number of handover processes. Moreover, when UEs move at high speed, the central controller has to compute the complex handover processes with the massive number of UEs in a short time. Therefore, if small cells join together to form temporarily bigger cells, the number of handover processes significantly declines.
Therefore, adjacent small cells have an incentive to join clusters to obtain a higher value for the objective instead of independent transmission. The coalition game. as a type of cooperative game where players tend to cooperate to achieve higher total utility, is appropriate to model this situation. The standard form of a coalition game is CG = PL, v , where PL = B active is the set of players (i.e., SCDs, which serve at least one SUEs), and v(C) ⊆ R C is the set of nontransferable utility functions of clusters [28]. The utility function of cluster C i is defined as where SIR u is the SIR of the links from cluster C i to SUE u and is expressed as where The solution of this coalition game is the set of coalitions of players or clusters of SCDs in UDNs [28]. In these clusters, there is no SCD that has an incentive to change clusters for achieving higher utility. This is the state where all clusters have reached equilibrium. The two operations that change the SCDs between clusters are split-merge and swap [29]. If SCD b ∈ C j leaves C j to join C i , then the operation is denoted as In addition, if SCD b ∈ C j swaps with b ∈ C i , then the operation is described as To obtain the solution to this coalition game, we use Algorithm 1. In this algorithm, we use the variable d max for determining the SUEs that are less than d max m away from SCD b. Then, the set of neighbor clusters C near that consists of these SUEs is determined. Therefore, different from [30], we only consider some close clusters with SCD b because the interference from b has a minor impact on far SUEs at very high frequency. In addition, variables g and g n are used for choosing the best option for SCD b to move to a new cluster. If the number of SCDs in each cluster increases, then the number of pairs of SCD and SUE with a distance larger than d max increases. Thus, we use N max for limiting the number of SCDs in each cluster.
Algorithm 1 Clustering SCDs using the coalition game. Assume b ∈ C j , C near is the set of clusters that have some SUEs covered with radius d max and center b, gap g = 0 7: if g n > g 15:

Subchannel Allocation (P2, P3) and Beamforming Design for MBSs (P4)
In this section, we propose methods to solve the SCA problems (6) and (7). Satisfactory subchannels are allocated to UEs with high channel gain while guaranteeing the constraints are met. In addition, the beamforming design for MBSs is also described.

Subchannel Allocation
(P2) For clusters, constraints (7b)-(7d) are the same as the constraints of an assignment problem where each task is assigned to only one worker to minimize the cost. One of the effective methods to solve assignment problems is the Hungarian method [31]. Compared with other methods, the Hungarian algorithm can solve large-scale assignment problems in polynomial time. Addtionally, it is a deterministic and global algorithm that guarantees the optimal solution is the best assignment. Regarding flexibility, the Hungarian method can be applied to solve assignment problems with equal or unequal numbers of resources and tasks, especially in the cases of variable numbers of network elements in clusters in our proposed SCA problem. Therefore, we use the Hungarian algorithm for making decisions in choosing the optimal subchannels for SUEs. Matrix VC maxprob , which consists of the value of VC u,n with all pairs of SUE u ∈ C i and subchannel n, is given in Table 2.
To obtain the solution to the optimization problem (7), we use the Hungarian method in Algorithm 2 [30]. Because of the limitation of the number of UEs in each cluster by N, the complexity of Algorithm 2 is O(N 4 ) [31]. Cover all zeros in the matrix with minimum number of vertical/horizontal lines nl min . 9: if nl min < N then 10: Subtract the lowest number from all elements that are not covered by any line. 11: Add this lowest number to elements that are crossed by any two lines. 12: Until nl min = N. 13: Obtain the optimal association from zeros in the processed matrix. 14: Output: The association matrix X = {x i,k | i = 1, . . . , , U C i and k = 1, . . . , N}.
(P3) For macrocells, constraint (6b) requires one MUE to be allocated only one subchannel. Moreover, the number of MUEs in a subchannel is not important due to the mMIMO MBSs. Therefore, the solution to the optimization problem (6) can be easily obtained by choosing the maximum value of VM m,n , n = 1, . . . , N for each MUE m.

Beamforming Design
According to the model in Figure 2, after SCA, we need to design the beamforming vectors for MBSs. Assuming that the CSI is known at the MBSs, we introduce two types of appropriate beamforming vectors for two cases.

•
If there are multiple MUEs that are served by MBS f in subchannel n, then the beamforming vector is one vector in the null space of the channel matrix from MBS f to all of its MUEs except MUE m, and is expressed as Therefore, we can neglect the cross-talk interference in (2). • If only one MUE m is served by MBS f in subchannel n, then according to the zeroforcing beamforming, the beamforming vector is as follows In both cases, the magnitude of the beamforming vector is normalized to 1 (||w m, f || = 1).
Thus, the PA for MBSs and SCDs only depends on p

Power Allocation for Macrocells and Small Cells (P5, P6)
Power allocation (PA) is a promising approach for interference mitigation in UDNs. In this section, we propose both a distributed method and a centralized method to obtain the solutions to the optimization problems (8) and (9). It is clear that the powers of MBSs and SCDs are variables in both optimization problems (8) and (9). This causes these optimization problems to be extremely complex and hard to solve. In addition, the channels from MBSs to MUEs are usually more strongly affected by environmental factors than the channel between SCDs and SUEs, while guaranteeing the QoS. Therefore, the optimization problem of MBSs has higher priority than that of SCDs. We assume that the interference from small-cell tiers to each MUE m is a constant at first as follows: Each active SCD uses the maximum power of P max SCD , which is divided equally to th subchannels, to serve theSUEs. Thus, the interference from the small-cell tier to MUEs is extremely high when first solving the PA optimization problem of MBSs. This helps the QoS constraints to remain satisfied after obtaining the solution to the PA optimization problem of SCDs with very low interference.

Centralized Power Allocation
In this subsection, we propose practical methods to solve centralized PA problems for macrocells and small cells. In these methods, we combine three techniques: changing variables, approximating objective functions, and the Dinkelbach algorithm. To this end, the very complex maximizing EE optimization problems especially in UDNs are converted into much easier convex problems in iterative algorithms. This is a promising approach to achieving real-time computing in multitier massive-cell networks such as UDNs.
The EE maximization PA problem for macrocells is expressed as follows where γ m min 2 T pm min /BW − 1 is the minimum SINR to guarantee the QoS of MUEs. By using the beamforming vectors in Section 4.2, the value of cross-talk interference I (n) cross-talk in the macrocells equals 0. I (n) m,C is fixed by using estimation (12). We can easily convert constraint (13d) into a linear inequality becuase each γ m in (13d) is an affineaffine fraction. Therefore, all constraints in the optimization problem (13) are linear. Let F , is a feasible point of optimization problem (13). Following the method in [32] and using in- According to the Dinkelbach algorithm, the suboptimization problem to find P (κ+1) F for iteration (κ + 1) is as follows Becauseã m > 0, the objective function in (14), which is the sum of pairs of convex functions, i.e., constant × variable and −constant/variable (constant > 0), is convex. Therefore, problem (14) is convex due to the convexity of the objective function and the affine functions at the constraints. The proof of the convergence is provided in Appendix B.1. Based on [33], the computational complexity of each iteration for solving It is given that the UE u is served by cluster C i , which combines SCDs b. The EE maximization PA problem for small cells is given by where γ s min 2 T ps min /BW − 1 is the minimum SINR that guarantees the QoS of the SUEs. Using the same technique as for the PA problem for macrocells, (15d) can be easily converted into a linear constraint. Therefore, the optimization problem (15) is a linear-constrained optimization problem. However, the objective function of the optimization problem (15) is more challenging than the optimization problem (13) of macrocells because there is more than one variable in the numerator of γ u . Let I u (P C ) = I C is a feasible point of the optimization problem (15). Following inequality (A3), we have . According to the Dinkelbach algorithm, the suboptimization problem to find P (κ+1) C for iteration (κ + 1) is as follows The same as for optimization problem (14), optimization problem (16) is also a convex problem. The proof of convergence is presented in Appendix B.2. Figure 3 shows a diagram for solving the centralized PA optimization problems, with denoting the tolerance value and n max being the maximum number of iterations. If the values of the objective functions in both (14) and (16) are lower than the tolerance value or n l p ≥ n max , with n l p denoting the number of iterations, then the convergence holds true. With the number of decision variables BUN and the number of constraints B + BUN + U, the computational complexity of each iteration for solving (14) is O (BUN) 2 (B + BUN + U) 2.5 + (B + BUN + U) 3.5 .

Distributed Power Allocation by the Stackelberg Game
The solution of centralized PA methods is the best option for PA to MBSs and SCDs. However, when the number of small cells increases, the execution time of these methods also significantly rises. Therefore, centralized PA methods cannot satisfy the strict time constraints of UDNs in B5G. In this subsection, we propose a distributed PA method based on a SG to considerably reduce the complexity of the original optimization problems while retaining an acceptable efficiency of the solution. A SG is a type of NCG where players have no incentive to collaborate. There are two levels of priority in a SG. In detail, players with a high level of priority are called leaders and the others are followers. When starting a SG, leaders choose optimal actions first. Then, followers see the actions of the leaders to choose their optimal actions. Thus, a SG is also called a two-stage NCG.
The distance from MBSs to MUEs is usually much larger than the distance between SCDs to SUEs. Therefore, the problems of reflection, diffraction, and scattering influence the received signal at MUEs more than at SUEs. Meanwhile, MBSs must guarantee serving all UEs that are out of the coverage of the small cells, with the constraint of data rate. Therefore, MBSs, as leaders, have a high level of priority, and SCDs are followers in this SG.

NCG for Leaders (MBSs)
We formulate the optimization problem (8) as an NCG, where NCG leader = PL F , P F , ul F , where PL F = F is the set of players (i.e., MBSs), P F is a matrix of PA for all MBSs, and ul F ( f ) is the utility function of the EE of MUEs served by MBS f . The solution to NCG leader is a Nash equilibrium where no MBS has an incentive to change its power, and the power of each MBS is the best response to the power of other MBSs [30]. By using an NCG, the optimization problem (8) can be divided into multiple suboptimization problems for MBSs to obtain the best response for each player. In each subOP, we find an optimal solution for an MBS while the power of the other MBSs is fixed. A subOP for MBS f is formulated as The objective function (17a) is a concave-convex fraction. We can use the Dinkelbach method to solve the optimization problem (17). In detail, we introduce a new variable λ f , and the optimization problem (17) can be rewritten as The objective function in (18) is a concave function in the maximization problem, so the optimization problem (18) is a convex problem. Therefore, the optimization problem (18) can be effectively solved by using CVX or CVXPY [34,35]. The optimal solution P * f is also the optimal solution to optimization problem (18) We use Algorithm 3 to find the solution for the NCG NCG leader . With the number of decision variables M f N and the number of constraints 1 + M f N + M f , the computational complexity of each iteration for solving (18) It is clear that the complexity of (18) is much lower than that of (14) because M f < M and only one MBS is considered.

NCG for Followers (Clusters)
We formulate optimization problem (9) as an NCG NCG f ollower = PL C , P C , ul C , where PL C = C is the set of players (i.e., clusters), P C is a matrix of PA for all clusters of SCDs, and ul C (C i ) is the utility function of the EE of cluster C i . The PA problem for SCDs in (9) can be divided into multiple suboptimization problems by using NCG with players as clusters. A suboptimization problem of cluster C i is as follows To obtain the best response of cluster C i , the power of the other clusters is fixed. This best response is the solution to optimization problem (19). We use the Dinkelbach method to solve the optimization problem (19) because its objective function is a concave-convex fraction. In detail, we introduce a new variable λ C i , and the optimization problem (19) can be rewritten as max The objective function in (20) is a concave function, so the optimization problem (20) is a convex optimization problem. Therefore, the optimization problem (20) can be effectively solved by using CVX or CVXPY [34,35]. The optimal solution P * C i is also the optimal solution to the optimization problem (20) if UT C i (P * C i , λ * C i ) = 0. We use Algorithm 3 to find the solution to the NCG NCG f ollower . With the number of decision variables U C i C i N and the number of constraints C i + U C i C i N + U C i , the computational complexity of each iteration for solving (20) The complexity of (20) is much lower than that of (16) since U C i < U and C i < B.

Simulation Results
Next, we investigated the performance of the proposed methods in a two-tier UDN. Three hexagon macrocells with a coverage radius of 100 m were deployed next to each other. In each macrocell, one MBS equipped with 128 antennas was located in the center of the cell. In addition, many small cells that had a radius of 20 m were randomly distributed by PPP. A small cell had an SCD equipped with one omnidirectional antenna. The UEs in the coverage of any small cell were served by the cluster that contained that small cell, and the others were served by theMBSs. The parameters in the simulations are given in Table 3 [36]. Additionally, we introduced two more clustering methods for comparison. In the traditional clustering method, each SUE is served by the nearest SCD. Meanwhile, with random clustering, each SUE randomly connected to an SCD, which had a distance to this SUE of less than d = 3 × the radius of small cells. On the other hand, equal power allocation, where the power of P max MBS /2 of each MBS is equally divided among its MUEs and the power P max SCD of each active SCD is equally divided among the subchannels, was also used for comparison with the EE maximization PA methods. To prove the performance improvement, the proposed methods were compared with conventional methods based on multiple aspects such as convergence speed of PA methods, energy efficiency, sum rate, and transmit power usage. The acronyms for the combinations of clustering and PA methods are given in Table 4. To solve the optimization problems and build simulating programs, CVXPY and PYTHON were used [35]. The computing platform was a PC with CPU@3.7 GHz and 32 GB RAM memory.  Figure 4 illustrates an example of a triple-macrocell UDN that combines 100 SUEs and 24 MUEs with a density of SCDs of 2000 SCDs/km 2 . In the first close-up in Figure 4a, two adjacent small cells use two same subchannels to serve their SUEs. This causes high interference for these SUEs. Meanwhile, after using the coalition game, these two small cells form a cluster to cooperatively transmit (the first close-up in Figure 4b). Therefore, there is no interference in a cluster. Another advantage of clustering that is shown in the second close-up is to enhance the quality of the desired signal. The SUE (red) subchannel is served by two small cells tgat have almost the same distance to it instead of being served by a single cell. When the density of SCDs is very high, these situations occur many times. Thus, clustering SCDs based on coalition game, and cooperative transmission help to efficiently mitigate intercell interference and improve the power of the received signal.

Convergence Speed of PA Methods
In this subsection, we consider the convergence of CGCO and CGDO in the network that has an SCD density of 2000 SCDs/km 2 , 100 SUEs, 3 MBSs, and 24 MUEs. The changes and differences in the quantity of the EE of these methods, along with the iterations in the macrocell tier and small-cell tier, are shown in Figure 5a,b, respectively. It is given that one iteration of GT methods is complete after updating the best responses of all players. At the convergence, the EE of the CGDO is always lower than that of the CGCO. In detail, the EE of the CGDO is lower than that of CGCO in the macrocell tier and small-cell tier by 0.65% and 15.39%, respectively. However, the CGDO almost converges to the equilibrium after only one loop for updating all the best responses of the players, while it takes three and four iterations in the macrocell tier and small-cell tier, respectively, for the CGCO to achieve an EE of more than 95% at the convergence. On the other hand, in Figure 5b, the EE of the CGDO after one iteration is greater than the value of the convergence. This is reasonable because, in game theory, the value of the utility function at the Nash equilibrium may not be the best solution.

Performance Analysis of CGCO, CGDO, and TCDO
We evaluated the performance of the proposed methods regarding the EE, the sum rate, the total transmit power, and the execution time in different scenarios with the number of network elements that are given in Table 5, and the number of SUEs equaled 60% of the number of SCDs.  Figure 6 displays data on the changes in the quantity of the EE of CGCO, CGDO, TCDO, RCDO, and TCEP in the small-cell tier when the density of SCDs increases from 500 to 1000 SCDs/km 2 . To overview, the EE of all methods witnesses a decrease with different levels since the denser the deployment of small cells is, the higher interference is. The EE of CGCO is greatest in both the five methods since it uses the centralized PA method. In detail, the EE of CGDO equals 74.09% of the EE of CGCO on average. Meanwhile, the EE of CGDO is greater from 11.15% to 17.23% than the one of TCDO in all testing scenarios. In addition, CGCO, CGDO, and TCDO outperform RCDO and TCEP in terms of the EE. Therefore, using the coalition game for clustering small cells before SCA and PA based on distributed method helps improve the efficiency of energy usage in UDNs compared with the traditional clustering and the random clustering. Along with the EE, the computing execution time is an important factor to evaluate the performance of different methods. Table 6 describes the data in detail regarding the time required for the PA of CGCO, CGDO, and TCDO. Thanks to (12), the intensity of interference from small cells is pre-estimated in the PA optimization problems of the macrocell tier. Therefore, the execution time of the PA methods in the macrocell tier does not depend on the density of SCDs. On the other hand, the execution time for the PA of the MBSs of CGDO and TCDO is higher than the one of CGCO when the number of MUEs is small. The reason is that optimization problem (18) with the objective, a logarithmic function, takes a longer time to solve than optimization problem (16) with the objective as the sum of affine functions and rational functions when the number of variables is slightly different. However, the execution time in macrocells is less important than that in small cells with much more variables in the OPs. Considering the small-cell tier, when the density of the SCDs changes from 500 to 1000 SCDs/km 2 , the execution time of the CGCO significantly increases by 1.7 times, while those of the distributed methods CGDO and TCDO do not depend on the number of SCDs. Additionally, the execution time of CGCO is always longer, from 1.95 to 3.49 times that of CGDO. In large-scale UDNs, quick adaptation is key. If a very complex algorithm is used and the execution time is long, when the solution is obtained, it is not optimal anymore because the channel and user locations have changed too much. Therefore, the trade-off between EE and a quick execution time using our proposed PA method is essential and effective in UDNs.  Figure 7 shows the average of the sum rate and the total power consumption of CGCO, CGDO, and TCDO in the macrocell tier of UDNs. The level of the sum-rate growth of the three methods gradually decreases when the number of small cells increases in a certain area. The sum rate of these methods increases and reaches the maximum at a density of 900 SCDs/km 2 . Then, it gradually decreases when densities are greater than 900 SCDs/km 2 despite the increased power usage because MUEs experience extreme interference from the small-cell tier. To clearly identify this trend, we added two more scenarios with densities of 1100 and 1200 SCDs/km 2 , and we set the number of MUEs to 24. Additionally, in comparison with traditional clustering, SCDs were used in more subchannels for transmitting the signals to both their primal SUEs and the other SUEs in their clusters. This causes more interference from small cells to macrocells. Therefore, MBSs with CGDO use more power to achieve a higher sum rate than those of TCDO. Figure 8 describes the sum rate and the total transmit power of the three methods, namely CGCO, CGDO, and TCDO, in the small-cell tiers of UDNs with different network sizes. When the density of the SCDs increases, the energy consumption considerably increases, while the rate of increase in the sum rate declines. Due to using centralized optimization, the data rate of CGCO is always larger than those of the others, while the total transmit power of CGCO is the lowest for the three methods. The power consumption of CGDO is lower by 9.51% to 15.19% than that of TCDO. On the other hand, when the deployment of SCDs becomes denser, to achieve the same data rate, the SCDs use much more power. For example, in Figure 8a, the sum rates, which CGDO and TCDO achieve with small cells with densities of 800 and 900 SCDs/km 2 , are nearly equal. Meanwhile, the power consumption of CGDO and TCDO significantly increases in this range of density, as shown in Figure 8b. Moreover, the data rate of CGDO is slightly higher than that of TCDO. To summarize, the clustering method based on the coalition game helps to achieve a higher data rate than the traditional clustering method with lower power consumption.

Execution Time of Coalition Game (Algorithm 1) and the Hungarian Method (Algorithm 2)
Using the same simulation parameters as in Section 6.2, we evaluate the execution time of Algorithms 1 and 2 in Table 7. It is clear that the execution times of the coalition game and the Hungarian method increase with the increasing density of SCDs for the different levels. In detail, from 500 to 1000 SCDs/km 2 , the execution time of the coalition game increases 3.8 times and that of the Hungarian method increases 3.1 times. The computational complexity of the Hungarian method is O(N 4 ) in the worst case for the SCA in each cluster. Thus, its computational time has a linear relationship with the number of clusters in the small-cell tier.

Conclusions
We investigated the problem of resource allocation for UDNs with two tiers consisting of multiple MBSs and massive randomly distributed drones in this study. We proposed a paradigm that combines three processes: clustering using the coalition game, SCA using the Hungarian method, and PA using the SG as a distributed optimization method. Our simulations proved that the execution time of the centralized method for PA is much higher, by 1.95 to 3.49 times, than that of our distributed method, while the network EE of the proposed method using clustering based on the coalition game is improved by 11.15% to 17.23% compared with that of the distributed one with traditional clustering. Therefore, the centralized method is more suitable for UDNs with a small number of network elements, while the distributed method can be efficiently employed in large-scale UDNs where many MBSs and drones communicate simultaneously. Moreover, real-time computing plays an important role in drone-aided UDNs due to the mobility of both base stations (i.e., drones) and users. The execution time of the coalition game and the Hungarian method showed that they can be applied to real-time systems in UDNs.
For future work, with a massive number of users, data traffic is extremely high in drone-aided networks. To avoid network data congestion, drones equipped with caching storage, which can prestore common data packets, are a promising solution. However, the optimization problem for caching is integer programming with a large number of decision variables. Therefore, low-complexity algorithms for solving caching problems with the objective of minimizing the risk of network congestion are essential to be considered as an extension of this method. Additionally, to guarantee the quality of the signal received by users, drone-to-user communications need to be line of sight. Thus, the trajectory design for multiple drones in order to adapt to the movements of users is also an open topic to be addressed in the future.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.

Conflicts of Interest:
The authors declare no conflicts of interest.

Appendix B. Proving the Convergence
Appendix B.1 Let P (κ+1) F be the solution of the convex optimization problem (14). Thus, the value of objective function in (14) at P (κ+1) F is greater than the value of the other feasible points (i.e., P (κ) F ) as [32] ϕ M (P (κ+1) F , P (κ) From (A1), the equality is achieved when x =x and y =ȳ. Therefore, we havẽ With (A4) and (A5), we have Therefore, the first loop in Figure 3 converges at least to a local maximal point of the optimization problem (13).