Yiying Zhang

Yue Huang

Chao Huang

Hailong Huang

Anh-Tu Nguyen

Joint Optimization of Deployment and Flight Planning of Multi-UAVs for Long-distance Data Collection from Large-scale IoT Devices

Keywords: Data collection, Internet of Things, multi-UAVs, flight planning, differential evolution

Internet of Things (IoT) devices have been widely deployed to build smart cities. How to efficiently collect data from large-scale IoT devices is a valuable and challenging research topic. Benefiting from agility, flexibility, and deployability, an unmanned aerial vehicle (UAV) has great potential to be an aerial base station. However, given the limited battery capacity, the flight time of a UAV is limited. This paper focuses on using multi-UAVs to execute long-distance data collection from largescale IoT devices. We design a multi-UAVs-assisted large-scale IoT data collection system. The core facilities of this system are the data center and charging stations, which are equipped with a limited number of charging piles to provide charging services for UAVs. To ensure the efficient operation of the system, the problem of deployment and flight planning of UAVs is formulated as a joint optimization problem. To solve the problem, a population-based optimization algorithm with a three-layer structure, namely EDDE-DPDE, is proposed. It includes two core components: elitedriven differential evolution (EDDE) and differential evolution with a dynamic population (DPDE), which are two variants of differential evolution. Thanks to ideas of reusing elite individuals and historical information, the proposed EDDE-DPDE shows an improvement of at least 11.11% compared with four powerful algorithms in terms of average travel time.

I. INTRODUCTION

I N recent years, the concept of the smart city is put forward to provide a comfortable life to citizens by fully utilizing modern technology, such as super-computing systems, Internet of Things (IoT), big data analysis, and unmanned aerial vehicles (UAVs). The major function of IoT in building smart cities is to collect and send data with the assistance of internetconnected devices [START_REF] Alsamhi | Predictive estimation of optimal signal strength from drones over IoT frameworks in smart cities[END_REF]. For some stand-alone IoT devices, their energy supply is constrained. Thus, the efficient usage of limited energy is of great concern. Y. Huang is with Data61 Commonwealth Scientific and Industrial Research Organisation, Melbourne, Australia (Email: yuehuangcsiro@gmail.com).

H. Huang is with the Department of Aeronautical and Aviation Engineering and The Research Institute for Research Institute for Sports Science and Technology (RISports), the Hong Kong Polytechnic University, Hong Kong. (E-mail: hailong.huang@polyu.edu.hk).

A. Nguyen is with the LAMIH laboratory, UMR CNRS 8201, Universit Polytechnique Hauts-de-France, 59300 Valenciennes, France, and also with the INSA Hauts-de-France, 59300 Valenciennes, France (e-mail: nguyen.trananhtu@gmail.com).

UAVs can play many roles in building smart cities, such as target monitoring, logistics distribution, topographic mapping, and power line inspection [START_REF] Raja | Efficient and secured swarm pattern multi-UAV communication[END_REF]- [START_REF] Wang | Taking drones to the next level: Cooperative distributed unmannedaerial-vehicular networks for small and mini drones[END_REF]. In addition, a UAV can also be regarded as an aerial base station to collect data from IoT devices [START_REF] Yang | Multi-UAV-enabled load-balance mobile-edge computing for IoT networks[END_REF]- [START_REF] Li | Rechargeable multi-UAV aided seamless coverage for QoS-guaranteed IoT networks[END_REF]. A lot of effort has been conducted for improving the efficiency of using UAVs to collect data from IoT devices [START_REF] Wang | Energyefficient data collection and device positioning in UAV-assisted IoT[END_REF]- [START_REF] Dong | Joint optimization of deployment and trajectory in UAV and IRS-assisted IoT data collection system[END_REF]. However, they do not cover how to use multi-UAVs simultaneously to carry out long-distance data collection tasks for large-scale IoT devices. There are two challenges associated with this topic:

• With the continuous advancement of smart city construction, it is a trend that large-scale IoT devices are deployed in all aspects of the city. It is a promising approach to deploy multi-UAVs to improve collection efficiency. Specifically, multiple stop points are required for a UAV to collect data from large-scale IoT devices. How to determine the optimal number and locations of stop points for a UAV is a challenge. • In a real collection scenario, a UAV needs to fly from a starting point to a destination point. Then, it executes the collection task at the destination point, after which, it returns to the starting point. Then, how to support longdistance flights and design efficient flight planning (FP) for multi-UAVs is another challenge.

Motivated by these challenges, this paper aims at solving the joint optimization problem of deployment and FP for a multi-UAVs-assisted large-scale IoT data collection system. This problem is to minimize the average travel time of all UAVs to complete data collection. To solve the problem, a new population-based optimization algorithm named EDDE-DPDE is proposed. The main contributions are as follows: TABLE I: The summary of some key references.

Reference

Deployment optimization Flight planning

[9] × ✓ [10] ✓ × [11] × ✓ [12] × ✓ [13] ✓ × [14] ✓ × [15] × ✓ [16] × ✓ [17] × ✓ [18] × ✓ [19] × ✓ [20]
✓ ✓ [START_REF] Dong | Joint optimization of deployment and trajectory in UAV and IRS-assisted IoT data collection system[END_REF] ✓ ✓ a three-layer structure. EDDE is first used to search the optimal FP that all UAVs fly from the DC to their destination charging stations (DCSs). Then, DPDE is employed to determine the optimal deployment of all UAVs during data collection. Lastly, EDDE is applied again to search the optimal FP that all UAVs fly from their DCSs to the DC. • This paper designs a new variant of differential evolution (DE), namely EDDE, to determine the optimal FP for all UAVs. The FP problem is converted into a largescale optimization problem with the help of the priority sequence (PS) generated by a priority-based encoding mechanism. The basic idea of EDDE is to enhance the global search ability of DE by making full use of the obtained best FP, which is reflected in three aspects: 1) an elite archive used to save some promising individuals is created, and a new mutation operator is designed by reusing these promising individuals, 2) a new crossover operator is introduced based on the obtained best PF, and 3) a local escape operator with random numbers produced by opposition-based learning is presented, which is a further optimization for the obtained best FP. • This paper presents another variant of DE, namely DPDE, to determine the optimal deployment of all UAVs during data collection. The basic idea of DPDE is to adaptively search the optimal number and locations of stop points of UAVs, which is reflected in two aspects: 1) three mutation strategies are designed by reusing historical population information, and 2) a random crossover operator is designed based on random numbers with a uniform distribution between 0 and 1. • This paper simulates a realistic scenario of using multi-UAVs to execute long-distance data collection tasks with the multi-UAVs-assisted large-scale IoT data collection system. Specifically, the nearest center rule is made to assign the data collection task and determine the DCS for each UAV. The impact of EDDE and DPDE on the performance of EDDE-DPDE is investigated by comparing with four powerful algorithms.

The rest of this paper includes five sections. Related work is introduced in Section II. Section III and Section IV present the system model and problem formulation, respectively. Section V describes the proposed algorithms. Experimental results are discussed in Section VI. Section VII concludes this paper.

II. RELATED WORK

In recent years, using the UAV to collect data from IoT devices has attracted the interest of many scholars as shown in Table I. Wang et al. [9] propose a novel UAV-assisted IoT network to provide energy-efficient data collection service. Huang et al. [START_REF] Huang | Differential evolution with a variable population size for deployment optimization in a UAVassisted IoT data collection system[END_REF] design a differential evolution with a variable population size to optimize the deployment of a UAV in a UAV-assisted IoT data collection system. Samir et al. [START_REF] Samir | UAV trajectory planning for data collection from time-constrained IoT devices[END_REF] use a UAV to collect data from time-constrained IoT devices. Hu et al. [START_REF] Hu | Aoi-minimal trajectory planning and data collection in UAV-assisted wireless powered IoT networks[END_REF] investigate the UAV-assisted wireless-powered IoT system. Li et al. [START_REF] Li | Data collection maximization in IoT-sensor networks via an energy-constrained UAV[END_REF] employ an energy-constrained UAV to collect data of IoT devices in a sparse IoT-sensor network. Ghdiri et al. [START_REF] Ghdiri | Offline and online UAV-enabled data collection in time-constrained IoT networks[END_REF] adopt UAVs to collect data from timeconstrained sensor nodes. Liu et al. [START_REF] Liu | DRL-UTPS: DRL-based trajectory planning for unmanned aerial vehicles for data collection in dynamic IoT network[END_REF] aim at considering a scenario of using UAVs to collect data from a dynamic IoT network. Wang et al. [START_REF] Wang | Learningbased UAV path planning for data collection with integrated collision avoidance[END_REF] pay attention to employing UAVs to collect data from distributed IoT nodes. In [START_REF] Na | UAV-assisted wireless powered Internet of Things: Joint trajectory optimization and resource allocation[END_REF], a UAVassisted multicarrier wireless-powered communication model is proposed for IoT scenarios. In [START_REF] Tran | UAV relay-assisted emergency communications in IoT networks: Resource allocation and trajectory optimization[END_REF], resource allocation and trajectory design for UAV-assisted full-duplex IoT networks with the emergency communication system are investigated. In [START_REF] Na | UAV-supported clustered noma for 6G-enabled Internet of Things: Trajectory planning and resource allocation[END_REF], a UAV-supported clustered non-orthogonal multiple access system is designed for 6G-enabled IoT. The mentioned references focus on improving the efficiency of data collection by optimizing the deployment or FP of UAVs and do not address the FP and deployment of UAVs at the same time.

In addition, although the FP and deployment problems of UAVs are addressed simultaneously in [START_REF] Han | Joint deployment optimization and flight trajectory planning for UAV assisted IoT data collection: A bilevel optimization approach[END_REF], [START_REF] Dong | Joint optimization of deployment and trajectory in UAV and IRS-assisted IoT data collection system[END_REF], they do not cover large-scale IoT devices and long-distance flights of UAVs. Specifically, the considered deployment and FP problems have their characteristics:

• The deployment problem. Considering the large workload, multi-UAVs are employed. In addition, a UAV needs to collect data at multiple stop points and add its energy multiple times due to the heavy workload. Thus, the considered deployment problem is to determine the optimal number and locations of stop points for multi-UAVs, which takes the least average travel time with the guarantee of adding energy for multi-UAVs.

• The FP problem. The considered FP problem is to determine the optimal FP for multi-UAVs that can take the least average travel time on the round trip paths between the DC and their DCSs. Note that, the distance between the DC and their DCSs is usually out of their maximum range. Thus, how to simultaneously support multi-UAVs to carry out long-distance flights is the main challenge of the considered FP problem. The deployment and FP problems of UAVs can be regarded as optimization problems. The superiority of population-based metaheuristic algorithms on such problems has been proven [START_REF] Wang | Energyefficient data collection and device positioning in UAV-assisted IoT[END_REF], [START_REF] Huang | Differential evolution with a variable population size for deployment optimization in a UAVassisted IoT data collection system[END_REF], [START_REF] Dong | Joint optimization of deployment and trajectory in UAV and IRS-assisted IoT data collection system[END_REF]- [START_REF] Pan | Joint power and 3D trajectory optimization for UAV-enabled wireless powered communication networks with obstacles[END_REF]. The basic idea of the population-based metaheuristic algorithm is as follows. Each individual in the population represents a solution to the solved problem. All individuals share the obtained valuable information, which generate the next generations by the defined iteration rules inspired by natural phenomenon. The iteration rules have two remarkable characteristics: 1) they consider both local Fig. 1. The operation of the built multi-UAVs-assisted large-scale IoT data collection system ("connected line" is that the connected two facilities are reachable; "charging line" is that a UAV returns to its DCS to recharge; "data collection line" is that a UAV returns to the operating point for data collection).

search and global search, and 2) they refer to random numbers that follow uniform distribution or normal distribution. Thus, population-based metaheuristic algorithms have more chances to find better solutions than numerical algorithms in solving complex engineering problems [START_REF] Shaikh | A review on swarm intelligence and evolutionary algorithms for solving the traffic signal control problem[END_REF]. At present, two main methods support the long-distance flight of a UAV, which are the UAV-vehicle model and the UAV-station model [START_REF] Huang | Drone stations-aided beyondbattery-lifetime flight planning for parcel delivery[END_REF]. UAV-vehicle model is the cooperation between the UAV and the vehicle. The vehicle can add the energy of a UAV by replacing the battery or providing charging service [START_REF] Pan | Extending delivery range and decelerating battery aging of logistics UAVs using public buses[END_REF]. UAV-station model can add the energy of a UAV by the deployed charging stations on the ground [START_REF] Huang | Drone stations-aided beyondbattery-lifetime flight planning for parcel delivery[END_REF]. Given the uncertainty of the traffic, the UAV-station model outperforms the UAV-vehicle model in terms of reliability.

Motivated by the mentioned research, this paper builds a multi-UAVs-assisted large-scale IoT data collection system based on the UAV-station model and designs a new populationbased metaheuristic algorithm based on the characteristics of the considered deployment and FP problems.

III. SYSTEM MODEL

On the premise of being close to a real scenario, a multi-UAVs-assisted large-scale IoT data collection system is built as shown in Fig. 1. It consists of a DC, some CSs, some UAVs, and large-scale IoT devices. The DC plays two roles: 1) processing the collected data, and 2) providing some services to a UAV, such as charging and maintenance. The DC and CSs are equipped with a limited number of charging piles to provide charging services to UAVs, which are called nodes that form a CSN to support long-distance flights of multi-UAVs. A UAV is a platform for collecting data from IoT devices, which has three working stages: 1) it flies from the DC to its operating area with the help of CSN as shown in Fig. 1(a), 2) it executes collection task at its operating area as shown in Fig. 1(b), and 3) it returns to the DC with the help of CSN after completing its task as shown in Fig. 1(c).

We take UAV 1 in Fig. 1 as an example to describe the operation of this system. In the first stage, UAV 1 flies from the DC to DCS 1 with the help of CSN. DCS 1 is the closest CS to the operating area of UAV 1 among all CSs. In the second stage, UAV 1 executes the data collection task. During data collection, UAV 1 can return to DCS 1 for charging. In the third stage, UAV 1 flies from DCS 1 to the DC. Note that, given the limited CSs, UAV 1 needs to compete with other UAVs for charging resources, such as CS 5 in Fig. 1(a) and CS 1 and CS 5 in Fig. 1(c). Clearly, as the number of UAVs increases, UAV 1 will bear greater competitive pressure, which also applies to other UAVs.

IV. PROBLEM STATEMENT

This section first introduces preliminary knowledge. Then, the considered deployment problem and FP problem are stated. In addition, to convert the considered practical problems into mathematical models, the following assumptions are made:

• All UAVs are identical, are fully charged at the DC and do not collide with obstacles during the task. • A UAV has three flight states: take-off, landing, and normal flight. Specifically, a UAV flies between two nodes in a straight line with a constant speed and altitude, which takes off and lands vertically above the take-off point and landing point, respectively. • The factors on the energy consumption of a UAV involve take-off, landing, charging, normal flight, and data collection.

• The factors on the travel time of a UAV involve take-off, landing, normal flight, charging, wait, and data collection.

• After a UAV establishes a communication link with an

IoT device, the communication link will not be disturbed by other factors. Remark 1. This work focuses on the high-level flight planning problem of multiple UAVs rather than the low-level motion control problem. So, this work does not present the UAV kinematics and dynamics models. In the considered context, to conduct a mission, a UAV normally takes off, does a level flight to the destination and then lands. During these processes, classic algorithms such as PID can easily achieve the expected motion control.

A. Preliminary knowledge

This section introduces the charging function, distance relationships between nodes and status chart of charging piles.

1) Charging function: A node in the CSN can provide charging service to UAVs. A UAV is typically equipped with series-connected lithium-ion polymer battery cells [START_REF] Kim | Optrone: Maximizing performance and energy resources of drone batteries[END_REF], whose charging process follows a "slow-fast-slow" curve [START_REF] Valle | An area and powerefficient analog Li-ion battery charger circuit[END_REF]. Function (1) describes the charging process:

C = C 0 1 + exp(6 -t/5) , (1)
where C 0 is the maximum energy capacity of a battery, C is the residual energy of a battery, and t is charging time.

Algorithm 1 Calculate T i,j and t w,k,i

Input:

T i,j , t v,k,i , t w,k,i , t w,k,i,j Output: t w,k,i , T i,j 1: Initialize flag bit ft by ft = 0 and t w,k,i,j by t w,k,i,j = 0; 2: for j = 1 to nc do 3:

if t v,k,i > t e,i,j then 4: ft = 1, t w,k,i = 0, t s,i,j =t v,k,i , t e,i,j = t v,k,i + t c,k,i ; 5:

Break; 6: else 7: t w,k,i,j = t e,i,j -t v,k,i ; 8:

end if 9: end for 10: if ft == 0 then 11:

Search the charging pile j * corresponding to the minimum wait time;

12:

t w,k,i = t w,k,i,m , t s,i,j * = t v,k,i , t e,i,j * = t v,k,i + t c,k,i + t w,k,i . 13: end if
2) Distance relationships between nodes: When the distance between two nodes in the CSN is beyond the maximum range of a UAV, the UAV cannot fly between the two nodes. Let n d , n s , and S denote the number of DC, the number of CS, and the sequence of numbers of all nodes, respectively. S can be expressed by:

S = [1, 2, . . . , n d , n d + 1, n d + 2, . . . , n d + n s ]. (2)
Let n n and D max denote the total number of nodes in the CSN and the maximum range of a UAV, respectively. If the distance between two nodes is less than D max , the two nodes are connected; otherwise, the two nodes are not connected. From (2), the distance relationship among nodes can be expressed by a 0-1 matrix with n rows and n columns, which is called the connected matrix M CSN . Specifically, in M CSN , "0" denotes that it is not reachable; "1" denotes that it is reachable.

3) Status chart of charging piles: Let n c denote the number of charging piles at one node. When the number of UAVs (which need to be charged) is larger than n c , some UAVs need to wait until idle charging piles appear. This paper uses a status chart to record the status of charging piles, which can be described by:

T i,j = [t s,i,j t e,i,j ], i ∈ [1, n n ], j ∈ [1, n c ], (3)
where T i,j , t s,i,j , and t e,i,j are the occupied period, start time, and end time of charging pile j at node i, respectively. Let t v,k,i , t w,k,i , and t w,k,i,j denote the time of UAV k arriving at node i, the wait time of UAV k at node i, and the wait time of UAV k at the jth charging pile of node i, respectively. T i,j and t w,k,i can be computed by Algorithm 1. Here, k ∈ [1, n u ] and n u is the number of UAVs.

B. The FP problem

The considered FP of a UAV includes: 1) the UAV flies from the DC to its DCS as shown in Fig. 1

L a,k = [L a,k,1 , L a,k,2 , . . . , L a,k,l a,k ], l a,k ≥ 3, (4)
L c,k = [L c,k,1 , L c,k,2 , . . . , L c,k,l c,k ], l c,k ≥ 3, (5)
where l a,k and l c,k are the lengths of L a,k and L c,k , respectively.

The travel time consumed on L a,k consists of four parts:

• Take-off and landing time t tol,k . From (4), UAV k takes off l a,k -1 times and lands l a,k -1 times. Thus, t tol,k can be obtained by:

t tol,k = t tl (l a,k -1), (6)
where t tl is the consumed total time by UAV k on taking off once and landing once. • Normal flight time t nf,k . Let v denote the constant speed of UAV k on the normal flight. t nf,k is computed by:

t nf,k = L a,k-1 s=1 d L a,k,s ,L a,k,s+1 v , (7)
where d L a,k,s ,L a,k,s+1 is the Euclidean distance between node L a,k,s and node L a,k,s+1 .

• Charging time t ct,k . Let t k,L a,k,h denote the charging time of UAV k at node L a,k,h
, where h is an integer between 1 and l a,k . Let C r,k,L a,k,h and C n,k,L a,k,h are the rest energy and required energy of UAV k at node L a,k,h , respectively. From (1), to get t k,L a,k,h , C r,k,L a,k,h and C n,k,L a,k,h need to be computed. To improve the flight efficiency, C r,k,L a,k,h is the minimum energy that can allow UAV k to fly from node L a,k,h to node L a,k,h+1 . In addition, the penalty energy C p is introduced to prevent a UAV from running out of energy. Thus, C n,k,L a,k,h can be achieved by:

C p = ηC 0 , (8)
C n,k,L a,k,h = P 0 d L a,k,h ,L a,k,h+1 v + C tl + C p , (9)
C n,k,L a,k,h = C 0 , if C n,k,L a,k,h > C 0 , ( 10
)
where η is the penalty factor that is set to 0.1 in this paper, C tl is the consumed total energy on take-off once and landing once, and P 0 is the rated power of UAV k. Note that, ( 10) is to make C n,k,L a,k,h obtained from ( 9) not more than the maximum energy capacity of a UAV. To compute C r,k,L a,k,h , the real penalty energy C np,k,L a,k,h is defined by:

C np,k,L a,k,h = C 0 -C n,k,L a,k,h , if C r,k,L a,k,h > C 0 ; C p , if C n,k,L a,k,h ≤ C 0 . (11
) According to [START_REF] Samir | UAV trajectory planning for data collection from time-constrained IoT devices[END_REF], C r,k,L a,k,h can be computed by: [START_REF] Li | Rechargeable multi-UAV aided seamless coverage for QoS-guaranteed IoT networks[END_REF][START_REF] Wang | Energyefficient data collection and device positioning in UAV-assisted IoT[END_REF][START_REF] Huang | Differential evolution with a variable population size for deployment optimization in a UAVassisted IoT data collection system[END_REF] and C r,k,L a,k,h obtained from [START_REF] Samir | UAV trajectory planning for data collection from time-constrained IoT devices[END_REF][START_REF] Hu | Aoi-minimal trajectory planning and data collection in UAV-assisted wireless powered IoT networks[END_REF] into (1)), which can be described by:

C r,k,L a,k,h = C 0 - P0d L a,k,h ,L a,k,h+1 v -C tl , if h = 2; C np,k,L a,k,h , if h ∈ (2, l a,k ). ( 12
) Let t ct,k,L a,k,h denote the charging time of UAV k at node L a,k,h . t ct,k,L a,k,h can be computed by substituting C n,k,L a,k,h obtained from
t ct,k,L a,k,h = 5 ln( C 0 C r,k,L a,k,h -1) -5 ln( C 0 C n,k,L a,k,h -1).
(13) Thus, t tc,k can be expressed by:

t ct,k = L a,k-1 h=2 t ct,k,L a,k,h . (14)
• Wait time t wt,k . From Algorithm 1, t wt,k is obtained by:

t wt,k = L a,k-1 h=2 t w,k,L a,k,h . (15)
Let t fdc,k and t fcd,k denote the travel time on L a,k and L c,k , respectively. Based on ( 6), ( 7), [START_REF] Ghdiri | Offline and online UAV-enabled data collection in time-constrained IoT networks[END_REF], and (15), t fdc,k can be computed by:

t fdc,k = t tl,k + t nf,k + t ct,k + t wt,k . (16)
Similarly, t fcd,k also can be obtained according to [START_REF] Alsamhi | Predictive estimation of optimal signal strength from drones over IoT frameworks in smart cities[END_REF][START_REF] Raja | Efficient and secured swarm pattern multi-UAV communication[END_REF][START_REF] Tang | An improved UAV-PHD filter-based trajectory tracking algorithm for multi-UAVs in future 5G IoT scenarios[END_REF][START_REF] Andreou | UAV trajectory optimisation in smart cities using modified A* algorithm combined with haversine and vincenty formulas[END_REF][START_REF] Wang | Taking drones to the next level: Cooperative distributed unmannedaerial-vehicular networks for small and mini drones[END_REF][START_REF] Yang | Multi-UAV-enabled load-balance mobile-edge computing for IoT networks[END_REF][START_REF] Wang | Joint UAV hovering altitude and power control for space-air-ground IoT networks[END_REF][START_REF] Li | Rechargeable multi-UAV aided seamless coverage for QoS-guaranteed IoT networks[END_REF][START_REF] Wang | Energyefficient data collection and device positioning in UAV-assisted IoT[END_REF][START_REF] Huang | Differential evolution with a variable population size for deployment optimization in a UAVassisted IoT data collection system[END_REF][START_REF] Samir | UAV trajectory planning for data collection from time-constrained IoT devices[END_REF][START_REF] Hu | Aoi-minimal trajectory planning and data collection in UAV-assisted wireless powered IoT networks[END_REF][START_REF] Li | Data collection maximization in IoT-sensor networks via an energy-constrained UAV[END_REF][START_REF] Ghdiri | Offline and online UAV-enabled data collection in time-constrained IoT networks[END_REF][START_REF] Liu | DRL-UTPS: DRL-based trajectory planning for unmanned aerial vehicles for data collection in dynamic IoT network[END_REF][START_REF] Wang | Learningbased UAV path planning for data collection with integrated collision avoidance[END_REF].

C. The deployment problem

As shown in Fig. 1(b), the consumed time on the deployment of UAV k is associated with two aspects:

1) Data collection time: Let n IoT,k denote the number of IoT devices assigned to UAV k and (x i , y i , 0) denote the location of IoT device i (i ∈ [1, n IoT,k ]), where x i and y i are the coordinate values in the x-axis and y-axis of the location of IoT device i, respectively. Let n sp,k denote the number of stop points of UAV k to complete the assigned collection task and (X k,j , Y k,j , H k,j ) denote the location of stop point j (j ∈ [1, n sp,k ]), where X k,i and Y k,i , and H k,i are the coordinate values in the x-axis, y-axis, and z-axis of the location of stop point j, respectively.

Like in [START_REF] Huang | Differential evolution with a variable population size for deployment optimization in a UAVassisted IoT data collection system[END_REF], a binary variable c s,i,j is used to record the connection status between IoT device i and stop point j. If c s,i,j is equal to 0, the connection is closed; otherwise, the connection is open. In addition, IoT device i sends its data to the nearest stop point. So, c s,i,j can be written by:

c s,i,j = 1, if j = arg min d ij j∈1,2,...,n sp,k , i ∈ [1, n IoT,k ]; 0, otherwise. (17) n sp,k j=1 c s,i,j = 1, i ∈ [1, n IoT,k ], (18)
n IoT,k i=1 n sp,k j=1 c s,i,j = n IoT,k , (19)
where d i,j is the Euclidean distance between between IoT device i and stop point j. [START_REF] Tran | UAV relay-assisted emergency communications in IoT networks: Resource allocation and trajectory optimization[END_REF] is to ensure that IoT device i only sends data to UAV k once. ( 19) is to ensure that all IoT devices can send data to UAV k. In addition, considering the bandwidth, the number of IoT devices simultaneously sending data to UAV k does not exceed n sd , which can be denoted by:

n IoT,k i=1 c s,i,j ≤ n sd , j ∈ [1, n sp,k ]. ( 20
)
The channel modeling is same with that of [START_REF] Huang | Differential evolution with a variable population size for deployment optimization in a UAVassisted IoT data collection system[END_REF]. In this channel modeling, an IoT device sends its data to a UAV by a wireless channel. The data rate R i,j can be computed by [START_REF] Huang | Differential evolution with a variable population size for deployment optimization in a UAVassisted IoT data collection system[END_REF]:

R ij = B log 1+ P i G ij δ 2 2 = B log 1+ P i G 0 d -2 ij δ 2 2 , ( 21
)
where B is the system bandwidth, δ 2 is the white Gaussian noise power, G i,j is the channel gain between UAV k at stop point j and IoT device i, G 0 is the channel power gain at the reference distance d 0 = 1m, and P i is the transmitting power between IoT device i and UAV k. Thus, when IoT device i sends D i amount of data to UAV k at the stop point j, the required transmission time t tt,i can be obtained by:

t tt,i = D i R ij . ( 22
)
Let t ht,k,j and C ht,k,j denote the hover time and consumed energy of UAV k at stop point j, respectively. According to [START_REF] Han | Joint deployment optimization and flight trajectory planning for UAV assisted IoT data collection: A bilevel optimization approach[END_REF], t ht,k,j and C ht,k,j can be represented by:

t ht,k,j = max(t tt,i ), i ∈ [1, n IoT,k ], (23)
C ht,k,j = t ht,k,j P hp , (24)
where P hp is the hover power of UAV k.

2) Round-trip time: Looking at Fig. 1(b), UAV k needs to return to its DCS to charge its battery during data collection. Let D cs,k , L sp,k , T ht,k , C he,k , and t ctt,k denote the location of the DCS of UAV k, the set of locations of all stop points of UAV k, the set of of the consumed hover time at all stop points of UAV k, the set of of the consumed hover energy at all stop points of UAV k, and the travel time of UAV k during data collection, respectively. Specifically, L sp,k , T ht,k , and C he,k can be denoted by

L sp,k ={(X k,1 , Y k,1 , H k,1 ), (X k,2 , Y k,2 , H k,2 ), . . ., X n sd,k , Y n sd,k , H n sd,k }, T ht,k = {t ht,k,1 , t ht,
k,2 , . . . , t ht,k,nsd }, and C ht,k = {C ht,k,1 , C ht,k,2 , . . . , C ht,k,nsd }, respectively. The computing method of t ctt,k is shown in Algorithm 2, which is described as follows. Firstly, the distances between stop points and the DCS of UAV k are obtained. Then, all stop points are sorted in order from closest to farthest according to the distances. Lastly, UAV k collects the data sequentially according to the order of the sorted stop points. Once the energy of UAV k cannot complete the collection task at the next stop point, it will return to its DCS for charging.

D. Objective function

From the above description, the characteristics of the three stages shown in Fig. 1 can be summarized as follows: travel time, they may need to compete with each other for charging resources in stage 2.

• Stage 3 has the same characteristics as stage 1. Note that, the travel time of a UAV in stage 1 and stage 2 directly affects its travel time in stage 3. Given these characteristics, the considered deployment and FP problem can be seen as a three-layer optimization problem:

• The first layer is to determine the optimal FP of UAVs in stage 1, which is measured by the average travel time of all UAVs. Let L * a denote the optimal L a . Let t adc denote the average travel time of all UAVs in stage 1. To find L * a , the objective function can be defined by:

Minimize t adc = f (L a,k ) = 1 n u nu k=1 t fdc,k , (25)
subject to

L a,k ∈ M CSN , k ∈ [1, n u ]. (26)
• The second layer determines the optimal deployment of UAVs in stage 2, which aims to search for the optimal average travel time of all UAVs. Let t act denote the consumed average time, which can be defined by:

Minimize t act = f (D cs,k , L sp,k , T ht,k , C he,k ) = 1 n u nu k=1 t ctt,k , (27)
subject to ( 17), ( 18), ( 19), [START_REF] Han | Joint deployment optimization and flight trajectory planning for UAV assisted IoT data collection: A bilevel optimization approach[END_REF], and

X min,k,j ≤ X k,j ≤ X max,k,j , j ∈ [1, n sp,k ], (28)
Y min,k,j ≤ Y k,j ≤ Y max,k,j , j ∈ [1, n sp,k ], (29)
Z min,k,j ≤ Z k,j ≤ Z max,k,j , j ∈ [1, n sp,k ], (30)
where X min,k,j , Y min,k,j , and Z min,k,j are the lower limits of the coordinate values of UAV k at stop point j in the a-axis, y-axis, and z-axis, respectively; X max,k,j , Y max,k,j , and Z max,k,j are the upper limits of the coordinate values of UAV k at stop point j in the x-axis, y-axis, and z-axis, respectively.

• The third layer is to determine the optimal FP of all UAVs in stage 3, which is measured by the average travel time of all UAVs. Let L * c denote the optimal L c . Let t acd denote the average travel time of all UAVs in stage 3. To find L * c , the objective function can be defined by:

Minimize t acd = f (L c,k ) = 1 n u nu k=1 t fcd,k , (31)
subject to L c,k ∈ M CSN , k ∈ [1, n u ]. ( 32
)
Let T att denote the average travel time of all UAVs in the three stages. The considered problem can be formulated as:

Minimize T att = t adc + t act + t acd , (33)
subject to ( 17), ( 18), ( 19), ( 20), ( 26), ( 28), ( 29), [START_REF] Rahnamayan | Oppositionbased differential evolution[END_REF], and (32).

C TR = C R -C sp,k,i ; 15: if C TR < 0 then 16: C R = C R - (D spc,k,I d,i-1 +D spc,k,I d,i )P 0 v -C tl +C 0 -C he,
C R = C R -D I P 0 v -C he,k,I d,i ; 26: t ctt,k = t ctt,k + D I v + T ht,k,I d,i . 27:
end if 28:

end if 29: end for Algorithm 3 The method of constructing path by the PBEM

Input:

S (the priority sequence), V p (the set of nodes corresponding to the constructed path), m N (the number of nodes in the network), ms (the source node), m d (the destination node) Output:

V p 1: Initialize V p by V p = [ms] and update S by S(ms) = -∞; 2: for i = 1 to m N do 3:

Select node m k that has the highest priority among the nodes having direct links with node V p(i); 4:

Update

V p by V p = [V p m k ] and S by S(m k ) = -∞; 5: if m k == m d then 6:
Break.

7:

end if 8: end for

V. METHODOLOGY

To solve [START_REF] Ren | A scatter learning particle swarm optimization algorithm for multimodal problems[END_REF], EDDE-DPDE is proposed. This section first introduces DE and then EDDE-DPDE is presented.

A. DE

DE is a popular algorithm, which consists of three parts: 1) Mutation: The mutation operator aims at generating new individuals, which can be expressed by:

v i = x a + M f (x b -x c ), i = 1, 2, . . . , n p , ( 34
)
where n p is the population size, v i is the mutation vector of individual i, M f is the mutation factor, and x a , x b , x c Compute Mp by (38); 8:

for j = 1 to np do 9:

Initialize the flag bit fp by fp = 0; 10:

Generate a, b, and c randomly; 11:

Perform the mutation operator by (37); 12:

Perform the crossover operator by (39); 13:

Compute La corresponding to v i by PBEM; 14:

Evaluate La by (1-16), and ( 25), M CSN , and F fma ; 15:

Perform the selection crossover by [START_REF] Friedman | The use of ranks to avoid the assumption of normality implicit in the analysis of variance[END_REF] and update fp (if v i is better than X dc,i , fp = 1); end for 29: end for are three selected randomly individuals from the population X = {x 1 , x 2 , . . . , x np }, where a, b, and c are three selected integers randomly between 1 and n p .

2) Crossover: The crossover operator is to generate the trial vector z i of v i , which is represented by:

z i,j = v i,j , if d ≤ C f ||j = e; x i,j , otherwise, (35)
where z i,j is the jth variable of z i , v i,j is the jth variable of v i , x i,j is the jth variable of individual i, d is a random number between 0 and 1, C f is the crossover factor, and e is a random integer between 1 and n p . In addition, j meets j ∈ [1, D], where D is the number of variables.

3) Selection operator: The selection operator is to find a better solution between z i and x i , which can be denoted by:

x i = z i , if f (z i ) ≤ f (x i ); x i , otherwise. (36)

B. The proposed EDDE-DPDE

EDDE-DPDE is based on the proposed EDDE and DPDE. Next, the proposed EDDE and DPDE are first introduced. Then, EDDE-DPDE is presented.

1) EDDE: EDDE is proposed to search L * a and L * c from M CSN . To achieve this, this paper adopts a priority-based encoding mechanism (PBEM) [START_REF] Mohemmed | Solving shortest path problem using particle swarm optimization[END_REF] as shown in Algorithm 3, whose basic idea is to use a PS to generate a L a,k or L c,k . In terms of population, each individual consists of n u modules and each module denotes a PS (the length of a PS is equal to the number of nodes in the CSN) of a UAV, which can be expressed as follows: X DC = {X dc,1 , X dc,2 , . . . , X dc,np } and X dc,i = {x dc,i,1 , x dc,i,2 , . . . , x dc,i,nu }, where X DC is a population with n p individuals, x dc,i,k (k ∈ [1, n u ]) is the kth module of individual i in the X dc , and x dc,i,k,j (k ∈ [1, n u ] and j ∈ [1, n n ]) is the jth element of the kth module of individual i. Next, taking X DC as an example to introduce EDDE. EDDE consists of three core components:

• Mutation operator driven by the elite archive. To improve the directionality of the mutation, an elite archive X E is created to save some promising individuals, whose length is equal to n e . That is, X E = {X e,1 , X e,2 , . . . , X e,np }.

The designed mutation operator can be expressed by:

v i = X e,a + φ 1 (X dc,b -X dc,c ) + φ 2 (X e,a -M p ), (37)
where φ 1 and φ 2 are two random numbers between 0 and 1 (φ 1 and φ 2 meet φ 1 + φ 2 = 1), and M p is the mean position of X DC . M p is computed by:

M p = 1 n p np k=1 X dc,i . (38)
From ( 37) and ( 38), the created elite archive can guide the direction of mutation. In addition, the designed mutation operator also can be seen as a random mutation operator that uses random numbers to balance the local mutation factor (the second term on the right of ( 37)) and the global mutation factor (the third term on the right of ( 37)), which can help to keep the population diversity. • Enhanced crossover operator. An individual with a better fitness value usually contains more valuable information to search for the global optimal solution. To make full use of the obtained best individual x * dc , an enhanced crossover operator is introduced, which is expressed by:

z i,j =    v i,j , if d ≤ 1/3; x * dc,j , if d ≥ 2/3; x dc,i,j , otherwise, i ∈ [1, n p ], j ∈ [1, n u n n ].
(39) From (39), v i,j , x * dc,j , and x dc,i,j have the same importance to the designed crossover operator, which does not refer to C f that can help to the stability of EDDE.

• Local escape operator. To increase the chance of EDDE to escape from the local optima, a local escape operator is designed, whose basic idea is to use random numbers generated by opposition-based learning [START_REF] Rahnamayan | Oppositionbased differential evolution[END_REF] to replace some elements of some modules of x * dc . The local escape operator can be expressed by:

x * dc,j,o = -x * dc,j,o , j ∈ [1, ⌈θn u ⌉], o ∈ [1, ⌈ηn n ⌉], (
40) where θ and η are two random numbers between 0 and 1, j and o are two integers. Like [START_REF] Mohemmed | Solving shortest path problem using particle swarm optimization[END_REF], x dc,i,k,j meets x dc,i,k,j ∈ [-n n , n n ] in using PBEM to generate L a,k . Specifically, X DC is initialized by:

x dc,i,k,j = 2n n µ -n n , ( 41
)
where i ∈ [1, n p ], k ∈ [1, n u ], j ∈ [1, n n ],
and µ is a random number between 0 and 1. Algorithm 4 presents the computing method to obtain L * a and t * adc (the optimal t adc ) by EDDE.

Algorithm 5 Dynamic population mechanism Input: X D,k (the population of UAV k), n sp,k (the number of stop points of UAV k), z (trail vector) Output:

X * D,k (the optimal X D,k ) 1: for l = i to nsp do 2:

Get X ID by inserting randomly z i to X D,k ; 3:

Get X RD by using z i to replace randomly an element of X D,k ; 4:

Get X MD by removing randomly an element of X D,k ; 5:

Evaluate X ID , X RD , and X MD by [START_REF] Na | UAV-assisted wireless powered Internet of Things: Joint trajectory optimization and resource allocation[END_REF][START_REF] Tran | UAV relay-assisted emergency communications in IoT networks: Resource allocation and trajectory optimization[END_REF][START_REF] Na | UAV-supported clustered noma for 6G-enabled Internet of Things: Trajectory planning and resource allocation[END_REF][START_REF] Han | Joint deployment optimization and flight trajectory planning for UAV assisted IoT data collection: A bilevel optimization approach[END_REF][START_REF] Dong | Joint optimization of deployment and trajectory in UAV and IRS-assisted IoT data collection system[END_REF][START_REF] Wang | Joint deployment and task scheduling optimization for large-scale mobile users in multi-UAVenabled mobile edge computing[END_REF][START_REF] Pan | Joint power and 3D trajectory optimization for UAV-enabled wireless powered communication networks with obstacles[END_REF][START_REF] Shaikh | A review on swarm intelligence and evolutionary algorithms for solving the traffic signal control problem[END_REF] and Algorithm 2; 6:

Update X D,k by choosing the best of X ID , X RD , and X MD ; 7: end for 8: Get X * D,k by X * D,k =X D,k .

2) DPDE: DPDE is to compute the optimal number and locations of stop points of UAV k to complete the assigned data collection task. As done in [START_REF] Huang | Differential evolution with a variable population size for deployment optimization in a UAVassisted IoT data collection system[END_REF], X D,k in DPDE denotes an entire deployment of UAV k. X D,k can be expressed by

X D,k = {X d,k,1 , X d,k,2 , . . . , X d,k,n sp,k }, where X d,k,i = {(X k,i , Y k,i , H k,i )}, i ∈ [1, n sp,k ].

DPDE includes three parts:

• Mutation operator based on the historical population. As shown in [START_REF] Rao | Jaya: A simple and new optimization algorithm for solving constrained and unconstrained optimization problems[END_REF], there is only a mutation strategy, which does not help to keep the population diversity [START_REF] Rakhshani | Snap-drift cuckoo search: A novel cuckoo search optimization algorithm[END_REF]. Given this, a set of mutation operators with historical population is designed, which can be written by:

v i =                X d,k,a + ϕ 1 (X d,k,b -X d,k,c )+ ϕ 2 (X d,k,i -M d ), if κ ≤ 1 3 , X d,k,a + ϕ 1 (X h,k,g -X d,k,c )+ ϕ 2 (X d,k,i -M d ), if κ ≤ 2 3 , X d,k,a + ϕ 1 (X h,k,g -X h,k,r )+ ϕ 2 (X d,k,i -M d ), otherwise, (42)
where a, b, and c are three random numbers that meet a, b, and c ∈ [1, n sp,k ], g and r are two random integers, ϕ 1 and ϕ 2 are two random numbers between 0 and 1 (ϕ 1 and ϕ 2 meet ϕ 1 +ϕ 2 = 1), M d is the mean position of X D,k , X h,k,g and X h,k,r are the gth and the rth elements of X H,k , respectively. Let n h,k denote the length of X H,k . g and r meet g, r ∈ [1, n h,k ]. M d can be computed by:

M d = 1 n sp,k n sp,k k=1 X d,k,i . (43)
In (42), the three learning strategies have the same importance and are assigned the same selected probability. Like ( 37), (42) also can be seen as a random mutation operator, which employs random numbers to balance the local mutation factor (the second term on the right of each mutation strategy) and the global mutation factor (the third term on the right of each mutation strategy). This can help to increase the population diversity. • Random crossover operator. Using random numbers to replace crossover factor is the basic idea of the designed random crossover operator that can be expressed by:

z i,j = v i,j , if d 1 ≤ d 2 ||j = e; X d,k,i,j , otherwise. ( 44
)
• Dynamic population mechanism. To obtain the optimal number and locations of UAV k, a dynamic population mechanism [START_REF] Huang | Differential evolution with a variable population size for deployment optimization in a UAVassisted IoT data collection system[END_REF] is employed, which consists of one insert t * act (the optimal tact) 1: for k = 1 to nu do 2:

Initialize X D,k and X H,k by (28-30); 3:

Evaluate X D,k by [START_REF] Na | UAV-assisted wireless powered Internet of Things: Joint trajectory optimization and resource allocation[END_REF][START_REF] Tran | UAV relay-assisted emergency communications in IoT networks: Resource allocation and trajectory optimization[END_REF][START_REF] Na | UAV-supported clustered noma for 6G-enabled Internet of Things: Trajectory planning and resource allocation[END_REF][START_REF] Han | Joint deployment optimization and flight trajectory planning for UAV assisted IoT data collection: A bilevel optimization approach[END_REF][START_REF] Dong | Joint optimization of deployment and trajectory in UAV and IRS-assisted IoT data collection system[END_REF][START_REF] Wang | Joint deployment and task scheduling optimization for large-scale mobile users in multi-UAVenabled mobile edge computing[END_REF][START_REF] Pan | Joint power and 3D trajectory optimization for UAV-enabled wireless powered communication networks with obstacles[END_REF][START_REF] Shaikh | A review on swarm intelligence and evolutionary algorithms for solving the traffic signal control problem[END_REF] and Algorithm 2; 4:

Update X * D,k by X * D,k =X D,k ; 5:

for i = 1 to M Iter2 do 6:

Generate a random number ω between 0 and 1; 7:

if ω ≤ 0.5 then 8:

X H,k = X D,k ; 9:

end if 10:

Update X H,k by randomly sorting the elements in the X H,k ; 11:

Perform the mutation operator by ( 42) and ( 43 (the optimal t ctt,k ) by [START_REF] Na | UAV-assisted wireless powered Internet of Things: Joint trajectory optimization and resource allocation[END_REF][START_REF] Tran | UAV relay-assisted emergency communications in IoT networks: Resource allocation and trajectory optimization[END_REF][START_REF] Na | UAV-supported clustered noma for 6G-enabled Internet of Things: Trajectory planning and resource allocation[END_REF][START_REF] Han | Joint deployment optimization and flight trajectory planning for UAV assisted IoT data collection: A bilevel optimization approach[END_REF][START_REF] Dong | Joint optimization of deployment and trajectory in UAV and IRS-assisted IoT data collection system[END_REF][START_REF] Wang | Joint deployment and task scheduling optimization for large-scale mobile users in multi-UAVenabled mobile edge computing[END_REF][START_REF] Pan | Joint power and 3D trajectory optimization for UAV-enabled wireless powered communication networks with obstacles[END_REF][START_REF] Shaikh | A review on swarm intelligence and evolutionary algorithms for solving the traffic signal control problem[END_REF] and Algorithm 2; 16: end for 17: Get t * act by [START_REF] Kim | Optrone: Maximizing performance and energy resources of drone batteries[END_REF].

Algorithm 7 The implementation of EDDE-DPDE

Input:

np (population size), nu (the number of UAVs), nn (the number of nodes), M CSN (the connected matrix), M Iter1 (the maximum number of iterations in stage 1), M Iter2 (the maximum number of iterations in stage 2), F fma (the flight schedule in stage 1), F fmd (the flight schedule in stage 3), M Iter3 (the maximum number of iterations in stage 3) Output:

The optimal average travel time t if tatt ≤ t * att then 8:

t * att = tatt. 9:

end if 10: end for operator, one remove operator, and one replace operator. Algorithm 5 shows the flowchart of this mechanism. In [START_REF] Huang | Differential evolution with a variable population size for deployment optimization in a UAVassisted IoT data collection system[END_REF], n sp,k is initialized by n sp,k = n IoT,k . To accelerate the convergence, n sp,k is initialized by n sp,k = 0.5n IoT,k . The implementation of DPDE has been presented in Algorithm 6.

3) EDDE-DPDE: Algorithm 7 shows the implementation of EDDE-DPDE. According to Algorithm 7, EDDE-DPDE has a three-layer structure: 1) the first one (line 3 in Algorithm 7) is to determine the optimal FP in stage 1, 2) the second one (line 4 in Algorithm 7) is to determine the optimal deployment in stage 2, and 3) the third one (line 5 in Algorithm 7) is to determine the optimal FP in stage 3.

VI. RESULTS AND COMPARISONS

A. Experiment Preparation Fig. 2 shows the considered scenario of using multi-UAVs for long-distance data collection from large-scale IoT devices according to the geographic information of Hong Kong. The considered scenario consists of a DC, 20 CSs, and a data collection zone (which is surrounded by a black box). The method of determining DCSs is as follows. That is, the data collection zone can be regarded as a rectangle with a length of 6 Km and a width of 3 Km. Specifically, as shown in Fig. 3, the data collection zone can be divided into 18 data collection units (DCUs), and each DCU can be seen as a square with a side length of 1 Km. Each DCU corresponds to a UAV, whose DCS is the CS closest to the center of the DCU. Thus, 18 UAVs are required, whose DCSs are listed in Table II.

In addition, the parameters of a UAV and IoT devices are extracted from [START_REF] Huang | Differential evolution with a variable population size for deployment optimization in a UAVassisted IoT data collection system[END_REF] and [START_REF] Huang | Drone stations-aided beyondbattery-lifetime flight planning for parcel delivery[END_REF], which are as follows: C 0 = 3.20 × 10 5 J, C tl = 1.25 × 10 4 J, t tl = 50s, v = 4m/s, P 0 = 300W, L max = 4.27Km, the number of IoT devices in each DCU is set to 300, the amount of data of each IoT device is between 500 MB and 1500 MB, G 0 = -30dB, δ 2 = -250dB, n sd = 10, B = 1MHz, P hp = 1000W, and P i = 0.1W.

To verify the performance of EDDE-DPDE, it is compared with DE-DEVIPS, BA-DEVIPS, PSO-DEVIPS, and JAYA-DEVIPS as shown in Table III. In Table III, DE, bat algorithm (BA) [START_REF] Yang | Chapter 10 -bat algorithms[END_REF], particle swarm optimization (PSO) [START_REF] Ren | A scatter learning particle swarm optimization algorithm for multimodal problems[END_REF], and Jaya algorithm (JAYA) [START_REF] Rao | Jaya: A simple and new optimization algorithm for solving constrained and unconstrained optimization problems[END_REF] are four popular population-based metaheuristic algorithms. Like DPDE, differential evolution with a variable population size (DEPS) [START_REF] Huang | Differential evolution with a variable population size for deployment optimization in a UAVassisted IoT data collection system[END_REF] is a also dynamic population differential evolution. In addition, the population size is set to 50 for EDDE, DE, BA, PSO, and JAYA. From Algorithm 4, each individual in EDDE is evaluated twice in one loop. To make a fair comparison, the maximum number of iterations is set to 250 and 500 for EDDE and the compared algorithms (i.e., DE, BA, PSO, and JAYA), respectively. The maximum number of iterations of DPDE and DEPS is set to 8,000. Other parameters of the compared algorithms are from the corresponding references. n c is set to 3 at each node and the number of independent runs for each algorithm is set to 30. Besides, the following quality indicators are employed:

• Value-based indicator. This indicator is to evaluate the solution accuracy, which includes the best solution (BEST), the mean solution (MEAN), the worst solution (WORST), the standard deviation (STD), and the average improvement rate (AIR). In terms of BEST, MEAN, and WORST, the smaller the value, the higher the accuracy of the solution obtained by an algorithm. A smaller STD means that an algorithm has stronger solution stability. Specifically, BEST, MEAN, WORST, and STD can be obtained by the statistical results of 30 independent runs. In addition, AIR AB (the AIR of Algorithm A relative to Algorithm B) can be defined by:

AIR AB = M B -M A M A × 100%, (45)
where M A and M B are the MEAN of Algorithm A and Algorithm B, respectively. From (45), the larger the value of AIR AB , the greater the advantage of Algorithm A over Algorithm B. Tables IV-VII present the results of the value-based indicators, and the best results are in bold. • Significance-based indicator. This indicator is to determine whether there are significant differences between the results obtained by the applied algorithms. To achieve this, a non-parametric Wilcoxon signed-rank test with a significance level α = 0.05 [START_REF] Derrac | A practical tutorial on the use of nonparametric statistical tests as a methodology for comparing evolutionary and swarm intelligence algorithms[END_REF] has been conducted for the mean results of 30 runs for each algorithm. The results produced by Wilcoxon signed-rank test have been displayed in Table VIII. In Table VIII, R + is the sum of positive rank, R -is the sum of negative rank, P is the probability, and '+' means the null hypothesis is rejected (P ≤ 0.05) with 95% certainty (the results of EDDE-DPDE and the compared algorithm have significant differences, and EDDE-DPDE shows better performance). • Rank-based indicator. This indicator is to rank all algorithms. Friedman test [START_REF] Friedman | The use of ranks to avoid the assumption of normality implicit in the analysis of variance[END_REF] is a very powerful nonparametric statistical test, which is often used to compare the differences among algorithms [START_REF] Liu | Adaptive particle swarm optimizer combining hierarchical learning with variable population[END_REF], [START_REF] Pan | Improved meta-heuristics for solving distributed lot-streaming permutation flow shop scheduling problems[END_REF]. In this paper, the mean results of the applied algorithms are subjected to the Friedman test, and the Friedman mean rank (FMR) is recorded. Specifically, an algorithm with a smaller FMR means it has a better performance. The results of FMR have been shown in Table IX. IX, all the applied algorithms can be sorted in the following order (from best to worst): EDDE-DPDE, DE-DEPS, JAYA-DEPS, PSO-DEPS, and BA-DEPS. EDDE-DPDE is the best of all algorithms.

When the UAV-station model is employed to add the energy of a UAV, charging time will take up a large proportion of the total travel time. As discussed above, EDDE-DPDE shows obvious advantages over the compared algorithms in terms of saving the charging time, which proves the great potential of EDDE-DPDE to be applied to determine the flight planning of a UAV to carry out a long-distance flight task under the UAV-station model. Like charging time, the straight flight time also takes up a large proportion of the total travel time under the UAV-station model. From the above discussion, EDDE-DPDE shows better comprehensive performance than DE-DEPS, BA-DEPS, PSO-DEPS, and JAYA-DEPS in terms of reducing the average straight flight time. This supports the superiority of EDDE-DPDE in determining the flight planning of a UAV to execute a long-distance flight task under the UAV-station model.

VII. CONCLUSION

This paper aims to use multi-UAVs to collect data from large-scale IoT devices. We design a multi-UAVs-assisted large-scale IoT data collection system that can support longdistance data collection. To make the system work efficiently, a new population-based optimization algorithm, namely EDDE-DPDE, is reported to solve the joint optimization problem of FP and deployment for multi-UAVs. EDDE-DPDE has a three-layer structure. The first and third layers are EDDE, which is to determine the optimal FP of multi-UAVs between the DC and their DCSs. The second layer is DPDE, which is to optimize the deployment of multi-UAVs in the data collection zone. The EDDE is based on three improved strategies: mutation operator driven by the created elite archive, enhanced crossover operator, and local escape operator. The DPDE is based on two improved strategies: mutation operator with historical population and random crossover operator. To investigate the performance of EDDE-DPDE, a data collection scenario based on real geographic information is designed. Numerical results demonstrate that EDDE-DPDE shows an improvement of at least 11.11% compared with four powerful algorithms in terms of average travel time, which proves the validity of the improved strategies.

In this paper, data collection units are divided based on the area of the data collection zone. In future research, we will pay attention to designing an adaptive task allocation method based on the data amount and location distribution of IoT devices to further improve the work efficiency of this system.

  This work was supported by Department General Research Grant [P0040253] and funding from the Research Institute for Sports Science and Technology [P0043566]. (Corresponding author: Chao Huang.) Y. Zhang and C. Huang are with the Department of Industrial and Systems Engineering, the Hong Kong Polytechnic University, Hong Kong. (E-mail: {yiyingchreo.zhang, hchao.huang}@polyu.edu.hk).

3 :

 3 Stage 1: UAVs fly from the DC to their DCSs. UAVs fly from their DCSs to the DC.

  (a), and 2) the UAV flies from its DCS to the DC as shown in Fig. 1(c). According to Fig. 1, when a UAV flies from one node to another, its FP can be seen as a sequence consisting of some nodes. Let L a = {L a,1 , L a,2 , . . . , L a,nu } and L c = {L c,1 , L c,2 , . . . , L c,nu } denote the flights of all UAVs in stage 1 and stage 3, respectively. Specifically, L a,k and L c,k are the flight of UAV k in the stage 1 and stage 3, respectively. L a,k and L c,k can be denoted by:

•

  In stage 1, a UAV flies from the DC to its DCS, whose flight is constrained by its battery capability and a limited number of charging piles (it needs to compete with other UAVs for charging resources) at one node. • In stage 2, a UAV needs to be recharged multiple times to complete the data collection task at multiple stop points. Therefore, once DCSs are identified, the impact of stage 1 on stage 2 includes: 1) stage 1 determines the start time of stage 2 and 2) when multiple UAVs share the same DCS and their flights determined in stage 1 consume similar This article has been accepted for publication in IEEE Internet of Things Journal. This is the author's version which has not been fully edited and content may change prior to final publication. Citation information: DOI 10.1109/JIOT.2023.3285942 © 2023 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission. See https://www.ieee.org/publications/rights/index.html for more information. Authorized licensed use limited to: UNIVERSITE POLYTECHNIQUE HAUTS DE FRANCE. Downloaded on June 20,2023 at 20:04:56 UTC from IEEE Xplore. Restrictions apply.

Algorithm 2

 2 Calculate t ctt,k Input: D cs,k , L sp,k , T ht,k , and C he,k Output: t ctt,k 1: Initialize t ctt,k by t ctt,k = 0; 2: Get D spc,k (Euclidean distance between D cs,k and L sp,k ); 3: Sort D spc,k from small to large and get the index I d ; 4: for i = 1 to n sd do 5: if i == 1 then 6: C sp,k,i = C he,I d,k,i + 2D spc,k,I d,i D I (Euclidean distance of D spc,k,I d,i and D spc,k,I d,i-1 ); 9: C sp,k,i = D I P 0 v + C he,k,I d,i + D spc,k,I d,i end for 12: C R = C 0 ; 13: for i = 1 to n sd do 14:

Algorithm 4 1 :

 41 The implementation of EDDE Input: np (population size), nu (the number of UAVs), nn (the number of nodes), F fma (the flight schedule in stage 1), M CSN (the connected matrix), and N Iter1 (the maximum number of iterations in stage 1) Output: L * a and t * adc Initialize the population X DC by (41); 2: Compute La corresponding to each individual in the X DC by PBEM; 3: Evaluate each La by (1-16), and (25), M CSN , and F fma ; 4: Find x * dc , L * a , and t * adc ; 5: Initialize X E by X E = x * dc and L E (the length of X E ) by L E = 1; 6: for l = 1 to N Iter1 do 7:

Algorithm 6

 6 The implementation of DPDE Input:nu (the number of UAVs), M Iter2 (the maximum number of iterations in stage 2) Output:

Fig. 2 .

 2 Fig. 2. The considered scenario of using multi-UAVs for longdistance data collection from large-scale IoT devices.

Fig. 3 .

 3 Fig. 3. The divided data collection units. The data collection zone can be described by four pairs of longitude and latitude: (22.247781 • N, 114.162139 • E), (22.24778 • N, 114.21976 • E), (22.27475 • N, 114.21976 • E), and (22.27475 • N, 114.16214• E). That is, the data collection zone can be regarded as a rectangle with a length of 6 Km and a width of 3 Km. Specifically, as shown in Fig.3, the data collection zone can be divided into 18 data collection units (DCUs), and each DCU can be seen as a square with a side length of 1 Km. Each DCU corresponds to a UAV, whose DCS is the CS closest to the center of the DCU. Thus, 18 UAVs are required, whose DCSs are listed in TableII.In addition, the parameters of a UAV and IoT devices are extracted from[START_REF] Huang | Differential evolution with a variable population size for deployment optimization in a UAVassisted IoT data collection system[END_REF] and[START_REF] Huang | Drone stations-aided beyondbattery-lifetime flight planning for parcel delivery[END_REF], which are as follows: C 0 = 3.20 × 10 5 J, C tl = 1.25 × 10 4 J, t tl = 50s, v = 4m/s, P 0 = 300W, L max = 4.27Km, the number of IoT devices in each DCU is set to 300, the amount of data of each IoT device is between 500 MB and 1500 MB, G 0 = -30dB, δ 2 = -250dB, n sd = 10, B = 1MHz, P hp = 1000W, and P i = 0.1W.To verify the performance of EDDE-DPDE, it is compared with DE-DEVIPS, BA-DEVIPS, PSO-DEVIPS, and JAYA-DEVIPS as shown in TableIII. In TableIII, DE, bat algorithm (BA)[START_REF] Yang | Chapter 10 -bat algorithms[END_REF], particle swarm optimization (PSO)[START_REF] Ren | A scatter learning particle swarm optimization algorithm for multimodal problems[END_REF], and Jaya algorithm (JAYA)[START_REF] Rao | Jaya: A simple and new optimization algorithm for solving constrained and unconstrained optimization problems[END_REF] are four popular population-based metaheuristic algorithms. Like DPDE, differential evolution with a variable population size (DEPS)[START_REF] Huang | Differential evolution with a variable population size for deployment optimization in a UAVassisted IoT data collection system[END_REF] is a also dynamic population differential evolution. In addition, the population size is set to 50 for EDDE, DE, BA, PSO, and JAYA. From Algorithm 4, each individual in EDDE is evaluated twice in one loop. To make a fair comparison, the maximum number of iterations is set to 250 and 500 for EDDE and the compared algorithms (i.e., DE, BA, PSO, and JAYA), respectively. The maximum number of iterations of DPDE and DEPS is set to 8,000. Other parameters of the compared algorithms are from the corresponding references. n c is set to 3 at each node and

Fig. 4 .

 4 Fig. 4. 12 optimal flights from six UAVs obtained by EDDE-DPDE (the direction of the arrow is the movement direction of the UAV between two nodes).

  the second best BEST, MEAN, WORST, and STD. PSO-DEPS is the worst of all algorithms in terms of BEST, MEAN, and WORST. Specifically, in terms of AIR, EDDE-DPDE is better than DE-DEPS, BA-DEPS, PSO-DEPS, and JAYA-DEPS by 11.11%, 15.68%, 15.69%, and 14.53%, respectively. This indicates that EDDE-DPDE can provide a more efficient solution to the considered joint optimization problem compared with DE-DEPS, BA-DEPS, PSO-DEPS, and JAYA-DEPS. In addition, as shown in the fourth row of Table VIII, the test results of the Wilcoxon signed-rank test also support that EDDE-DPDE outperforms DE-DEPS, JAYA-DEPS, PSO-DEPS, and BA-DEPS. According to the obtained FMR in the second row of Table IX, all the applied algorithms can be sorted in the following order (from best to worst): EDDE-DPDE, DE-DEPS, JAYA-DEPS, PSO-DEPS, and BA-DEPS. EDDE-DPDE is the best of all algorithms. Besides, Fig. 4 presents 12 optimal flights from six UAVs that EDDE-DPDE searches from all candidate flights under the premise of consuming the least average travel time. By observing Fig. 4, the optimal flights of six UAVs obtained by EDDE-DPDE in stage 1 and stage 3 are not much different, which can be explained as follows: • In stage 1, UAVs fly from the DC to six CSs at the same time and starting point. Thus, although UAVs fly towards six CSs, UAVs will face great competitive pressure for charging resources in the early stage of the flight. • In stage 3, UAVs fly toward the DC. Note that, there is a certain time interval for the departure time of UAVs in stage 3 due to the spent different times by each UAV in stage 1 and stage 2, which can help to reduce the competitive pressure of UAVs for charging resources. • The considered objective function focuses on the average travel time of all UAVs rather than a single UAV. The above discussion supports the superiority of EDDE-DPDE in solving the considered joint optimization problem and the effectiveness of the improved strategies. C. Case 2: comparison on average waiting time Average waiting time is the consumed average time by all UAVs waiting for charging in the process of executing the task. Specifically, an efficient algorithm should reduce the consumed average waiting time as much as possible on the premise of completing the task. Table V shows the statistical results of the average waiting time achieved by five algorithms. Looking at Table V, although DE-DEPS can get the optimal BEST, DE-DEPS cannot compete with EDDE-DPDE on MEAN, WORST, and STD. DE-DEPS and EDDE-DPDE can achieve similar solutions on BEST. Note that, DE-DEPS is the worst of all algorithms in terms of STD, which means the stability of DE-DEPS is the worst of all algorithms. In addition, EDDE-DPDE shows significant advantages over BA-DEPS, PSO-DEPS, and JAYA-DEPS on BEST, MEAN, WORST, and STD. Specifically, according to AIR, EDDE is better than DE-DEPS, BA-DEPS, PSO-DEPS, and JAYA-DEPS by 30.84%, 44.88%, 24.59%, and 28.84%, respectively. That is, EDDE-DEPS is significantly superior to DE-DEPS, BA-DEPS, PSO-DEPS, and JAYA-DEPS in terms of saving the average waiting time. Besides, according to the fifth row of Table VIII, the test results from the Wilcoxon signed-rank test also prove that the compared algorithms are inferior to EDDE-DPDE. In terms of the obtained FMR presented in the third row of Table IX, all the applied algorithms can be sorted in the following order (from best to worst): EDDE-DPDE, PSO-DEPS, JAYA-DEPS, DE-DEPS, and BA-DEPS. EDDE-DPDE is the best of all algorithms.From the above discussion, EDDE-DPDE shows obvious comprehensive advantages over DE-DEPS, BA-DEPS, PSO-DEPS, and JAYA-DEPS in terms of reducing average waiting time. This advantage is very helpful for EDDE-DPDE to be used for determining the flight planning of a UAV to execute a long-distance flight task under the UAV-station model.D. Case 3: comparison on average charging timeAverage charging time is the consumed average time by all UAVs charging in the process of executing the task. An efficient algorithm should reduce average charging time as much as possible on the premise of completing the task.

E. Case 4 :

 4 comparison on average straight flight time Average straight flight time is the consumed average time by all UAVs on the straight flight during the task. An efficient algorithm should reduce average straight flight time as much as possible on the premise of completing the task. The statistical results of average straight flight time from five algorithms have been shown in Table VII. As presented in Table VII, in terms of BEST, MEAN, WORST, and STD, EDDE-DPDE is the best of all algorithms, which proves its excellent global search ability and stability. DE-DEPS shows slight advantages over BA-DEPS, PSO-DEPS, and JAYA-DEPS on BEST, MEAN, WORST, and STD. From AIR, EDDE-DPDE is better than DE-DEPS, BA-DEPS, PSO-DEPS, and JAYA-DEPS by 3.89%, 3.90%, 3.90%, and 3.90%, respectively. In addition, as displayed in the seventh row of Table VIII, the test results of the Wilcoxon signed-rank test support that DE-DEPS, BA-DEPS, PSO-DEPS, and JAYA-DEPS cannot compete with EDDE-DPDE. According to the

  att 1: Initialize t * att by t * att = +∞; 2: for l = i to M Iter do

	3: 4: 5: 6:	Get t * adc by Algorithm 3; Get t * act by Algorithm 5; Get t * acd by Algorithm 3; Compute tatt by (33);
	7:

TABLE II :

 II The task scheduler of UAVs.

	No.	DCS	No.	DCS	No.	DCS
	UAV 1	CS 11	UAV 7	CS 15	UAV 13	CS 19
	UAV 2	CS 10	UAV 8	CS 14	UAV 14	CS 18
	UAV 3	CS 10	UAV 9	CS 14	UAV 15	CS 18
	UAV 4	CS 11	UAV 10	CS 15	UAV 16	CS 19
	UAV 5	CS 11	UAV 11	CS 15	UAV 17	CS 19
	UAV 6	CS 10	UAV 12	CS 14	UAV 18	CS 18

TABLE III :

 III The structures of the applied algorithms.

	Algorithm EDDE-DPDE DE-DEPS BA-DEPS PSO-DEPS JAYA-DEPS
	Layer-1	EDDE	DE	BA	PSO	JAYA
	Layer-2	DPDE	DEVIPS DEVIPS DEVIPS	DEVIPS
	Layer-3	EDDE	DE	BA	PSO	JAYA

TABLE IV :

 IV The statistical results of average travel time from 30 independent runs (the unit is seconds).

	Indicator DE-DEPS BA-DEPS PSO-DEPS JAYA-DEPS EDDE-DPDE
	BEST	40932.72 42652.12 42912.45	42587.22	38083.76
	MEAN 44737.83 46577.89 46583.57	46115.18	40265.08
	WORST 47568.01 50585.71 51197.28	50809.4	42427.76
	STD	1641.803 2136.288 2126.387	2134.692	1134.90
	AIR	11.11%	15.68%	15.69%	14.53%	-

TABLE V :

 V The statistical results of average waiting time from 30 independent runs (the unit is seconds). Average travel time is the consumed average time by all UAVs to complete the task. Clearly, the smaller the consumed average travel time, the more efficient the algorithm.The statistical results of average travel time obtained by five algorithms have been presented in TableIV. By observing TableIV, EDDE-DPDE is the best of all algorithms in terms of BEST, MEAN, WORST, and STD, which means that EDDE-DPDE shows stronger global search ability and stability. DE-DEPS also can provide a competitive solution, which achieves

	Indicator DE-DEPS BA-DEPS PSO-DEPS JAYA-DEPS EDDE-DPDE
	BEST	1064.42	1442.26	1147.91	1201.17	1083.79
	MEAN 1808.08	2002.03	1721.78	1780.47	1381.87
	WORST 2717.99	2763.24	2488.08	2654.08	1775.76
	STD	395.86	376.41	370.45	357.14	161.64
	AIR	30.84%	44.88%	24.59%	28.84%	-
	B. Case 1: comparison on average travel time

TABLE VI :

 VI The statistical results of average charging time from 30 independent runs (the unit is seconds).

	Indicator DE-DEPS BA-DEPS PSO-DEPS JAYA-DEPS EDDE-DPDE
	BEST	23139.60 23970.37 24122.31	23882.87	21580.08
	MEAN 25411.36 26414.71 26488.58	26193.71	23748.68
	WORST 27188.22 28897.30 29255.35	29312.66	2441.39
	STD	1017.73	1263.38	1246.43	1308.29	717.72
	AIR	7.00%	11.23%	11.54%	10.29%	-

TABLE VII :

 VII The statistical results of average straight flight time from 30 independent runs (the unit is seconds).

	Indicator DE-DEPS BA-DEPS PSO-DEPS JAYA-DEPS EDDE-DPDE
	BEST	5484.75	5484.87	5485.67	5485.44	5269.16
	MEAN 5848.37	5848.98	5849.14	5848.95	5629.52
	WORST 6249.64	6250.62	6251.02	6251.01	6033.56
	STD	855.74	856.31	856.25	856.52	830.26
	AIR	3.89%	3.90%	3.90%	3.90%	-

  Table VI presents the statistical results of average charging time from five algorithms. From Table VI, EDDE-DPDE can get the optimal BEST, MEAN, WORST, and STD, which means EDDE-DPDE outperforms DE-DEPS, BA-DEPS, PSO-DEPS, and JAYA-DEPS in terms of global search ability and stability. DE-DEPS also shows excellent performance, which is second only to EDDE-DPDE in terms of BEST, MEAN, WORST, and STD. According to AIR, EDDE-DPDE is better than DE-DEPS, BA-DEPS, PSO-DEPS, and JAYA-DEPS by 7.00%, 11.23%, 11.54%, and 10.29%, respectively. In addition, looking at the sixth row of Table VIII, the test results of the Wilcoxon signed-rank test determine that EDDE-DPDE can provide better solutions than DE-DEPS, BA-DEPS, PSO-DEPS, and JAYA-DEPS. According to the obtained FMR shown in the fourth row of Table

TABLE VIII :

 VIII The results of Wilcoxon signed-rank test.

		EDDE-DPDE vs.
	No.	DE-DEPS		BA-DEPS		PSO-DEPS		JAYA-DEPS
		R +	R -	P	R +	R -	P	R +	R -	P	R +	R -	P
	Case 1	465	0	1.73E -6 +	465	0	1.73E -6 +	465	0	1.73E -6 +	465	0	1.73E -6 +
	Case 2	453	12	5.75E -6 +	465	0	1.73E -6 +	432	33	4.07E -5 +	450	15	7.69E -6 +
	Case 3	465	0	1.73E -6 +	465	0	1.73E -6 +	465	0	1.73E -6 +	465	0	1.73E -6 +
	Case 4	465	0	1.73E -6 +	465	0	1.73E -6 +	465	0	1.73E -6 +	465	0	1.73E -6 +

TABLE IX :

 IX The obtained FMR by executing Friedman test. FMR listed in the fifth row of TableIX, all algorithms can be sorted in the following order (from best to worst): EDDE-DPDE, DE-DEPS, JAYA-DEPS, BA-DEPS, and PSO-DEPS. That is, EDDE-DPDE is the best of all algorithms.

	No.	DE-DEPS BA-DEPS PSO-DEPS JAYA-DEPS EDDE-DPDE
	Case 1 2.20	4.20	4.00	3.60	1.00
	Case 2 3.37	4.13	2.87	3.17	1.47
	Case 3 2.23	4.27	4.07	3.43	1.00
	Case 4 2.33	3.90	4.10	3.67	1.00

obtained

This article has been accepted for publication in IEEE Internet of Things Journal. This is the author's version which has not been fully edited and content may change prior to final publication. Citation information: DOI 10.1109/JIOT.2023.3285942 © 2023 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission. See https://www.ieee.org/publications/rights/index.html for more information. Authorized licensed use limited to: UNIVERSITE POLYTECHNIQUE HAUTS DE FRANCE. Downloaded on June 20,2023 at 20:04:56 UTC from IEEE Xplore. Restrictions apply.

This article has been accepted for publication in IEEE Internet of Things Journal. This is the author's version which has not been fully edited and content may change prior to final publication. Citation information: DOI 10.1109/JIOT.2023.3285942 © 2023 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission. See https://www.ieee.org/publications/rights/index.html for more information.