# Distributed Collaborative 3D-Deployment of UAV Base Stations for On-Demand Coverage

Deployment of unmanned aerial vehicles (UAVs) performing as flying aerial base stations (BSs) has a great potential of adaptively serving ground users during temporary events, such as major disasters and massive events. However, planning an efficient, dynamic, and 3D deployment of UAVs in adaptation to dynamically and spatially varying ground users is a highly complicated problem due to the complexity in air-to-ground channels and interference among UAVs. In this paper, we propose a novel distributed 3D deployment method for UAV-BSs in a downlink network for on-demand coverage. Our method consists mainly of the following two parts: sensing-aided crowd density estimation and distributed push-sum algorithm. The first part estimates the ground user density from its observation through on-ground sensors, thereby allowing us to avoid the computationally intensive process of obtaining the positions of all the ground users. On the basis of the estimated user density, in the second part, each UAV dynamically updates its 3D position in collaboration with its neighboring UAVs for maximizing the total coverage. We prove the convergence of our distributed algorithm by employing a distributed push-sum algorithm framework. Simulation results demonstrate that our method can improve the overall coverage with a limited number of ground sensors. We also demonstrate that our method can be applied to a dynamic network in which the density of ground users varies temporally.

## Authors

• 6 publications
• 2 publications
06/12/2018

### Downlink Analysis in Unmanned Aerial Vehicle (UAV) Assisted Cellular Networks with Clustered Users

The use of unmanned aerial vehicles (UAVs) operating as aerial base stat...
04/22/2018

### Cellular Connectivity for UAVs: Network Modeling, Performance Analysis and Design Guidelines

The growing use of aerial user equipments (UEs) in various applications ...
06/01/2021

### Energy-aware placement optimization of UAV base stations via decentralized multi-agent Q-learning

Unmanned aerial vehicles serving as aerial base stations (UAV-BSs) can b...
12/04/2018

### Coverage and Rate Analysis for Unmanned Aerial Vehicle Base Stations with LoS/NLoS Propagation

The use of unmanned aerial vehicle base stations (UAV-BSs) as airborne b...
05/21/2018

### Joint Configuration of Transmission Direction and Altitude in UAV-based Two-Way Communication

When considering unidirectional communication for unmanned aerial vehicl...
10/15/2018

### Learn to Fly: A Distributed Mechanism for Joint 3D Placement and Users Association in UAVs-assisted Networks

In this paper, we study the joint 3D placement of unmanned aerial vehicl...
11/18/2019

### Optimal deployments of UAVs with directional antennas for a power-efficient coverage

To provide a reliable wireless uplink for users in a given ground area, ...
##### This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

## I Introduction

Unmanned aerial vehicle (UAV)-enabled networks have been gaining substantial attention because of their wide variety of applications including surveillance, military, and rescue operations [1]. In particular, deployment of UAVs performing as flying aerial base stations (BSs) to support existing cellular networks is a key application [2]

. Typically, UAV-BSs can establish a line-of-sight (LoS) connection to on-ground users with high probability owing to their high altitudes. Thus, the user coverage can be improved significantly in an efficient manner

[3]. In addition, since UAVs have high autonomous mobility, UAV-BSs can provide connections to on-ground users in disaster areas (e.g., flood- and earthquake-affected areas), or rural areas more robustly and cost-effectively than ground BSs of cellular networks [4]. Furthermore, rapid and flexible deployment of UAV-BSs can enable them to respond to the occurrence of hotspots in sports events and open-air concerts and can help achieve on-demand coverage for these.

Despite the potential benefits of UAV-BSs, determining the most effective deployment of UAVs poses several research challenges [2]

. First, owing to the flexible mobility of UAVs, the 3D deployment problem of UAVs exhibits a higher degree of freedom and is more complicated than that of ground BSs. Furthermore, air-to-ground (A2G) channels have different characteristics from a terrestrial one because they are highly dependent on the altitudes of the UAVs and the blockage effect of obstacles (e.g., buildings)

[5, 6]. Moreover, the deployment of multiple UAV-BSs induces an inter-cell interference problem, which may degrade the communication quality of users. Since the interference received from each UAV depends additionally on the A2G channel condition, the overall characteristics of interference and the signal-to-interference-plus-noise-ratio (SINR) of users are also exceedingly complicated. As a result, the optimal 3D deployment of UAVs considering the above factors is a significantly challenging problem.

Owing to the importance of the UAV deployment problem, there has been an increasing number of studies that address the aforementioned challenges [8, 7, 2]. However, there are few studies that consider the spatial and temporal variations of ground users. For example, the works [5, 8] proposed optimal UAV deployment methods that maximize the coverage area. However, in general, the density of users is spatially varied for several reasons, such as the higher densities at stations and sports events. Therefore, to maximize the user coverage, it is insufficient to maximize the coverage area, and the spatial inhomogeneous density of users must be considered.

Furthermore, the density of ground users varies temporally because they may leave and join the network at any time and move dynamically at each moment. Although several studies

[7, 4] proposed optimal UAV deployment methods assuming that the specific positions of all the users are known, these positions may change temporally, and tracking all the user movements is unrealistic. In addition, these methods were fundamentally centralized approaches, which are unsuitable for a more dynamic network. This is because the computational expense is likely to escalate with an increase in the number of users or in the total area. Therefore, it is crucial to develop an efficient 3D deployment method for UAV-BSs that can respond to the spatial and temporal dynamics of users.

In this paper, we propose a novel distributed 3D deployment method for UAV-BSs in a downlink network. In our method, each UAV dynamically updates its 3D position by collaborating with neighboring UAVs so that the overall coverage of users is maximized in an on-demand manner. The distributed nature and incremental updates of our method enable it to be applied to spatially and temporally varying networks.

Our proposed method consists mainly of two key parts that address the above challenges in spatially and temporally varying ground users: i) sensing-aided crowd density estimation; and ii) distributed push-sum algorithm. The first part is used to estimate the density of users based on the information from ground sensors deployed in the network. As it is computationally intensive to obtain all the positions of ground users at each moment in time, we assume that ground users are distributed by a spatial point process. However, its intensity (density) function is also generally unavailable. For this problem, we assume that the ground sensors can detect users nearby them, and estimate the number of users, e.g., by a video surveillance system with human detection/tracking methods [9] and Wi-Fi access points with received signal strength indicator (RSSI)- or channel state information (CSI)-based sensing methods [10, 11]. Since the intensity of users is typically considered to be spatially correlated (e.g., road system), we infer the entire intensity function from periodically gathered sensing results.

In the second step, the UAVs dynamically optimize their positions in a distributed manner. Fundamentally, the challenge posed by distributed 3D UAV deployment is additional to the present problems such as the A2G channel characteristics. This is because the SINR of a user served by a UAV is affected by interference from other UAVs. Thus, to control the SINR of the user, each UAV needs to consider the positions of all the other UAVs and update its position so that the SINR of all users are improved, and not only its serving users is imporved. To address this problem, we apply a distributed push-sum algorithm framework [12, 13]. With this framework, UAVs incrementally optimize their 3D positions by collaborating with their neighbor UAVs. We also prove that the algorithm converges to a consensus among the UAVs. To the best of our knowledge, this is the first study to develop a distributed 3D UAV deployment method with guaranteed convergence. Regarding the performance metric of a UAV deployment, we adopt the total coverage of users, i.e., the total expected number of users whose SINR exceeds a certain threshold. By expressing the coverage probability of a user theoretically and analyzing the total coverage, we apply this performance metric to the distributed push-sum algorithm framework. Furthermore, we evaluate our method with extensive simulations. The results reveal that our method can successfully improve the overall coverage with a limited number of ground sensors. We also apply our method to a dynamic network scenario and demonstrate that it can respond to moving hotspots and provide coverage to users in an on-demand manner.

The remaining part of this paper is organized as follows. Section II summarizes related studies. In Section III, we describe our network model in detail. Section IV presents the proposed method. In Section V, we provide several simulation results. Finally, Section VI concludes our paper.

## Ii Related Work

Since UAV-BSs can be potentially used for demanding future networks, several deployment methods have been proposed. In particular, several studies considered an optimal 3D-deployment of a single UAV. Al-Hourani et al. [5] presented the optimal altitude of a UAV that maximizes the coverage area. They also developed a closed-form expression of the LoS probability of a UAV based on a probabilistic LoS model provided by the International Telecommunication Union (ITU-R) [14]. Bor-Yaliniz et al. [4] proposed an efficient deployment method that maximizes the number of covered users. Alzenad et al. [15] further considered a scenario where users have heterogeneous quality of service (QoS) requirements.

In contrast to the single-UAV cases, multiple-UAV scenarios lead to an inter-cell interference problem. Mozaffari et al. [16] considered the interference among two UAVs and proposed an optimal placement method that maximizes the covered area. Moreover, the authors extended this method to a multiple-UAV scenario in [8]. Kalantari et al. [7]

proposed a meta-heuristic-based optimal 3D placement algorithm that maximizes the user coverage using the minimum number of UAVs. Furthermore, the work

[17] proposed a Q-learning-based 3D deployment method for maximizing the data rate of ground users. However, the above methods are fundamentally centralized and offline approaches. Since the 3D placement optimization problem of UAVs has a high degree of freedom and the computational cost increases with the number of users or deployment area, they are not applicable to dynamic scenarios in which the density of users changes temporally.

Despite the importance of the dynamic scenarios, dynamic and distributed UAV-BSs deployment problems have not been extensively studied. Several recently proposed dynamic deployment methods focused on constructing a desired formation for relay networks [18] or target tracking [19, 20]. Thus, they are not for UAV-BS networks, i.e., did not consider ground users or backbone networks. Zhao et al. [21] proposed a heuristic-based distributed motion-control method for UAV-BSs. However, they focused on the connectivity among UAVs and only UAVs’ horizontal movements. They did not consider the probabilistic characteristics of the SINR of users including A2G channel and interference. In contrast, we theoretically formulate the distributed UAV 3D deployment problem by considering interference among UAVs and the coverage probability of users. We also demonstrate that by controlling the UAV altitudes, the total coverage can be efficiently improved.

## Iii Model Description

### Iii-a Network modeling

In this paper, we consider a downlink network in which UAVs act as flying BSs and aim to enhance the communication quality of users. We first provide an overview of the network model. As shown in Fig. 1, the network consists mainly of four types of elements: UAVs, on-ground user equipments (UEs), on-ground sensors (GSs), and a remote gateway (RG). UAVs provide internet connections to the UEs by connecting to the RG, which also provides UAVs connections to back-bone networks (e.g., satellites and larger UAVs [21]). GSs are assumed to be capable of detecting UEs nearby the GSs and estimating their (average) number (e.g., Wi-Fi access points with RSSI- or CSI-based sensing methods [10, 11]). The sensing results observed at GSs are periodically reported to a central server through a wired or wireless link. Then, the central server estimates the entire intensity of the UEs and distributes it among the UAVs through the RG.

We now describe our model in detail. Table I summarizes the notations used in this paper. We assume that time is slotted as . Each time slot corresponds to the time at which the UAV positions are updated. We assume that the UAVs can move anywhere in a certain closed convex set where and that the UEs are in a closed convex set on the ground. Here, represents the projection from the 3D space to the plane (ground). Let () denote the position of the -th UAV at time . Here, denotes the total number of UAVs. We also write . Furthermore, we assume that the UEs are distributed according to a certain (inhomogeneous) spatial point process with an intensity function . Although depends on due to the dynamics of UEs, we omit the index for simplicity.

Next, we describe the activities of GSs. Let denote the positions of GSs. Here, represents the total number of the GSs. The GSs are assumed to be capable of detecting and estimating the number of UEs around them by sensing and report their sensing results periodically at time . Here, denotes the sensing period and is assumed to be larger than UAVs’ update interval (i.e., ) due to the computational time of sensing or crowd density estimation. Furthermore, let denote the sensing result reported from the -th GS at time , i.e., the (average) expected number of UEs per unit space in the neighborhood of . Thus, it can be regarded as an observation of at time . We also express . Using , the crowd density estimation part infers the entire intensity function as (). Note that we can easily consider the scenario where UAVs are equipped with sensors by taking into account sensing results at ’s. However, in this study, we assume that only GSs can perform sensing, for simplicity.

### Iii-B Channel modeling and UAV-cell definition

We consider the following radio propagation model for UAV–UE channels. Due to the blockage effect in the A2G link, the link conditions can be divided into LoS and non-LoS (NLoS) conditions. We adopt a probabilistic LoS model similar to that in [5, 23], wherein a link condition (, ) is determined by the following probability:

 P(L;θi,y(t)) =11+b1exp(−b0(180πθi,y(t)−b1)), (1) P(N;θi,y(t)) =1−P(L;θi,y(t)), (2)

where and are constants dependent on the environments (e.g., urban and rural) and are provided in [5], and is the elevation angle from a UE at to the -th UAV. Thus, if we denote as the distance between the -th UAV and a UE at , . For simplicity, we assume that all the UAVs and UEs are equipped with omni-directional antennas. Although the channel condition is spatially and temporally correlated in general, for simplicity, we assume that they are independently determined at each time and position. We assume a distance-dependent path-loss model [23, 24] such that, for a distance ,

 ℓq(d)=βq(ε0+d)−αqfor q∈{L,N},

where is a path-loss exponent, denotes a constant to avoid a singularity at , and is a constant that depends on the environment. Furthermore, we model the small-scale fading effect in UAV-UE channels with the independent Nakagami- fading similar to that in [24]. Let denote the small-scale fading gain corresponding to the channel between the -th UAV and a UE at . Then,

is distributed according to the normalized gamma distribution with a parameter

. For simplicity, we assume to be an integer throughout this paper. In addition, we omit the effects of shadowing and Doppler shift for mathematical tractability. The transmission power of each UAV is identical and is assumed to be normalized to one.

We next define the serving area of each UAV, i.e., a UAV-cell. In this paper, we assume that each UAV has a capacity sufficient for serving UEs. Moreover, each UE is assumed to associate with the UAV that provides the strongest average signal power. By definition, the average signal power from the -th UAV to a UE at is expressed as

 Si,y(t)=∑q∈{L,N}P(q;θi,y(t))ℓq(di,y(t)). (3)

The UAV-cell is then defined as a signal-weighted Voronoi cell, i.e., for the -th UAV,

 Vi(\boldmathu(t))={\boldmathy∈Ω2∣Si,y(t)≥Sj,y(t),∀j∈U∖{i}}. (4)

From this expression, we can construct an associated signal-weighted Delaunay graph by choosing the set of vertices as and the set of edges as pairs of UAVs whose signal-weighted Voronoi cells are adjacent. Let denote the set of neighboring UAVs of the -th UAV on the graph . We assume that the -th UAV can communicate with only in our distributed collaborative deployment method.

### Iii-C Problem description

In this paper, we consider the total coverage of UEs, i.e., the total expected number of covered UEs, as the performance metric of UAV deployment. We assume that a UE is covered if the SINR at the UE exceeds a threshold . By considering a random channel condition and fading, the coverage probability of a UE at served by the -th UAV is expressed as

 Ci,y(\boldmathu(t)):=Ci,y(\boldmathu(t);Θ)=P(SINRi,y(t)>Θ). (5)

Here, is expressed as

 SINRi,y(t)=hi,y(t)ℓq(di,y(t))Ii,y(t)+σ, (6)

where is the total interference power defined as

 Ii,y(t)=∑j∈U∖{i}hj,y(t)ℓq(di,y(t)), (7)

and denotes the thermal noise power and is assumed to be constant. By considering the coverage probability, we can characterize the communication quality of UEs more flexibly than the existing performance metrics based on mean path-loss [4, 15] or average SINR [7, 21]. By taking the expectation of the number of covered UEs on and applying Campbell’s theorem [26], the total coverage can be expressed as

 ∫Ω2E[\mathbbold1(UE at % \boldmathy is covered)]λ(\boldmathy)d% \boldmathy (8)

where denotes the indicator function.

We are now prepared to state our optimization problem. Our objective is to determine an optimal 3D deployment of UAVs that maximizes (8) in a distributed manner. Since the true intensity cannot be obtained in general, we use an estimated intensity for the deployment optimization instead.

###### Problem 1

Assume that i) each UAV can communicate only with its neighborhood on and ii) the estimated intensity at () based on sensing result is given by . For each sensing period , determine the optimal UAV deployment that maximizes the estimated total coverage defined as

 F(\boldmathu∣\boldmathx(Tk)) =∑i∈U∫Vi(\boldmathu)Ci,y(\boldmathu)˜λk(\boldmathy)d\boldmathy. (9)

## Iv Distributed UAV 3D Deployment Method

In this section, we present the distributed collaborative UAV 3D deployment method. Our method mainly consists of two parts: i) sensing-aided crowd density estimation; and ii) distributed push-sum algorithm. At each sensing time , GSs send to a central server. The crowd estimation part then estimates . During the -th sensing period, UAVs update their 3D positions by using via the distributed push-sum algorithm, which is aimed at solving Problem 1.

This section is organized as follows. We first analyze our performance metrics, namely the coverage probability and total coverage in Section IV-A. We then describe the crowd density estimation part in Section IV-B and present the distributed push-sum algorithm part in Section IV-C. Finally, we prove the convergence of the distributed algorithm in Section IV-D.

### Iv-a Coverage analysis

Since our distributed optimization method is gradient-based, we analyze the total coverage and derive its gradient in this section. We omit the index in this section because a fixed time is being considered. We start by deriving the theoretical expression of the coverage probability of UEs. By using the expressions of (5) and (6), we can obtain the following result. A sketch of the proof is provided in Appendix A.

###### Lemma 1

If the deployment of UAVs is , the coverage probability can be approximated by

 Ci,y(\boldmathu)≈˜Ci,y(% \boldmathu) ≜∑q0∈{L,N}P(q0;θi,y)mq0∑k=1(−1)k+1 ×(mqk)e−kσγq0,i,yLIi,y(kγq0,i,y), (10)

where with for , and is given by

 LIi,y(s)=∏j∈U∖{i}∑qj∈{N,L}P(qj;θj,y)(1+sℓqj(dj,y)mqj)−mqj.

On the basis of Lemma 1, we aim to solve Problem 1 through in what follows. In other words, we replace in the objective function in (9) by and consider the optimization problem of the following function:

 ˜F(\boldmathu∣\boldmathx(Tk)) =∑i∈U∫Vi(\boldmathu)˜Ci,y(\boldmathu)˜λk(% \boldmathy)d\boldmathy. (11)

Furthermore, from (10), we can calculate the derivatives of with respect to , i.e.,

 ˜\boldmathcji,y(\boldmathu)≜dd\boldmathuj˜Ci,y(% \boldmathu),i,j∈U.

Since we can conveniently obtain from Lemma 1, we omit their explicit expressions due to space limitations.

Next, we demonstrate that the derivative of with respect to can be approximated through . A sketch of the proof is provided in Appendix B. In what follows, we simply write for readability.

###### Lemma 2

If the estimated intensity function is (),

 d˜F(\boldmathu)d\boldmathu% j≈∑i∈U∫Vi(\boldmathu)˜\boldmathcji,y(\boldmathu)˜λ(\boldmathy)d\boldmathy,j∈U. (12)

According to Lemma 2, the gradient of exhibits the following useful (approximate) separation property, which enables us to optimize in a distributed manner:

 ∇˜F(\boldmathu) ≈∑i∈Uˆ\boldmathfi(\boldmathu), (13)

where, for ,

 ˆ\boldmathfi(\boldmathu) ≜[∫Vi(\boldmathu)˜\boldmathcji,y(\boldmathu)˜λ(\boldmathy)d\boldmathy]⊤j∈U =∫Vi(\boldmathu)∇˜Ci,y(\boldmathu)˜λ(\boldmathy)d\boldmathy. (14)

The above relationship indicates that there exists such that and .

### Iv-B Crowd density estimation

We next describe the crowd density estimation. The aim of this part is to estimate the intensity function from the sensing results at sensing time (), i.e., to obtain from . Typically, the density of UEs is spatially correlated owing to geographical factors, such as roads and attractions. Thus, by modeling this spatial correlation, we infer the entire intensity function with a limited number of GSs.

For a fixed sensing time, we consider a spatial Gaussian process (GP) prior for . Here, the constant expresses the mean of and corresponds to a prior expectation of . Meanwhile, is the kernel function of and represents the spatial correlation of . In this paper, we adopt a well-accepted Gaussian kernel:

 k0(\boldmathy,\boldmathy′)=A0exp(−∥\boldmathy−\boldmathy′∥2A1),\boldmathy,\boldmathy′∈R2, (15)

where denotes the Euclidean norm and and are constants. We assume that the parameters of GP (i.e., , , and ) are obtained from statistical data. GPs are widely used in modeling various spatial data [27] ranging from geology to environmental sciences. According to an established property of a GP, if we regard observations as realizations of () at

, they are distributed with a multivariate Gaussian distribution. Furthermore, the posterior distribution of

given also becomes a Gaussian (e.g., [27]): as follows:

 λ(\boldmathy)∣\boldmathx(Tk),{\boldmathgn;n=1,…,N} ∼N(μy(\boldmathx(Tk)),k0(\boldmathy,\boldmathy)−\boldmathk⊤y% \boldmathK−1g\boldmathky), (16)

where

 μy(\boldmathx(Tk)) =μ+\boldmathk⊤y\boldmathK−1g(\boldmathx(Tk)−μ\boldmathe), (17) \boldmathky =[k0(\boldmathg1,\boldmathy),…,k0(\boldmathgN,\boldmathy)]⊤, n,n′ =k0(\boldmathgn,\boldmathg′n),n,n′=1,…,N,

and

denotes a vector of ones with an appropriate dimension. As a result, UAVs can estimate the entire intensity from the observations

by setting

 ˜λk(\boldmathy)=μy(\boldmathx(Tk)),\boldmathy∈Ω2. (18)

By substituting (18) into (11), in (11) is now rewritten as

 ˜F(\boldmathu∣\boldmathx(Tk)) =∑i∈U∫Vi(\boldmathu)˜Ci,y(\boldmathu)μy(\boldmathx(Tk))d\boldmathy. (19)

### Iv-C Distributed push-sum algorithm

In this section, we describe the distributed push-sum algorithm. In this part, UAVs update their positions and maximize (19) in a distributed manner. Our objective function (19) depends on the positions of all the UAVs. Therefore, it is not trivial to solve the problem in a distributed manner. However, the useful separation property of in (13) enables us to apply a distributed push-sum framework. This framework was first proposed by [12] for average computation. It was recently extended to a time varying network by [28], and to a non-convex optimization problem by Tatarenko and Touri [13].

To state the algorithm, we first introduce an extension of that is defined on and satisfies

 (20)

where denotes the complement of an open set disjoint with . In addition, is a constant and denotes a unit vector directed from toward its closest point on the boundary of . In the next section, we prove the existence of that satisfies several conditions required for the convergence of the distributed push-sum algorithm.

We now describe the algorithm in detail. The pseudo-code is given in Algorithm 1. Each UAV (i.e., the -th UAV) maintains a -dimensional vector , which represents the pseudo positions of the UAVs. Here, pseudo indicates that this location information is local for the UAV and does not always represent the exact positions of the other UAVs (only expresses the exact position of the UAV ). This is because each UAV can only communicate with its neighboring UAVs. They also maintain auxiliary -dimensional vectors , and , and a scholar variable . In each time step, each UAV updates its local variables according to the following rules:

 \boldmathw[i](t+1) :=∑j∈Ni(t)\boldmathξ[j](t)dj(t), ϕ[i](t+1):=∑j∈Ni(t)ϕ[j](t)dj(t), \boldmathu[i](t+1) :=\boldmathw[i](t+1)ϕ[i](t+1), (21) \boldmathξ[i](t+1) :=\boldmathw[i](t+1) +a(t+1)(ˆ\boldmathf∗i(% \boldmathu[i](t+1))+\boldmathκ[i](t+1)),

where

are i.i.d. random vectors whose entries are independent random variables with zero mean and unit variance for all

. Moreover, is a positive non-increasing step-size sequence such that and is the node degree of the UAV on .

The update rules in (21) can be interpreted simply as follows. Each UAV aims to reach a consensus among the UAVs by exchanging its local information with its neighbors and calculating their weighted sum. Furthermore, by combining with a gradient descent based on its local information (i.e., ), the consensus point is steered towards a critical point of the objective function. In addition, the random perturbation ensures the convergence to a local optimum [13].

### Iv-D Proof of convergence

In this section, we prove that the distributed push-sum algorithm given in (21) converges to a local optimum of the objective function among all the UAVs by employing the framework in [13]. To do this, we first prove the existence of an extension in (20) and the antiderivative of satisfying several regularity properties.

###### Lemma 3

There exist functions and on such that and in and

1. is bounded on for all ;

2. is Lipchitz continuous on for all ;

3. ; and

4. .

###### Proof.

Due to space limitations, we provide only a summary of the proof. In what follows, we fix arbitrary. We can confirm that are all bounded on . We can also confirm that is bounded on (see (17)) and is a closed set. It then follows from (14) and (18) that is bounded on . Furthermore, we can prove that is continuously differentiable in and the derivatives are bounded on , which indicates that is Lipschitz continuous with a certain Lipschitz constant . Thus, if we extend the domain of to similarly as in (20), is bounded and Lipschitz continuous with on . It thus follows from Tietze’s extension theorem that there exists a bounded and Lipschitz continuous function on with such that in . These results indicate that the statements (i) and (ii) hold. Furthermore, by using and (13), we can construct a continuously differentiable function on such that in and the statements (iii) and (iv) hold, which completes the proof. ∎

Combining Lemma 3 with the facts that is differentiable in and the sequence is strongly connected, i.e., the union of their edge sets is strongly connected, we can apply Theorem 5 in [13] to and in Lemma 3 and prove the convergence of the algorithm.

###### Theorem 1

Each and the average state variable in the distributed push-sum algorithm in (21) converge to a point in the set of local maxima of from any initial state.

## V Simulation Results

We evaluate our distributed 3D UAV deployment method by conducting numerical simulations. In each simulation, we assume that UAVs are deployed over a ( m) area () with a maximum altitude 1500 m. We adopted values similar to those in [25, 24] for the parameters of the channel model associated with the and conditions. We also considered sub-urban, urban, and dense-urban scenarios by selecting the values of and presented in [4]. The other parameters are listed in Table II unless otherwise stated. In addition, we randomly generated spatially correlated intensity by the following procedure. We first divided the whole area into small grids, and then sampled the values of in each grid from a multivariate Gaussian distribution whose covariance matrix is determined by the kernel function in (15).

### V-a Effectiveness of distributed push sum algorithm

We first focus on the performance of the distributed push-sum algorithm part and demonstrate its effectiveness. We also investigate the impacts of several parameters on the performance and the resulting UAV deployment. To achieve this, we focused on a sensing period, i.e., a fixed , and assumed that all the UAVs knew the true , i.e., (). We then set artificial hotspots, in which the intensity of UEs is (see Fig. 3 (a)). In addition, the initial positions of the UAVs were uniformly arranged with the equal intervals.

Fig. 3 shows the 3D placement of the UAVs at in the urban scenario. We also plotted in Fig. 3 (c) the case where the altitudes of the UAVs were fixed at the default altitude by controlling only their horizontal movements. The triangles in the graphs express the positions of the UAVs, and their altitudes are represented by their color depth. Each cell corresponds to each UAV-cell. The background color depth expresses the intensity of UEs, and the dark ellipses correspond to the artificially added hotspots. We also illustrate the process of the total coverage (i.e., the objective function) in Fig. 3. Fig. 3 (a) and (b) reveal that the 3D positions of the UAVs and associated shapes of the UAV-cells dynamically adjusted to the hotspots so that each hotspot was fitted in a UAV-cell to avoid interference. As a result, the total coverage was significantly improved from the initial deployment (approximately 35 % higher). On the other hand, the total coverage was less improved in the fixed altitude scenario. Fig. 3 (b) also shows that the altitudes of the UAVs serving hotspots were higher than those of the other UAVs. Moreover, the altitudes of the UAVs changed more aggressively than their horizontal positions. This is because as the altitude of a UAV increases, the LoS probability increases, and thus the UAV-cell and the coverage probability also increase. These results indicate that the total coverage is more sensitive to the vertical movements of the UAVs than to their horizontal ones. Thus, their 3D position must be considered for an efficient UAV deployment optimization.

Fig. 5 shows the process of the total coverage in sub-urban, urban, and dense urban scenarios. Although the total coverage increased in all cases, their improvements varied depending on the scenarios. The reason for this can be considered as follows. Since the coverage area of a UAV is larger in the sub-urban scenario than the other scenarios [5, 4], the whole area was sufficiently covered even with an initial placement. Thus, the benefit of our method was smaller in the sub-urban scenario than the urban one. On the other hand, in the dense-urban scenario, ground UEs could not be covered by UAV-BSs due to the small coverage area.

We next investigate the impacts of the SINR threshold . To see the difference more clearly, we calculated the ratio of the total coverage at each step to that at . Fig. 5 expresses the processes of the ratio with different . The graph shows that the improvement with higher becomes larger, whereas that with dB was highly limited. This is because with low , the total coverage was already sufficient even in the initial step, and thus the objective function did not exhibit much potential for improvement. Meanwhile, when dB, the total coverage became approximately two times of that in the initial step. This result indicates that our optimization method is highly efficient for a high SINR regime. We also illustrate the 3D deployment of the UAVs in Fig. 6 in each case. The deployment with dB was similar to that with dB, whereas the altitudes of the UAVs with dB were much lower than the other cases.

Finally, we demonstrate the impact of the initial positions of the UAVs. Fig. 7 shows the 3D deployment of the UAVs at when their initial positions were gathered at the center. The figures show that as the time elapsed, the positions of the UAVs adjusted to the hotspots and the UAVs maintained certain distances among themselves to reduce interference. In addition, the altitudes of the UAVs serving the hotspots became high again, whereas the resulting deployment was different from that in the previous scenario.

### V-B Effect of crowd density estimation

We next focus on the effect of the crowd density estimation part. Unlike the previous scenario, we assumed that the intensity function is available only at randomly generated GSs (grids where GSs exist). The entire intensity was calculated by the crowd density estimation method described in Section IV-B. Since the number of GSs is likely to impact the performance, we evaluated it by considering the GS deployment ratio . This is defined as the ratio of the grids where GSs exist to the total grids. We regarded each grid as the sensing range of a GS111In our example, the area of each grid is 100 m (2500 grids). and assumed that the intensity in each GS’s grid was obtained correctly. To evaluate the impact of , we consider the coverage gain, which was calculated as the difference between the total coverage at and . We first calculated the coverage gain in the case with complete information of true intensity, (i.e., ) as a baseline. We then calculated the ratio of the coverage gain with each to the baseline. Fig. 8 illustrates the result of the ratio of coverage gain. For each , we generated 10 samples of the GS placement. The figure shows that as expected, as increases, the performance of of our method also increases. This is because the proposed algorithm optimizes the positions of UAVs by using (see (11)). Thus, if the intensity function was incorrectly estimated, the proposed method failed to update the UAV placement effectively for improving the objective function. However, the simulation results demonstrate that even a small (e.g., , 3 GS/km) achieved 98% of the baseline in the urban scenario. This indicates that our method could improve the total coverage even with a limited number of GSs.

### V-C Dynamic network scenario

We next apply our algorithm to a more dynamic network, in which the positions of hotspots change dynamically. Unlike the previous examples, we consider multiple sensing periods and assume that hotspots move at each sensing time. We first artificially set hotspots on similarly as those in Fig. 3. We then moved each of them in a predefined direction at with . The initial deployment of UAVs was also the same as Fig. 3 (a). Fig. 9 illustrates the 3D deployment of UAVs at the end of each sensing period (e.g., (b) shows , the end of the first sensing period ()). We also show the process of the total coverage in Fig. 10. As the figures show, at each sensing time, the total coverage temporarily decreased due to the movement of the hotspots. However, as shown in Fig. 9, the UAVs adjusted their positions to the hotspot locations and successfully increased the total coverage. For example, in Fig. 9 (b) and (d), (i.e., ), the left hotspot moved to a position on the borders of UAV-cells. However, the UAVs changed their 3D positions and the shapes of the UAV-cells to reduce interference. Note that, at (b), the right hotspot is on the border of UAV cells. This is because our algorithm could not sufficiently optimize the UAV placement during a sensing period. However, the total coverage was successfully improved as shown in Fig. 10. Consequently, our method could respond to hotspots and provide coverage to UEs in a dynamic network in an on-demand manner.

## Vi Conclusion

In this paper, we proposed a distributed UAV-BS 3D deployment method for on-demand coverage. Since the specific positions of ground users may not be obtained in real networks, we proposed the sensing-aided crowd density estimation method. Furthermore, by adopting total coverage as the performance metric, we developed a push-sum-based distributed UAV 3D deployment algorithm. Several simulation results revealed that our method can improve the overall coverage of ground users with the aid of a limited number of sensors and can be applied to a dynamic network.

In our problem setting, we assumed that the UAVs exhibit identical transmission power. Thus, the joint optimization of 3D UAV deployment and energy consumption are topics for our future work. Furthermore, other physical constraints, such as the capacity of UAVs [7], and obstacle avoidance [21] should be considered while applying our method to a real network. Moreover, simulations considering more realistic user density, such as road constraints, remain for future work.

## Acknowledgment

This work was supported by JSPS KAKENHI Grant Numbers JP19K14981 and JP18K13777.

## Appendix A: Proof of Lemma 1

By conditioning on the status of the desired channel between and , can be rewritten as

 Ci,y(\boldmathu)=∑q0∈{L,N}P(q0;θi,y)P(SINRi,y>Θ∣q0). (22)

Since small-scale fading gain is distributed according to the normalized gamma distribution with the parameter , by using (6), we obtain, for

 P(SINRi,y>Θ∣q0)=P(hi,y>Θ(Ii,y+σ)ℓq0(di,y)∣∣∣q0) (a)≈1−EIi,y[(1−exp(−ηq0Θ(Ii,y+σ)ℓq0(di,y)))mq0] =mq0∑k=1(−1)k+1(mq0k)e−kσγq0,i,yEIi,y[e−γq0,i,yIi,y] (b)=mq0∑k=1(−1)k+1(mq0k)e−kσγq0,i,yLIi,y(kγq0,i,y), (23)

where we approximate the tail probability of the normalized gamma distribution in (a) by using the same method as in [25, 24], and we use the Laplace transform of in (b). Since the small-scale fading gain and channel condition are i.i.d. for each UAV–UE channel, can be calculated by

 LIi,y(s) =E⎡⎣exp⎛⎝−s∑j∈U∖{i}hj,yℓqj(dj,y)⎞⎠⎤⎦ =∏j∈U∖{i}E[exp(−shj,yℓq(∥\boldmathuj−\boldmathy∥))] =∏j∈U∖{i}∑qj∈{L,N}P(qj;θj,y)E[e−shj,yℓqj(dj,y)∣qj], =∏j∈U∖{i}∑q∈{LoS,NLoS}P(q;θj,y)(1+sℓq(dj,y)mq)−mq,

where the last equality follows from the Laplace transform of a gamma distribution with a parameter . Finally, combining the above with (23) and substituting it into (22), we obtain (10).

## Appendix B: Proof of Lemma 2

In this appendix, we provide an outline of the proof due to space limitations. First, for each UAV , the boundary of is determined by only UAV and its neighbors according to the definition of . Thus, by differentiating with respect to , we obtain

 d˜F(\boldmathu)d%\boldmath$u$j =∑i∈Nj∪{j}dd\boldmathuj∫Vi(\boldmathu