一类非自治系统混沌细胞模型研究

卢元元; 丘水生; 周小安

一类非自治系统混沌细胞模型研究

卢元元^{①;丘水生},
丘水生,
周小安

计量
- 文章访问数: 2055
- HTML全文浏览量: 125
- PDF下载量: 787
- 被引次数: 1
出版历程
- 收稿日期: 2000-02-03
- 修回日期: 2000-10-09
- 刊出日期: 2002-04-19

Study on chaotic cell model in a kind of nonautonomous systems

摘要

摘要: 该文研究非自治系统混沌吸引子,通过对Genesio系统深入的数学分析和数字仿真,得出一类非自治系统的混沌细胞模型。
- 非自治混沌系统; 鞍焦平衡点; 半稳定极限环; 混沌模型
Abstract: This paper studies chaotic attractors in nonautonomous nonlinear systems. Based on the mathematical analysis and digital simulation of the Genesio system in detail, a cell model of chaotic attractor in a kind of nonautonomous systems is obtained.

HTML全文

1. 引言

设备到设备(Device to Device, D2D)通信无需安装任何基础设施就能支持各种高速服务，已经成为第5代移动通信技术(5th Generation mobile communication technology, 5G)的一个重要功能^[1]。D2D用户复用蜂窝系统资源已经受到了极大关注^[2,3]。无人机拥有强视距通信链路和自适应高度等优势，能够为高密度城市区域中的无线通信用户提供服务。目前对无人机通信的研究还处于快速发展阶段^[4,5]，特别是其与地面D2D通信网络的共存引起了大量关注^[1-3,6,7]。无人机的灵活性为与地面设备共享同一频段提供了机会，而频谱共享也带来了一些有趣挑战，如功率控制和轨迹设计等^[6]。因此，D2D辅助的无人机通信(D2D Assisted Unmanned Aerial Vehicle communication, DAUAV)系统迫切需要先进频谱共享技术来提高频谱效率。

由于授权频谱资源的紧张，目前高通、三星和威瑞森等公司已经考虑调用免授权频谱来增加可利用的带宽。第3代合作伙伴计划提出的工作在免授权频段的5G空中接口(5G New Radio, 5G NR)能将5G NR扩展到免授权频谱来增加可利用的带宽。因此，本文提出一个基于全频谱共享的DAUAV系统。在这个系统中，无人机能通过控制无人机3维轨迹和发射功率管理不同用户之间的干扰，从而使得上行蜂窝用户、D2D用户和WiFi用户能使用相同的免授权信道，同时多个不同下行蜂窝用户能够使用相同的授权信道。

一些研究人员研究了无人机通信系统中授权频谱共享方案的性能。文献[8]采用主用户相同频段为地面次用户提供数据传输，通过联合优化无人机3维飞行轨迹和发射功率，在不干扰主用户的前提下实现了次用户总吞吐量最大化。文献[9]考虑多无人机同时进行上行和下行传输，通过联合优化无人机3维轨迹、通信调度和传输功率最大化系统的总吞吐量。文献[10,11]在保证主用户通信质量的前提下实现了无人机认知安全通信系统的安全速率最大化。文献[12]研究了能量采集的认知移动中继网络并实现了次用户吞吐量最大化。

在DAUAV系统中，D2D用户和蜂窝用户大多采用频谱共享方式进行共存。文献[13]通过联合优化功率和无人机飞行轨迹，在保证D2D用户最低通信速率和无人机传输信息约束的同时，最大化无人机传输给目的节点的总吞吐量。文献[3]在确保D2D用户服务质量的前提下实现了最小蜂窝用户的吞吐量的最大化。文献[1]在基于能量收集的无人机D2D通信网络中提出了一种鲁棒性资源分配算法。文献[7]通过联合优化功率和信道分配在保证蜂窝用户最小传输速率和D2D用户最小安全速率前提下实现了D2D用户能效的最大化。

最近，无人机免授权通信在缓存、虚拟现实和灾难恢复等应用中获得了越来越多的关注^[14-19]。文献[15-17]提出静止无人机基站(Base Station, BS)使用占空比(Duty Cycle, DC)方法为蜂窝用户提供服务；文献[18,19]提出基于先听后说(Listen Before Talk, LBT)机制的无人机2维飞行方案。文献[14]提出使用授权和免授权频谱为无人机用户提供服务，并优化子信道分配和功率实现了无人机用户上行和速率最大化。

由以上讨论可知，关于DAUAV研究存在以下问题：(1)只考虑了授权频谱共享；(2)鲜少研究蜂窝用户的通信性能^[3,6,13]；(3)忽略了无人机飞行能耗；(4)鲜少考虑无人机3维飞行轨迹。在无人机免授权通信系统中，大多数研究考虑静止无人机BS场景下基于DC和LBT方案的资源优化。虽然文献[14]考虑了全频谱共享方案，但是没有考虑到D2D用户，同时也忽略了无人机飞行轨迹和飞行能耗对频谱利用效率的影响。

因此，本文首先提出一种全频谱共享方法，即无人机能够通过控制上行蜂窝通信用户(Uplink Cellular Communication Users, UCCUs)和D2D用户的发射功率，在不影响WiFi设备正常传输的前提下使用免授权频谱，同时无人机也能够在不影响其他下行蜂窝通信用户(Downlink Cellular Communication Users, DCCUs)的前提下使用授权频谱。然后基于提出的全频谱共享方法，在无人机电池能量约束下联合优化无人机3维轨迹和发射功率，获得了蜂窝用户总吞吐量最大化。本文主要贡献总结如下：

(1)首次提出一种基于干扰控制的全频谱共享方法：在联合优化3维轨迹与UCCUs,D2D用户和无人机BS发射功率的前提下，UCCUs,D2D用户和多个WiFi用户能同时使用免授权信道，以及无人机能够重复使用授权信道；(2)采用分解的方法解决多变量耦合的非凸优化问题。本文将原问题分解为独立的两个子问题，即无人机3维轨迹 $\left\{ {{\boldsymbol{Q}}_{\text{u}}}[n], H[n] \right\}$ 联合优化子问题和UCCUs,D2D用户和无人机BS的发射功率 $\left\{ {{p_k}[n],{p_v}[n],\;{p_{\text{u}}}[n]} \right\}$ 联合优化子问题，然后采用连续凸逼近方法将两个非凸优化子问题进行非凸转凸。

2. 系统模型

2.1 场景描述

如图1所示，旋翼无人机携带BS调用全频谱为城市高密度区域的用户提供通信服务，其在一个任务周期T内从初始位置进行3维飞行到结束位置。本文将无人机的任务周期T划分为N个相等的时隙，无人机覆盖区域内存在K个UCCUs，随机地分布在WiFi 接入点(WiFi AP, WAP)的覆盖区域以外。UCCUs使用与WiFi设备相同的免授权信道，并采用频分多址(Frequency Division Multiple Access, FDMA)技术与无人机BS进行通信。同时，还存在V对D2D用户随机分布在无人机的覆盖区域内，每对D2D用户中的D2D发射端(D2D Transmitter, DT)与D2D接收端(D2D Receiver, DR)之间能够通过复用WiFi设备的免授权频谱直接进行数据通信，并在无人机BS控制下进行连接及资源分配。为了尽可能地提高免授权频谱效率，利用D2D用户传输距离短和发射功率低的特点，每对D2D用户都能够复用所有UCCUs的免授权频谱。此外，无人机覆盖区域内存在S个互不重叠的WAP，每个WAP为L个WiFi用户提供服务。在本文中只考虑使用相同频段的WAP，每个WAP之间不存在干扰。在无人机覆盖区域内，无人机能够复用P个DCCUs的频谱资源为区域内一个高速率需求的下行蜂窝用户(Downlink cellular User with High Rate demands, DUHR)提供通信服务。因此需要控制无人机BS发射功率使其对区域内所有DCCUs的干扰在阈值以下。无人机在飞行时一直覆盖着区域内所有用户，并不间断提供服务。关键符号变量总结如表1所示。

图 1 基于全频谱共享的DAUAV系统场景图

下载: 全尺寸图片幻灯片

表 1 关键的符号变量

符号变量	描述	符号变量	描述	符号变量	描述
T	飞行周期	${E_{\max }}$	无人机最大能量	${{\boldsymbol{L}}_{ { \rm{D} }_v^{\text{T} } } }$	DT的位置
N	飞行时隙数	$I_v^{ {\text{WiFi-R} } }$	WiFi对DR总干扰	${{\boldsymbol{L}}_{ {\text{W} }_l^s} }$	WiFi设备的位置
K	UCCUs的个数	${{\boldsymbol{Q}}_v}[n]$	D2D的吞吐量需求	${V_{\text{h}}}[n]$	3维飞行速度
L	WiFi的用户数	${p_s}$	WiFi的发射功率	$\beta$	地面路径损耗因子
S	WAP的个数	${{\boldsymbol{L}}_p}$	DCCU p的位置	${p_k}[n]$	UCCU k的发射功率
P	DCCUs的个数	$\sigma _{ \rm{L} }^2$	总授权干扰功率	${Q_k}$	UCCU k的速率需求
V	D2D的对数	${T_{\text{r}}}$	最大干扰阈值	$I_v^{ {\text{D2D-T} } }[n]$	DT对BS的总干扰
$I_v^{ {\text{CU-R} } }[n]$	UCCU对DR干扰	${L_{ { { \rm{O} }_k} } }$	UCCUs的位置	$I_l^{{\text{WiFi}}}[n]$	WiFi对BS的干扰
${{\boldsymbol{L}}_{ {\text{D} }_v^{\text{R} } } }$	DR的位置	${p_v}[n]$	DT v的发射功率	${p_{\text{u}}}[n]$	无人机的发射功率
${{\boldsymbol{Q}}_{\text{u} } }[n]$	无人机的2维坐标	$I_v^{{\text{D2D}}}[n]$	其他DT对DR干扰	$H[n]$	无人机的飞行高度

下载: 导出CSV

| 显示表格

2.2 授权和免授权频谱共享方案

如图2所示，在整个周期T内，K个UCCUs通过FDMA技术与无人机进行通信，而V对D2D用户可复用免授权总带宽 ${B_{\text{U}}}$ 直接进行数据传输，第 $k$ 个UCCU和第v对D2D用户在第n个时隙中分配的带宽分别为 ${B_k}[n]$ 和 ${B_v}[n]$ ，即 ${B_k}[n]{\text{ = }}{{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K}$ , ${B_v}[n]{\text{ = }}{B_{\text{U}}}$ 。同时无人机在飞行过程中能够复用地面无线通信的P个DCCUs的总授权频谱 ${B_{\text{L}}}$ 为区域内一个DUHR提供服务，并控制无人机BS对所有DCCUs的干扰^[8,10-12]。

图 2 基于全频谱共享的DAUAV系统方案图

下载: 全尺寸图片幻灯片

2.3 无人机轨迹模型

在本文中，无人机每个飞行时隙长度 $\delta = {T \mathord{\left/ {\vphantom {T N}} \right. } N}$ 。在第n个时隙中，无人机3维坐标为 ${{\boldsymbol{Z}}_{\text{u}}}[n] = \left( {{{\boldsymbol{Q}}_{\text{u}}}[n],H[n]} \right)$ , ${{\boldsymbol{Q}}_{\text{u}}}[n] = (X[n],Y[n])$ 表示为无人机的2维位置坐标， $H[n]$ 表示无人机的高度。假设无人机的最大水平飞行速度为 ${V_{{{x}},\max }}$ ，最大垂直飞行速度为 ${V_{{{z}},\max }}$ 。为了保证地面蜂窝用户与无人机之间的视距通信链路，无人机的最小飞行高度被设置为 ${H_{\min }}$ ，最大飞行高度为 ${H_{\max }}$ ^[1,8,9,20]。因此无人机的3维轨迹应满足约束为

$\qquad \left\| {{{\boldsymbol{Q}}_{\text{u}}}[n + 1] - {{\boldsymbol{Q}}_{\text{u}}}[n]} \right\| \le {V_{{{x}},\max }}\delta ,\forall n$

(1)

$\qquad \left\| {H[n + 1] - H[n]} \right\| \le {V_{{{z}},\max }}\delta ,\forall n$

(2)

$\qquad {H_{\min }} \le H[n] \le {H_{\max }},\forall n$

(3)

无人机在 $T$ 内从起始点 ${{\boldsymbol{Q}}_{\text{I}}}$ 飞到结束点 ${{\boldsymbol{Q}}_{\text{F}}}$ 的过程中为D2D和蜂窝用户提供服务，因此应满足约束为

$H[0] = H[N],{{\boldsymbol{Q}}_{\text{u}}}[0] = {{\boldsymbol{Q}}_{\text{I}}},{{\boldsymbol{Q}}_{\text{u}}}[N] = {{\boldsymbol{Q}}_{\text{F}}}$

(4)

2.4 D2D用户和UCCUs的传输模型

假设每个WAP区域中均存在一个虚拟WiFi设备，以这个虚拟WiFi设备来代表所有WiFi设备对无人机BS, DR和UCCUs的干扰。因此，UCCUs, DT和WiFi设备到无人机的信道功率增益分别表示为

$\left. \begin{split} & {g_k}[n] = {{{\beta _0}} \mathord{\left/ {\vphantom {{{\beta _0}} {\left( {{{\left\| {{L_{{{\text{O}}_k}}} - {Q_{\text{u}}}[n]} \right\|}^2} + {H^2}[n]} \right)}}} \right. } {\left( {{{\left\| {{{\boldsymbol{L}}_{{{\text{O}}_k}}} - {{\boldsymbol{Q}}_{\text{u}}}[n]} \right\|}^2} + {H^2}[n]} \right)}}\\ & g_v^{\text{T}}[n] = {{{\beta _0}} \mathord{\left/ {\vphantom {{{\beta _0}} {\left( {{{\left\| {{{\boldsymbol{L}}_{{\text{D}}_v^{\text{T}}}} - {{\boldsymbol{Q}}_{\text{u}}}[n]} \right\|}^2} + {H^2}[n]} \right)}}} \right. } {\left( {{{\left\| {{{\boldsymbol{L}}_{{\text{D}}_v^{\text{T}}}} - {{\boldsymbol{Q}}_{\text{u}}}[n]} \right\|}^2} + {H^2}[n]} \right)}}\\ & g_l^s[n] = {{{\beta _0}} \mathord{\left/ {\vphantom {{{\beta _0}} {\left( {{{\left\| {{{\boldsymbol{Q}}_{\text{u}}}[n] - {{\boldsymbol{L}}_{{{\text{W}}^s}}}} \right\|}^2} + {H^2}[n]} \right)}}} \right. } {\left( {{{\left\| {{{\boldsymbol{Q}}_{\text{u}}}[n] - {{\boldsymbol{L}}_{{{\text{W}}^s}}}} \right\|}^2} + {H^2}[n]} \right)}} \end{split}\right\}$

(5)

其中， ${\beta _0}$ 为1 m距离处的信道功率增益， ${{\boldsymbol{L}}_{{{\text{W}}^s}}} = {1 \mathord{\left/ {\vphantom {1 {L\sum\nolimits_{l = 1}^L {{L_{{\text{W}}_l^s}}} }}} \right. } {L\displaystyle\sum\nolimits_{l = 1}^L {{{\boldsymbol{L}}_{{\text{W}}_l^s}}} }}$ 表示虚拟WiFi设备的位置坐标， ${{\boldsymbol{L}}_{{{\text{O}}_k}}}$ 和 ${{\boldsymbol{L}}_{{\text{D}}_v^{\text{T}}}}$ , ${{\boldsymbol{L}}_{{\text{W}}_l^s}}$ 分别表示为UCCUs, DT和第s个WAP区域中第l个WiFi设备的位置坐标。WiFi设备到DR, UCCUs到DR以及UCCUs到WiFi设备信道功率增益分别表示为

$\left. \begin{split} & g_v^{{\text{R}},s} = {{\left( {{\beta _0}\xi _v^{{\text{R}},s}} \right)} \mathord{\left/ {\vphantom {{\left( {{\beta _0}\xi _v^{{\text{R}},s}} \right)} {{{\left\| {{L_{{{\text{W}}^s}}} - {L_{{\text{D}}_v^{\text{R}}}}} \right\|}^\beta }}}} \right. } {{{\left\| {{{\boldsymbol{L}}_{{{\text{W}}^s}}} - {{\boldsymbol{L}}_{{\text{D}}_v^{\text{R}}}}} \right\|}^\beta }}}\\ & g_{k,v}^{\text{R}} = {{\left( {{\beta _0}\xi _{k,v}^{\text{R}}} \right)} \mathord{\left/ {\vphantom {{\left( {{\beta _0}\xi _{k,v}^{\text{R}}} \right)} {{{\left\| {{L_{{{\text{O}}_k}}} - {L_{{\text{D}}_v^{\text{R}}}}} \right\|}^\beta }}}} \right. } {{{\left\| {{{\boldsymbol{L}}_{{{\text{O}}_k}}} - {{\boldsymbol{L}}_{{\text{D}}_v^{\text{R}}}}} \right\|}^\beta }}}\\ & g_{k,l}^s = {{\left( {{\beta _0}\xi _{k,l}^s} \right)} \mathord{\left/ {\vphantom {{\left( {{\beta _0}\xi _{k,l}^s} \right)} {{{\left\| {{L_{{{\text{O}}_k}}} - {L_{{{\text{W}}^s}}}} \right\|}^\beta }}}} \right. } {{{\left\| {{{\boldsymbol{L}}_{{{\text{O}}_k}}} - {{\boldsymbol{L}}_{{{\text{W}}^s}}}} \right\|}^\beta }}} \end{split}\right\}$

(6)

其中， $\beta$ 表示地面路径损耗因子，随机变量 $\xi _v^{{\text{R}},s},\,\xi _{k,v}^{\text{R}}$ 和 $\xi _{k,l}^s \sim \mathbb{C}\mathbb{N}\left( {0,1} \right)$ 均表示小规模衰落，代表遵循零均值指数分布的瑞利衰落系数， ${{\boldsymbol{L}}_{{\text{D}}_v^{\text{R}}}}$ 为DR的坐标。DT到WiFi设备和DT到DR的信道功率增益分别表示为

$\left. \begin{split} & g_{v,l}^{{\text{T}},s} = {{\left( {{\beta _0}\xi _{v,l}^{{\text{T}},s}} \right)} \mathord{\left/ {\vphantom {{\left( {{\beta _0}\xi _{v,l}^{{\text{T}},s}} \right)} {{{\left\| {{L_{{\text{D}}_v^{\text{T}}}} - {L_{{{\text{W}}^s}}}} \right\|}^\beta }}}} \right. } {{{\left\| {{{\boldsymbol{L}}_{{\text{D}}_v^{\text{T}}}} - {{\boldsymbol{L}}_{{{\text{W}}^s}}}} \right\|}^\beta }}}\\ & g_v^{{\text{T-R}}} = {{\left( {{\beta _0}\xi _v^{{\text{T-R}}}} \right)} \mathord{\left/ {\vphantom {{\left( {{\beta _0}\xi _v^{{\text{T - R}}}} \right)} {{{\left\| {{L_{{\text{D}}_v^{\text{T}}}} - {L_{{\text{D}}_v^{\text{R}}}}} \right\|}^\beta }}}} \right. } {{{\left\| {{{\boldsymbol{L}}_{{\text{D}}_v^{\text{T}}}} - {{\boldsymbol{L}}_{{\text{D}}_v^{\text{R}}}}} \right\|}^\beta }}} \end{split}\right\}$

(7)

其中，随机变量 $\xi _{v,l}^{{\text{T}},s}$ 和 $\xi _v^{{\text{T-R}}}$ 表示小规模衰落。因此，第 $k$ 个UCCU在每个时隙的吞吐量为

$\begin{split} {R_k}[n] =& {{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K}\delta {\log _2}\left( 1 + {\left( {{p_k}[n]{g_k}[n]} \right)} \right.\\ & \left.\mathord{\left/ {\vphantom {{\left( {{p_k}[n]{g_k}[n]} \right)} {\left( {{N_0}{{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K} + I_l^{{\text{WiFi}}}[n] + I_v^{{\text{D2D-T}}}[n]} \right)}}} \right. } {\left( {{N_0}{{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K} + I_l^{{\text{WiFi}}}[n] + I_v^{{\text{D2D-T}}}[n]} \right)} \right) \end{split}$

(8)

其中， ${N_0}$ 表示噪声功率谱密度， ${p_k}[n]$ 表示在第n个时隙中第k个UCCU的发射功率， $I_l^{{\text{WiFi}}}[n] = \displaystyle\sum\nolimits_{s = 1}^S {{p_s}g_l^s[n]}$ 表示WiFi设备对无人机BS的总干扰， ${p_s}$ 表示第s个WAP区域中对无人机造成干扰的WiFi设备的发射功率， $I_v^{{\text{D2D-T}}}[n] = \displaystyle\sum\nolimits_{v = 1}^V {{p_v}[n]g_v^{\text{T}}[n]}$ 表示DT设备对无人机BS的总干扰， ${p_v}[n]$ 表示在第n个时隙中第v对D2D设备中的DT设备的发射功率。因此，第k个UCCU在飞行周期T中的最低平均吞吐量应满足

$\left( {{1 \mathord{\left/ {\vphantom {1 T}} \right. } T}} \right)\sum\limits_{n = 1}^N {{R_k}[n] \ge {Q_k},\forall k}$

(9)

其中， ${Q_k}$ 为UCCU在整个周期内的最低速率需求。因此第v对D2D用户在每个飞行时隙的总吞吐量为

$\begin{split} {R_v}[n] =& {B_{\text{U}}}\delta {\text{E}}\left[ {{\log }_2}\left( 1 + {\left( {{p_v}[n]g_v^{{\text{T-R}}}} \right)}\right.\right. \\ & \left.\left. \mathord{\left/ {\vphantom {{\left( {{p_v}[n]g_v^{{\text{T-R}}}} \right)} {\left( {{N_0}{B_{\text{U}}} + I_v^{{\text{WiFi - R}}} + I_v^{{\text{CU - R}}}[n] + I_v^{{\text{D2D}}}[n]} \right)}}} \right. }{\left( {{N_0}{B_{\text{U}}} + I_v^{{\text{WiFi-R}}} + I_v^{{\text{CU-R}}}[n] + I_v^{{\text{D2D}}}[n]} \right)} \right) \right] \end{split}$

(10)

其中， ${\text{E}}\left[ \cdot \right]$ 是对随机变量 $\xi _v^{{\text{R}},s}$ , $\xi _{k,v}^{\text{R}}$ 和 $\xi _v^{{\text{T-R}}}$ 做期望运算， $I_v^{{\text{WiFi-R}}} = \displaystyle\sum\nolimits_{s = 1}^S {{p_s}g_v^{{\text{R}},s}}$ 表示WiFi设备对第v个DR的总干扰， $I_v^{{\text{CU-R}}}[n] = \displaystyle\sum\nolimits_{k = 1}^K {{p_k}[n]g_{k,v}^{\text{R}}}$ 表示UCCUs对第v对D2D设备中DR的总干扰， $I_v^{{\text{D2D}}}[n] = \displaystyle\sum\nolimits_{j = 1,j \ne v}^V {{p_j}[n]{g_{j,v}}}$ 表示第v对D2D设备中DR受到其他对D2D用户中DT的干扰之和， ${g_{j,v}}$ 表示第j对D2D设备中的DT到第v对D2D设备中DR的信道增益。因此第v对D2D用户设备在飞行时隙的最低吞吐量满足约束

${R_v}[n] \ge {Q_v}[n],\forall v,n$

(11)

其中， ${Q_v}[n]$ 表示每对D2D用户在第n个时隙内的最低传输吞吐量需求。为了保证WiFi设备通信质量，需保证每个WiFi设备受到的平均干扰功率之和不超过阈值 ${T_{\text{r}}}$ ^[10-12]，即满足约束

$\begin{split} & \left( {{1 \mathord{\left/ {\vphantom {1 N}} \right. } N}} \right)\left( \sum\limits_{k = 1}^K {\sum\limits_{n = 1}^N {{p_k}[n]{\text{E}}\left[ {g_{k,l}^s} \right]} } \right.\\ & \quad \left.+ \sum\limits_{v = 1}^V {\sum\limits_{n = 1}^N {{p_v}[n]{\text{E}}\left[ {g_{v,l}^{{\text{T}},s}} \right]} } \right) \le {T_{\text{r}}},\forall s,l \end{split}$

(12)

同时，UCCU的最小和最大功率 $p_{\rm{K}}^{\min }$ 和 $p_{\rm{K}}^{\max }$ ，D2D设备的最小和最大功率 $p_{\rm{V}}^{\min }$ 和 $p_{\rm{V}}^{\max }$ 应满足约束

$p_{\rm{K}}^{\min } \le {p_k}[n] \le p_{\rm{K}}^{\max },p_{\rm{V}}^{\min } \le {p_v}[n] \le p_{\rm{V}}^{\max },\forall k,v,n$

(13)

2.5 授权信道中DUHR的传输模型

无人机在第n个时隙的发射功率为 ${p_{\text{u}}}[n]$ ，保证DCCU受到的干扰不超过阈值 ${T_{\text{r}}}$ ^[8,10-12]需满足约束

$\left( {{{{p_{\text{u}}}[n]{\beta _0}} \mathord{\left/ {\vphantom {{{p_{\text{u}}}[n]{\beta _0}} {\left( {{{\left\| {{Q_{\text{u}}}[n] - {L_p}} \right\|}^2} + {H^2}[n]} \right)}}} \right. } {\left( {{{\left\| {{{\boldsymbol{Q}}_{\text{u}}}[n] - {{\boldsymbol{L}}_p}} \right\|}^2} + {H^2}[n]} \right)}}} \right) \le {T_{\text{r}}},\forall p,n$

(14)

其中， ${{\boldsymbol{L}}_p} = \left( {{X_p},{Y_p}} \right)$ 表示第p个DCCU的位置坐标。无人机BS的发射功率需要满足约束

${p_{\text{u}}}[n] \le P_{{\text{uav}}}^{\max },\forall n$

(15)

其中 $P_{{\text{uav}}}^{\max }$ 为无人机BS的最大发射功率。若 $\sigma _{\text{L}}^2$ 为高斯白噪声和地面干扰总功率^[8,10,11,21]，DUHR的吞吐量为

$\begin{split} {R_{\text{c}}}[n] = & {B_{\text{L}}}\delta {\log _2}\left( 1 + {{p_{\text{u}}}[n]{\beta _0}} \right.\\ & \left.\mathord{\left/ {\vphantom {{{p_{\text{u}}}[n]{\beta _0}} {\left( {\left( {{{\left\| {{Q_{\text{u}}}[n] - {L_{\text{c}}}} \right\|}^2} + {H^2}[n]} \right)\sigma _{\text{L}}^2} \right)}}} \right. }{\left( {\left( {{{\left\| {{{\boldsymbol{Q}}_{\text{u}}}[n] - {{\boldsymbol{L}}_{\text{c}}}} \right\|}^2} + {H^2}[n]} \right)\sigma _{\text{L}}^2} \right)} \right),\forall n \end{split}$

(16)

2.6 无人机能量损耗模型

为了增强无人机续航能力，本文考虑到旋翼无人机的3维推进能耗^[22-25]，因此能够得到

$\begin{split} {E_{{\text{fly}}}} =& \sum\limits_{n = 1}^N \delta ( {P_0}\left( {1 + {{3V_{\text{h}}^2[n]} \mathord{\left/ {\vphantom {{3V_{\text{h}}^2[n]} {S_{{\text{tip}}}^2}}} \right. } {S_{{\text{tip}}}^2}}} \right) \\ & + {P_{{\text{ii}}}}{{\left( {\sqrt {1 + {{V_{\text{h}}^4[n]} \mathord{\left/ {\vphantom {{V_{\text{h}}^4[n]} {\left( {4v_0^2} \right)}}} \right. } {\left( {4v_0^2} \right)}}} - {{V_{\text{h}}^2[n]} \mathord{\left/ {\vphantom {{V_{\text{h}}^2[n]} {\left( {2v_0^2} \right)}}} \right. } {\left( {2v_0^2} \right)}}} \right)}^{1/2}} \\ & + 0.5{d_0}\rho XAV_{\text{h}}^3[n] + W{V_{{z}}}[n] )\\[-15pt] \end{split}$

(17)

其中， $V_{\text{h}}^2[n] = V_{{x}}^2[n] + V_{{z}}^2[n]$ , ${V_{{x}}}[n]$ 和 ${V_{{z}}}[n]$ 分别表示水平飞行和垂直飞行速度，即 ${V_{{z}}}[n]\delta = \left\| H[n + 1] - H[n] \right\|$ 和 ${V_{{x}}}[n]\delta = \left\| {{{\boldsymbol{Q}}_{\text{u}}}[n + 1] - {{\boldsymbol{Q}}_{\text{u}}}[n]} \right\|$ ，A表示无人机转子盘面积， $\rho$ 表示空气密度， ${S_{{\text{tip}}}}$ 表示转子叶片叶尖速度， ${v_0}$ 表示悬停时平均旋翼诱导速度， ${d_0}$ 为机身阻力比，X表示转子坚固性， ${P_0}$ 和 ${P_{{\text{ii}}}}$ 分别表示悬停状态叶片轮廓和感应功率，W表示无人机重量。若 ${E_{\max }}$ 和 $E_{\rm{comm}}=(T/N)\displaystyle\sum\nolimits_N^{n=1}\;\;{p_{\text{u}}}[n]$ 分别表示无人机最大能量和通信能耗，则应满足约束

${E_{{\text{fly}}}} + E_{{\rm{comm}}} \le {E_{\max }}$

(18)

2.7 优化问题阐述

本文基于提出的全频谱共享方法，考虑在无人机电池能量约束下，联合优化无人机水平轨迹 ${\boldsymbol{Q}} = \left\{ {{{\boldsymbol{Q}}_{\text{u}}}[n]} \right\}$ 以及垂直飞行轨迹 ${\boldsymbol{H}} = \left\{ {H[n]} \right\}$ ，无人机BS,D2D用户和UCCUs的发射功率 ${\boldsymbol{P}} = \left\{ {p_{\text{u}}}[n], {p_v}[n],{p_k}[n] \right\}$ ，实现了UCCUs和DUHR的总吞吐量最大化。因此，针对优化目标建立的优化问题P1能够被表述为

$\begin{array}{l}\text{P}1:\underset{{\boldsymbol{Q}},{\boldsymbol{H}},{\boldsymbol{P}}}{ \mathrm{max}}{\displaystyle \sum \limits_{k=1}^{K}{\displaystyle \sum \limits_{n=1}^{N}{R}_{k}[n]}}+{\displaystyle \sum\limits _{n=1}^{N}{R}_{\text{c}}[n]},\\ \text{s}\text{.t}\text{. }式\text{(1)}-式\text{(4)}、式\text{(9)}、式\text{(11)}-式\text{(15)}、式\text{(18)}\end{array}$

(19)

约束式(9)和式(11)分别表示UCCUs和D2D用户的最低平均通信速率需求。约束式(12)保证第s个WAP区域中每个WiFi设备受到的平均干扰功率不超过WiFi设备的预设阈值，约束式(14)表示DCCUs受到的干扰功率之和不超过预设阈值，约束式(18)表示无人机每个飞行周期的总能耗不能超过最大电池能量。

3. 模型求解

由于 ${R_k}[n]$ 中的优化变量 $\left\{ {{\boldsymbol{Q}}_{\text{u}}}[n],H[n],{p_k}[n], {p_v}[n] \right\}$ 相互耦合， ${R_{\text{c}}}[n]$ 中的优化变量 $\left\{ {{{\boldsymbol{Q}}_{\text{u}}}[n],H[n],{p_{\text{u}}}[n]} \right\}$ 也是相互耦合的，因此目标函数是一个多变量耦合的复杂非凹函数，同时约束式(9)、式(11)和式(14)均为非凸约束，而约束式(18)也是非凸的，因此原问题P1是一个联合优化 $\left\{ {{{\boldsymbol{Q}}_{\text{u}}}[n],H[n],{p_k}[n],{p_v}[n],{p_{\text{u}}}[n]} \right\}$ 的复杂非凸优化问题，无法直接用现有凸优化方法解决。因此，本文首先基于块坐标下降方法将优化问题分解为两个联合优化子问题^{[6,8-11,19,20]}：3维轨迹 ${{\boldsymbol{Z}}_{\text{u}}}[n] = \left\{ {{{\boldsymbol{Q}}_{\text{u}}}[n],H[n]} \right\}$ 优化问题P2、功率 $\left\{ {p_k}[n], {p_v}[n],{p_{\text{u}}}[n] \right\}$ 优化问题P3，两个子问题均为非凸问题，需采用连续凸逼近方法进行非凸转凸，然后每个子问题通过使用凸优化工具CVX求解。

3.1 无人机3维轨迹联合优化子问题

令3维轨迹 ${{\boldsymbol{Z}}_{\text{u}}}[n] = \left\{ {{{\boldsymbol{Q}}_{\text{u}}}[n],H[n]} \right\}$ ，本文可以将原问题P1重写为非凹优化子问题P2

$\begin{split} & \text{P2}:\text{ }\underset{{{\boldsymbol{Z}}}_{\text{u}}[n]}{\mathrm{max}}\text{ }{\displaystyle \sum \limits_{k=1}^{K}{\displaystyle \sum \limits_{n=1}^{N}{R}_{k}[n]}}+{\displaystyle \sum \limits_{n=1}^{N}{R}_{\text{c}}[n]},\\ & \qquad \text{s}\text{.t}\text{. }式\text{(1)} - 式\text{(4)}、式\text{(9)}、式\text{(14)}、式\text{(18)} \end{split}$

(20)

首先通过引入松弛变量 ${\tilde {\boldsymbol{R}}_k} = \left\{ {{R_{k,{\text{lb}}}}[n],\forall k,n} \right\}$ 和 ${\tilde {\boldsymbol{R}}_{\text{c}}} = \left\{ {{R_{{\text{c}},{\text{lb}}}}[n],\forall n} \right\}$ ，能得到约束

$\qquad {R_k}[n] \ge {R_{k,{\text{lb}}}}[n],\forall k,n$

(21)

$\qquad {R_{\text{c}}}[n] \ge {R_{{\text{c}},{\text{lb}}}}[n],\forall n$

(22)

$\qquad \left( {{1 \mathord{\left/ {\vphantom {1 T}} \right. } T}} \right)\sum\limits_{n = 1}^N {{R_{k,{\text{lb}}}}[n] \ge {Q_k},\forall k}$

(23)

在约束式(19)中， ${R_k}[n]$ 可转化为 ${R_k}[n] = {{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K\delta {\log _2} \left( {{N_0}{{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K} + {Y_k}[n]} \right) - {{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K}\delta {\log _2} \left( {{N_0}{{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K} + Z[n]} \right)$ ，而 ${Y_k}[n] = \displaystyle\sum\nolimits_{s = 1}^S {{p_s}g_l^s[n]} + {p_k}[n]{g_k}[n] + \displaystyle\sum\nolimits_{v = 1}^V {{p_v}[n]g_v^{\text{T}}[n]}$ , $Z[n] = \displaystyle\sum\nolimits_{s = 1}^S {{p_s}g_l^s[n]} + \displaystyle\sum\nolimits_{v = 1}^V {{p_v}[n]g_v^{\text{T}}[n]}$ ，接下来引入非负的松弛变量 $\tilde {\boldsymbol{B}} = \left\{ {{b^s}[n],\forall s,n} \right\}$ 和 $\tilde {\boldsymbol{C}} = \left\{ c_v^{\text{T}}[n], \forall v, \right\}$ 对 ${R_k}{{[n]}}$ 进行处理。因此， ${R_k}{{[n]}}$ 下界表达式能够表达为

$\begin{split} {\hat R_k}{{[n]}}=& \left( {{{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K}} \right)\delta \left( {t_k}[n] - {{\log }_2}( {N_0}{{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K} \right. \\ & + \sum\limits_{s = 1}^S {\left( {{{{\beta _0}{p_s}} \mathord{\left/ {\vphantom {{{\beta _0}{p_s}} {{b^s}[n]}}} \right. } {{b^s}[n]}}} \right)} \\ & \left.\left.+ \sum\limits_{v = 1}^V {\left( {{{{\beta _0}{p_v}[n]} \mathord{\left/ {\vphantom {{{\beta _0}{p_v}[n]} {c_v^{\text{T}}[n]}}} \right. } {c_v^{\text{T}}[n]}}} \right)} \right) \right) \end{split}$

(24)

其中， ${t_k}[n] = {\log _2}\left( {{N_0}{{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K} + {Y_k}[n]} \right)$ ，假设给定第i次迭代局部值 ${\boldsymbol{Q}}_{\text{u}}^i[n]$ 和 ${H^i}[n]$ ，则 ${t_k}[n]$ 下界 ${\hat t_k}[n]$ 可转化为

$\begin{split} {\hat t_k}[n] =& {\log _2}\left( {U[n]} \right) - \frac{{{\beta _0}{p_k}[n]\left( {{W_k}[n] - W_k^i[n]} \right)}}{{\ln 2\left( {U[n]} \right){{\left( {W_k^i[n]} \right)}^2}}} \\ & - \sum\limits_{s = 1}^S {\frac{{{\beta _0}{p_s}\left( {{I_s}[n] - I_s^i[n]} \right)}}{{\ln 2\left( {U[n]} \right){{\left( {I_s^i[n]} \right)}^2}}}} \\ & - \sum\limits_{v = 1}^V {\frac{{{\beta _0}{p_v}[n]\left( {{J_v}[n] - J_v^i[n]} \right)}}{{\ln 2\left( {U[n]} \right){{\left( {J_v^i[n]} \right)}^2}}}} \end{split}$

(25)

其中， $U[n] = Y_k^i[n] + {N_0}{{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K}$ , ${W_k}[n] = \left\| {\boldsymbol{{L}}_{{{\text{O}}_k}}} - {{\boldsymbol{Q}}_{\text{u}}}[n] \right\|^2 + H{[n]^2}$ , ${I_s}[n] = {\left\| {{{\boldsymbol{Q}}_{\text{u}}}[n] - {{\boldsymbol{L}}_{{{\text{W}}^s}}}} \right\|^2} + H{[n]^2}$ , ${J_v}[n] = {\left\| {{{\boldsymbol{L}}_{{\text{D}}_v^{\text{T}}}} - {{\boldsymbol{Q}}_{\text{u}}}[n]} \right\|^2} + H{[n]^2}$ 。 $Y_k^i[n],W_k^i[n],I_s^i[n], J_v^i[n]$ 分别是 ${Y_k}[n],{W_k}[n],{I_s}[n],{J_v}[n]$ 在第i次迭代的局部值，此时 ${{\boldsymbol{Q}}_{\text{u}}}[n] = {\boldsymbol{Q}}_{\text{u}}^i[n]$ 和 $H[n] = {H^i}[n]$ 。

经过上述的转换，最终约束式(21)和约束式(22)能够分别转换为约束式(26)和式(27)

$\begin{split} & {{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K}\delta \left( {{\hat t}_k}[n] - {{\log }_2}\left( {N_0}{{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K} + \sum\limits_{s = 1}^S {\left( {{{{\beta _0}{p_s}} \mathord{\left/ {\vphantom {{{\beta _0}{p_s}} {{b^s}[n]}}} \right. } {{b^s}[n]}}} \right)} \right.\right.\\ & \left.\left. + \sum\limits_{v = 1}^V {\left( {{{{\beta _0}{p_v}[n]} \mathord{\left/ {\vphantom {{{\beta _0}{p_v}[n]} {c_v^{\text{T}}[n]}}} \right. } {c_v^{\text{T}}[n]}}} \right)} \right) \right) \ge {R_{k,{\text{lb}}}}[n],\forall k,n \end{split}$

(26)

$\begin{split} & {B_{\text{L}}}\delta( {\log _2}\left( {1 + {{\left( {{\beta _0}{p_{\text{u}}}[n]} \right)} \mathord{\left/ {\vphantom {{\left( {{\beta _0}{p_{\text{u}}}[n]} \right)} {\left( {\sigma _{\text{L}}^2G[n]} \right)}}} \right. } {\left( {\sigma _{\text{L}}^2G[n]} \right)}}} \right) - \tilde R[n]\\ & \cdot\left( {{{\left\| {{{\boldsymbol{Q}}_{\text{u}}}[n] - {{\boldsymbol{L}}_{\text{c}}}} \right\|}^2} + {H^2}[n] - G[n]} \right)) \ge {R_{{\text{c}},{\text{lb}}}}[n],\forall n \end{split}$

(27)

其中， $\tilde R[n] = {{\beta _0}{p_{\text{u}}}[n]} \mathord{\left/ {\vphantom {{{\beta _0}{p_{\text{u}}}[n]} {\left( {\ln 2\left( {\sigma _{\text{L}}^2{{\left\| {{\boldsymbol{Q}}_{\text{u}}^i[n] - {{\boldsymbol{L}}_{\text{c}}}} \right\|}^2} + \sigma _{\text{L}}^2{H^i}{{[n]}^2} + {\beta _0}{p_{\text{u}}}[n]} \right)G[n]} \right)}}} \right. } \left(\ln 2\left( \sigma _{\text{L}}^2{{\left\| {{\boldsymbol{Q}}_{\text{u}}^i[n] - {{\boldsymbol{L}}_{\text{c}}}} \right\|}^2} + \sigma _{\text{L}}^2{H^i}{{[n]}^2} + {\beta _0}{p_{\text{u}}}[n] \right)G[n] \right)$ , $G[n] = {\left\| {{\boldsymbol{Q}}_{\text{u}}^i[n] - {{\boldsymbol{L}}_{\text{c}}}} \right\|^2} + {H^i}{[n]^2}$ 。

在约束式(14)中引入非负松弛变量 $\tilde {\boldsymbol{A}} = \left\{ {{a_p}[n],\forall p,n} \right\}$ ，因此约束式(14)能够被重写为约束式(28)和式(29)

${{{p_{\text{u}}}[n]{\beta _0}} \mathord{\left/ {\vphantom {{{p_{\text{u}}}[n]{\beta _0}} {{a_p}[n]}}} \right. } {{a_p}[n]}} \le {T_{\text{r}}},\forall p,n$

(28)

$\begin{split} & {\left\| {{\boldsymbol{Q}}_{\text{u}}^i[n] - {{\boldsymbol{L}}_p}} \right\|^2} + 2{\left( {{\boldsymbol{Q}}_{\text{u}}^i[n] - {{\boldsymbol{L}}_p}} \right)^{\text{T}}}\left( {{{\boldsymbol{Q}}_{\text{u}}}[n] - {\boldsymbol{Q}}_{\text{u}}^i[n]} \right) \\ & \quad + {H^i}[n]\left( {2H[n] - {H^i}[n]} \right) \ge {a_p}[n],\forall p,n \end{split}$

(29)

其中， ${( \cdot )^{\text{T}}}$ 代表矩阵转置。相似地，能够得到

$\begin{split} & {\left\| {{\boldsymbol{Q}}_{\text{u}}^i[n] - {{\boldsymbol{L}}_{{{\text{W}}^s}}}} \right\|^2} + 2{\left( {{\boldsymbol{Q}}_{\text{u}}^i[n] - {{\boldsymbol{L}}_{{{\text{W}}^s}}}} \right)^{\text{T}}}\\ & \qquad \cdot \left( {{{\boldsymbol{Q}}_{\text{u}}}[n] - {\boldsymbol{Q}}_{\text{u}}^i[n]} \right) + {H^i}[n]\left( {2H[n] - {H^i}[n]} \right) \\ & \qquad \ge {b^s}[n],\forall s,n, \\ & {\left\| {{\boldsymbol{Q}}_{\text{u}}^i[n] - {\boldsymbol{{L}}_{{\text{D}}_v^{\text{T}}}}} \right\|^2} + 2{\left( {{\boldsymbol{Q}}_{\text{u}}^i[n] - {{\boldsymbol{L}}_{{\text{D}}_v^{\text{T}}}}} \right)^{\text{T}}}\\ & \qquad\qquad\cdot\left( {{{\boldsymbol{Q}}_{\text{u}}}[n] - {\boldsymbol{Q}}_{\text{u}}^i[n]} \right) + {H^i}[n]\left( {2H[n] - {H^i}[n]} \right) \\ & \qquad\ge c_v^{\text{T}}[n],\forall v,n \\[-10pt] \end{split}$

(30)

引入 $\tilde {\boldsymbol{O}} = \left\{ {o[n],\forall n} \right\}$ , $o[n] \ge ( \sqrt {1 + {{{{\left( {V_{{x}}^2[n] + V_{{z}}^2[n]} \right)}^2}} \mathord{\left/ {\vphantom {{{{\left( {V_{\text{x}}^2[n] + V_{\text{z}}^2[n]} \right)}^2}} {\left( {4v_0^2} \right)}}} \right. } {\left( {4v_0^2} \right)}}} - {{\left( {V_{{x}}^2[n] + V_{{z}}^2[n]} \right)} \mathord{\left/ {\vphantom {{\left( {V_{\text{x}}^2[n] + V_{\text{z}}^2[n]} \right)} {\left( {2v_0^2} \right)}}} \right. } {\left( {2v_0^2} \right)}} )^{0.5}$ ，能够得到

$\left( {{1 \mathord{\left/ {\vphantom {1 {{o^2}[n]}}} \right. } {{o^2}[n]}}} \right) \le {o^2}[n] + {{\left( {V_{{x}}^2[n] + V_{{z}}^2[n]} \right)} \mathord{\left/ {\vphantom {{\left( {V_{\rm{x}}^2[n] + V_{\rm{z}}^2[n]} \right)} {v_0^2}}} \right. } {v_0^2}},\forall n$

(31)

其中，右边表达式关于 $o[n]$ , ${V_{{x}}}[n]$ 和 ${V_{{z}}}[n]$ 是联合凸的，假设 ${o_i}[n],V_{{x}}^i[n],V_{{z}}^i[n]$ 分别是 $o[n],{V_{{x}}}[n],{V_{{z}}}[n]$ 在第i次迭代的局部值，因此约束式(31)和式(18)可以分别转换为约束式(32)和式(33)

$\begin{split} \left( {{1 \mathord{\left/ {\vphantom {1 {{o^2}[n]}}} \right. } {{o^2}[n]}}} \right) \le & {o_i}[n]\left( {2o[n] - {o_i}[n]} \right) \\ &+ {{\left( {\left( {V_{{x}}^i[n]} \right)\left( {2{V_{{x}}}[n] - V_{{x}}^i[n]} \right)} \right)} \mathord{\left/ {\vphantom {{\left( {\left( {V_{\text{x}}^i[n]} \right)\left( {2{V_{\text{x}}}[n] - V_{\text{x}}^i[n]} \right)} \right)} {v_0^2}}} \right. } {v_0^2}} \\ & + {{\left( {\left( {V_{{z}}^i[n]} \right)\left( {2{V_{{z}}}[n] - V_{{z}}^i[n]} \right)} \right)} \mathord{\left/ {\vphantom {{\left( {\left( {V_{\text{z}}^i[n]} \right)\left( {2{V_{\text{z}}}[n] - V_{\text{z}}^i[n]} \right)} \right)} {v_0^2}}} \right. } {v_0^2}},\forall n \end{split}$

(32)

$\begin{split} & \sum\limits_{n = 1}^N \delta \left( {P_0}\left( {1 + {{3\left( {V_{{x}}^2[n] + V_{{z}}^2[n]} \right)} \mathord{\left/ {\vphantom {{3\left( {V_{\text{x}}^2[n] + V_{\text{z}}^2[n]} \right)} {S_{{\text{tip}}}^2}}} \right. } {S_{{\text{tip}}}^2}}} \right) + {P_{{\text{ii}}}}o[n]\right. \\ & \quad\left. + 0.5{d_0}\rho XA{{\left( {V_{{x}}^2[n] + V_{{z}}^2[n]} \right)}^{\frac{3}{2}}} + W{V_{{z}}}[n] \right)\\ & \quad+ {E_{{\text{comm}}}} \le {E_{{\text{max}}}}\\[-10pt] \end{split}$

(33)

因此，非凸子问题(P2)能够转换为以下标准凸问题P2.1，可直接使用标准凸优化工具CVX来解决

$\begin{split} & \text{P2}\text{.1}:\text{ }\underset{\tilde{{\boldsymbol{B}}},\tilde{{\boldsymbol{C}}},\tilde{{\boldsymbol{A}}},\tilde{{\boldsymbol{O}}},{\tilde{{\boldsymbol{R}}}}_{\text{c}},{\tilde{{\boldsymbol{R}}}}_{k},{{\boldsymbol{Z}}}_{\text{u}}[n]}{\mathrm{max}}{\displaystyle {\sum\limits }_{k=1}^{K}{\displaystyle {\sum\limits }_{n=1}^{N}{R}_{k,\text{lb}}[n]}}\\ & \qquad\;\; +{\displaystyle \sum\limits_{n=1}^{N}{R}_{\text{c},\text{lb}}[n]},\\ & \text{s}\text{.t}\text{. }式\text{(1)}-式\text{(4)}、式\text{(23)}、式\text{(26)}-式\text{(30)}、\\ & \qquad 式\text{(32)}、式\text{(33)}\\[-1pt] \end{split}$

(34)

3.2 发射功率联合优化子问题

对发射功率 $\left\{ {{p_{\text{u}}}[n],{p_k}[n],{p_v}[n]} \right\}$ 联合优化子问题P3进行求解，因此可以将问题P1改写为问题P3

$\begin{split} & \text{P}3:\underset{{\boldsymbol{P}}}{\text{ }\mathrm{max}}\text{ }{R}_{\text{total}}\left({p}_{\text{u}}[n],{p}_{v}[n],{p}_{k}[n]\right),\\ & \quad\text{s}\text{.t}\text{. }式\text{(9)}、式\text{(11)}-式\text{(15)}、式\text{(18)} \end{split}$

(35)

令 ${D_{{\text{con}}}} = \displaystyle\sum\nolimits_{s = 1}^S {{p_s}g_l^s[n]} + {p_k}[n]{g_k}[n]$ , ${E_{{\text{con}}}} = {N_0}{{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K} + \displaystyle\sum\nolimits_{s = 1}^S {{p_s}g_l^s[n]}$ ，而 $t_k^3[n] = {\log _2}( {E_{{\text{con}}}} + \displaystyle\sum\nolimits_{v = 1}^V {{p_v}[n]g_v^{\text{T}}[n]} )$ , $t_k^4[n] = {\log _2}( {N_0}{{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K}{\text{ + }}{D_{{\text{con}}}}+ \displaystyle\sum\nolimits_{v = 1}^V {{p_v}[n]g_v^{\text{T}}[n]} )$ ，因此能够得到不等式(36)

$\begin{split} {R_k}[n] &= {{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K}\delta t_k^4[n] - {{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K}\delta t_k^3[n] \\ & \ge {{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K}\delta t_k^4[n] - {{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K}\delta \bar t_k^3[n] \end{split}$

(36)

在不等式(36)中，假设 $\left\{ {{p_v}[n]} \right\}$ 第i次迭代局部值为 $\left\{ {{p_{v,i}}[n]} \right\}$ ，因此 $t_k^3[n]$ 的上界 $\bar t_k^3[n]$ 可表示为

$\begin{split} \bar t_k^3[n] \triangleq& {\log _2}\left( {{E_{{\text{con}}}} + F_n^i} \right)\\ & + \sum\limits_{v = 1}^V \left( {{\beta _0}\left( {{p_v}[n] - {p_{v,i}}[n]} \right)} \right.\\ & /\left.{\left( {{{\left\| {{{\boldsymbol{L}}_{{\text{D}}_v^{\text{T}}}} - {{\boldsymbol{Q}}_{\text{u}}}[n]} \right\|}^2} + {H^2}} \right)} \right) \\ & /{\left( {\ln 2\left( {{E_{{\text{con}}}} + F_n^i} \right)} \right)} \end{split}$

(37)

其中， $F_n^i = \displaystyle\sum\nolimits_{v = 1}^V {{{\left( {{\beta _0}{p_{v,i}}[n]} \right)} \Bigr/ {\left( {{{\left\| {{{\boldsymbol{L}}_{{\text{D}}_v^{\text{T}}}} - {{\boldsymbol{Q}}_{\text{u}}}[n]} \right\|}^2} + {H^2}} \right)}}}$ 。因此，约束式(9)能够被重写为

$\left( {{1 \mathord{\left/ {\vphantom {1 T}} \right. } T}} \right)\sum\limits_{n = 1}^N {\left( {{{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K}\delta t_k^4[n] - {{{B_{\text{U}}}} \mathord{\left/ {\vphantom {{{B_{\text{U}}}} K}} \right. } K}\delta \bar t_k^3[n]} \right) \ge {Q_k},\forall k}$

(38)

若 ${E_n} = {N_0}{B_{\text{U}}} + \displaystyle\sum\nolimits_{s = 1}^S \left( {{{{\beta _0}{p_s}} \Bigr/ {{{\left\| {{{\boldsymbol{L}}_{{{\text{W}}^s}}} - {{\boldsymbol{L}}_{{\text{D}}_v^{\text{R}}}}} \right\|}^\beta }}}} \right) + {{\left( {{\beta _0}{p_k}[n]} \right)} \Bigr/ {{{\left\| {{{\boldsymbol{L}}_{{{\text{O}}_k}}} - {{\boldsymbol{L}}_{{\text{D}}_v^{\text{R}}}}} \right\|}^\beta }}} + \displaystyle\sum\nolimits_{j = 1,j \ne v}^V {{p_j}[n]} {g_{j,v}}$ ，则 ${R_v}[n]$ 为^[9]

$\begin{split} {R_v}[n] \approx& {B_{\text{U}}}\delta {\log _2}\left( {{E_n} + {p_v}[n]{{{\beta _0}} \mathord{\left/ {\vphantom {{{\beta _0}} {{{\left\| {{L_{{\text{D}}_v^{\text{T}}}} - {L_{{\text{D}}_v^{\text{R}}}}} \right\|}^\beta }}}} \right. } {{{\left\| {{{\boldsymbol{L}}_{{\text{D}}_v^{\text{T}}}} - {{\boldsymbol{L}}_{{\text{D}}_v^{\text{R}}}}} \right\|}^\beta }}}} \right) \\ & - {B_{\text{U}}}\delta {\log _2}\left( {{E_n}} \right) = {B_{\text{U}}}\delta {R_{v,1}}[n] - {B_{\text{U}}}\delta {R_{v,2}}[n] \end{split}$

(39)

若 $\left\{ {{p_j}[n]} \right\}$ 第 $i$ 次迭代局部值为 $\left\{ {p_{j,i}}[n], j = 1, 2, \cdots, V,j \ne v \right\}$ ，则 ${R_{v,2}}[n]$ 上界被表示为

$\begin{split} {\bar R_{v,2}}[n] \triangleq & {\log _2}\left( {{N_0}{B_{\text{U}}} + {E_{n,i}}} \right) + \Biggr( {E_n} - {E_{n,i}} - {N_0}{B_{\text{U}}} ) \\ & / {\left( {\ln 2\left( {{N_0}{B_{\text{U}}} + {E_{n,i}}} \right)} \right)}\\[-10pt] \end{split}$

(40)

其中， ${E_{n,i}} = \displaystyle\sum\nolimits_{s = 1}^S {\left( {{{{\beta _0}{p_s}} \mathord{\left/ {\vphantom {{{\beta _0}{p_s}} {{{\left\| {{L_{{{\text{W}}^s}}} - {L_{{\text{D}}_v^{\text{R}}}}} \right\|}^\beta }}}} \right. } {{{\left\| {{{\boldsymbol{L}}_{{{\text{W}}^s}}} - {{\boldsymbol{L}}_{{\text{D}}_v^{\text{R}}}}} \right\|}^\beta }}}} \right)} + \displaystyle\sum\nolimits_{j = 1,j \ne v}^V {{p_{j,i}}[n]} {g_{j,v}} + {{\left( {{\beta _0}{p_{k,i}}[n]} \right)} \mathord{\left/ {\vphantom {{\left( {{\beta _0}{p_{k,i}}[n]} \right)} {{{\left\| {{L_{{{\text{O}}_k}}} - {L_{{\text{D}}_v^{\text{R}}}}} \right\|}^\beta }}}} \right. } {{{\left\| {{{\boldsymbol{L}}_{{{\text{O}}_k}}} - {{\boldsymbol{L}}_{{\text{D}}_v^{\text{R}}}}} \right\|}^\beta }}}$ 。即约束式(11)可重写为

${B_{\text{U}}}\delta {R_{v,1}}[n] - {B_{\text{U}}}\delta {\bar R_{v,2}}[n] \ge {Q_v}[n],\forall v,n$

(41)

因此，优化问题P3能够转换为以下标准凸优化问题P3.1，可直接使用标准凸优化工具CVX来解决

$\begin{split} & \text{P3}\text{.1}:\text{ }\underset{{\boldsymbol{P}}}{\mathrm{max}}\text{ }\left({B}_{\text{U}}/K\right)\delta {\displaystyle \sum\limits _{k=1}^{K}{\displaystyle {\sum \limits}_{n=1}^{N}\left({t}_{k}^{4}[n]-{\overline{t}}_{k}^{3}[n]\right)}}\\ & \qquad\qquad\quad +{\displaystyle {\sum\limits }_{n=1}^{N}{R}_{\text{c}}[n]},\\ & \qquad\text{s}\text{.t}\text{.}\; 式\text{(12)}-式\text{(15)}、式\text{(18)}、式\text{(38)}、式\text{(41)} \end{split}$

(42)

3.3 算法收敛性与复杂度分析

接下来讨论整体算法的收敛性问题。首先，令 ${\boldsymbol{P}} = \left\{ {{p_{\text{u}}}[n],{p_k}[n],{p_v}[n],\forall k,v,n} \right\}$ , ${\boldsymbol{Z}} = {{\boldsymbol{Z}}_{\text{u}}}[n] = \left\{ {{\boldsymbol{Q}}_{\text{u}}}[n], H[n],\forall n \right\}$ ，接下来定义 $L({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^i})$ 是第i轮迭代得到的蜂窝用户的总吞吐量， ${{L}}_{\rm{p}}^{{\text{lb}},i}({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^i})$ 表示问题P2.1关于 $({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^i})$ 的目标函数值，即为第i轮迭代变量 ${{\boldsymbol{P}}^i}$ 近似得到的总吞吐量下界，而 $L_{{z}}^{{\text{lb}},i}({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^i})$ 表示问题P3.1关于 $({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^i})$ 的目标函数值，即为第i轮迭代变量 ${{\boldsymbol{Z}}^i}$ 近似得到的总吞吐量下界。在算法1的步骤4中，可得不等式

算法1 DAUAV系统中的3维轨迹优化和功率控制算法
(1) 初始化最大误差 $\varepsilon$ 、迭代次数初值i=0、目标函数初值obj、最大迭代次数 $\alpha$ 以及 ${{\boldsymbol{Z}}^i}\left( n \right) = \left\{ {p_{k,n}^i,p_{v,n}^i,p_{ {\text{u} },n}^i,{ {\tilde {\boldsymbol{O}}}^i},{\boldsymbol{Q}}_n^i,H_n^i} \right\}$ ;
(2) ${\bf{while} }$ $i < \alpha$ ${\bf{do} }$
(3) 　 $i = i + 1,$
(4) 　对于给定的 ${{\boldsymbol{Z}}^i}\left( n \right)$ 求解3维轨迹优化子问题P2.1,得到当前解为 $\left\{ { {\boldsymbol{Q} }_n^{i + 1},H_n^{i + 1},{ {\tilde {\boldsymbol{A} } }^{i + 1} },{ {\tilde {\boldsymbol{B} } }^{i + 1} },{ {\tilde {\boldsymbol{C} } }^{i + 1} },{ {\tilde {\boldsymbol{O} } }^{i + 1} },\tilde {\boldsymbol{R} }_c^{i + 1},\tilde {\boldsymbol{R} }_k^{i + 1} } \right\}$ ;
(5) 　对于给定的 $\left\{ {{\boldsymbol{Q}}_n^{i + 1},H_n^{i + 1},p_{k,n}^i,p_{v,n}^i} \right\}$ 求解功率优化子问题P3.1,得到当前问题解 $\left\{ {p_{k,n}^{i + 1},p_{v,n}^{i + 1},p_{{\text{u}},n}^{i + 1}} \right\}$ ;
(6) 　执行 $\left\{ { {\boldsymbol{Q} }_n^i,H_n^i,p_{k,n}^i,p_{v,n}^i,p_{ {\text{u} },n}^i,{ {\tilde {\boldsymbol{A} } }^i},{ {\tilde {\boldsymbol{B} } }^i},{ {\tilde {\boldsymbol{C} } }^i},{ {\tilde {\boldsymbol{O} } }^i},\tilde {\boldsymbol{R} }_{\text{c} }^i,\tilde {\boldsymbol{R} }_k^i} \right\} = \left\{ { {\boldsymbol{Q} }_n^{i + 1},H_n^{i + 1},p_{k,n}^{i + 1},p_{v,n}^{i + 1},p_{ {\text{u} },n}^{i + 1},{ {\tilde {\boldsymbol{A} } }^{i + 1} },{ {\tilde {\boldsymbol{B} } }^{i + 1} },{ {\tilde {\boldsymbol{C} } }^{i + 1} },{ {\tilde {\boldsymbol{O} } }^{i + 1} },\tilde {\boldsymbol{R} }_{\text{c} }^{i + 1},\tilde {\boldsymbol{R} }_k^{i + 1} } \right\};$
(7) 　 ${\bf{if } }{\text{ abs} }\left( {L\left( {p_{k,n}^i,p_{v,n}^i,p_{ {\text{u} },n}^i,{\boldsymbol{Q} }_n^i,H_n^i} \right) - {\rm{obj}}} \right) \le \varepsilon {\text{ } }{\bf{then} }$ break;
(8) 　else ${{\boldsymbol{Z}}^i}\left( n \right) = {{\boldsymbol{Z}}^{i + 1} }\left( n \right)$ 且 ${\rm{obj} } = L\left( {p_{k,n}^i,p_{v,n}^i,p_{ {\text{u} },n}^i,{\boldsymbol{Q}}_n^i,H_n^i} \right)$ ;
(9)end while
(10)输出最优参数值为 $\left\{ {{\boldsymbol{Q}}_n^i,H_n^i,p_{k,n}^i,p_{v,n}^i,p_{ {\text{u} },n}^i} \right\}$ ，计算获得当前最大的总吞吐量为 $L\left( {p_{k,n}^i,p_{v,n}^i,p_{ {\text{u} },n}^i,{\boldsymbol{Q}}_n^i,H_n^i} \right)$ 。

下载: 导出CSV

| 显示表格

$\begin{split} L({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^i}) & = L_{{z}}^{{\text{lb}},i}({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^i}) \le {\text{ }}L_{{z}}^{{\text{lb}},i}({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^{i + 1}}) \\ & \le L({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^{i + 1}}) \end{split}$

(43)

在式(43)中，因为1阶泰勒展开在给定局部点 ${Z^i}$ 下具有相同的目标值，即 $L({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^i}) = L_{{z}}^{{\text{lb}},i}({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^i})$ 成立。因为第i轮迭代中目标函数max(·)优化后最优值一定大于等于优化前的值，因此 $L_{{z}}^{{\text{lb}},i}({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^i}) \le {\text{ }}L_{{z}}^{{\text{lb}},i}({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^{i + 1}})$ 成立。而问题P2在给定 ${Z^{i + 1}}$ 时得到的下界等于问题P2.1的函数值，所以满足 ${\text{ }}L_{{z}}^{{\text{lb}},i}({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^{i + 1}}) \le L({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^{i + 1}})$ 。式(43)表明问题P2的目标值可以保证在每次迭代之后不会降低。相似地，针对算法1的步骤5，可得到

$\begin{split} L({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^{i + 1}}) & = L_{\rm{p}}^{{\text{lb}},i}({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^{i + 1}}) \le {\text{ }}L_{\rm{p}}^{{\text{lb}},i}({{\boldsymbol{P}}^{i + 1}},{{\boldsymbol{Z}}^{i + 1}}) \\ & \le L({{\boldsymbol{P}}^{i + 1}},{{\boldsymbol{Z}}^{i + 1}})\\[-10pt] \end{split}$

(44)

因此，最终可以得到不等式 $L({{\boldsymbol{P}}^i},{{\boldsymbol{Z}}^i}) \le L({{\boldsymbol{P}}^{i + 1}},{{\boldsymbol{Z}}^{i + 1}})$ 。综上所述，原问题P1的目标值在算法1每次迭代之后都不会下降，目标函数在迭代过程中不断趋向于一个最大值，最终会实现收敛。因此可以证明算法1是收敛的。具体流程如算法1所示。在步骤4中，令R=K+V+S，问题P2.1涉及对数形式，在CVX中通过采用内点法解决凸问题P2.1的复杂度为 $O(\mu ({(R + K + 5)^3}{N^3})\cdot \log_2 ({\varepsilon ^{ - 1}}))$ ^[9,22]，其中 $RN + KN + 5N$ 表示问题P2.1中的变量总数， $\mu$ 表示算法1的迭代次数，迭代过程中的最大容忍误差为 $\varepsilon$ 。相似地，在步骤5中，通过采用内点法在CVX中解决凸问题P3.1的复杂度为 $O(\mu ({(V + K + 1)^3}{N^3})\log_2 ({\varepsilon ^{ - 1}}))$ ^[9,22]， $VN + KN + N$ 表示问题P3.1中的变量总数，因此的总复杂度为 $O(\mu ({(R + K + 5)^3}{N^3} + (V + K + 1)^3\cdot {N^3})\log_2 ({\varepsilon ^{ - 1}}))$ 。

4. 仿真部分

4.1 仿真参数设置

本节进行大量仿真来验证提出算法的有效性和可行性。假设D2D用户随机分布在(X,Y)=500 m $\times$ 500 m的区域内，UCCU与DR设备以及WiFi设备之间的最短距离均为21 m，每对D2D用户中DT与DR之间距离均小于30 m^[3,7,26,27]，而DT与其他对D2D用户中的DR以及WiFi设备的最短距离均为15 m，最小高度为90/100 m^[1]，最大高度为180/200 m^[1]，无人机初始和终点位置坐标被分别设置为(–250, –350,H_min)m和(250,150,H_min)m^[8,20]。假设 $\beta$ 为3.76^[28], ${Q_v}[n] = {B_{\text{U}}}{R_{\min }}$ ，其中R_min=1.5～1.9 bit/(slot·Hz)，bit/(slot·Hz)表示每时隙1 Hz带宽传输的比特数。 $\sigma _{\text{L}}^2$ 为–80 dBm^[8,10,11,21]，T_r被设置为–77～–70 dBm^[8,14,28]。 ${N_0}$ 为–169 dBm/Hz^[3,6]， ${\beta _0}$ 为–50 dB， ${B_{\text{U}}}$ 为20 MHz^[14-19]， $p_{\rm{K}}^{\min }$ 和 $p_{\rm{K}}^{\max }$ 分别为0 dBm和27 dBm^[14]， $p_{\rm{V}}^{\min }$ 和 $p_{\rm{V}}^{\max }$ 分别为0 dBm和20 dBm^[3,6]， ${p_s}$ 为17 dBm^[28]，无人机BS的最大功率为27 dBm^[8]。部分参数设置如表2所示。

表 2 部分仿真参数列表

参数	取值	参数	取值
DCCUs的个数P	6/8/10/14	UCCUs的个数K	3～6
WAP的个数 $S$	2～4	飞行时隙数N	30
D2D用户对数V	6～10对	授权频谱带宽 ${B_{\text{L}}}$	6 MHz
WAP的WiFi用户L	10	最大水平速度	18 m/s
无人机飞行周期T	43～49 s	最大垂直速度	8 m/s
无人机最大能量 ${E_{\max }}$	13 kJ	UCCUs平均速率 ${Q_k}$	0.8 bit/(s·Hz)

下载: 导出CSV

| 显示表格

4.2 全频谱共享方案的结果分析

为了证明全频谱共享方案的有效性，将设置3种方案与提出方案进行对比：(1)基于全频谱共享的3维轨迹和功率优化(Three dimensional Trajectory optimization and Power optimization based on Full Spectrum Sharing, TTP-FSS)方案：基于全频谱共享方案，联合优化3维轨迹，UCCUs, D2D和无人机的发射功率，实现蜂窝用户的总吞吐量最大化；(2)基于全频谱共享的固定无人机最低高度(Fixed the Lowest Altitude based on Full Spectrum Sharing, FLA-FSS)方案：基于全频谱共享方案，联合优化无人机2维轨迹，UCCUs和D2D以及DUHR的功率，固定高度为最低飞行高度，实现总吞吐量最大化；(3)基于固定无人机直线轨迹的全频谱共享(Fixed Straight Trajectory based on Full Spectrum Sharing, FST-FSS)方案：基于全频谱共享方案，固定飞行轨迹为直线轨迹，联合优化UCCUs, D2D和DUHR的功率，实现总吞吐量最大化；(4)基于固定可行功率的全频谱共享(Fixed Feasible Power based on Full Spectrum Sharing, FFP-FSS)方案：基于全频谱共享方案，固定满足要求的功率集合，优化3维轨迹，实现总吞吐量最大化。

4.2.1 收敛性分析

显示了当K为3， $\beta$ 为3.76，T为48 s， ${E_{\max }}$ 为13 kJ， ${H_{\min }}$ 为100 m以及R_min为1.9 bit/(slot·Hz)时 ${T_{\text{r}}}$ 对提出方案与FLA-FSS方案的收敛性影响。显示了当 ${Q_k}$ 为0.4 bit/(s·Hz)， $\beta$ 为3.76，T为47 s， ${H_{\min }}$ 为90 m以及R_min为1.5 bit/(slot·Hz)时不同的K, V, P对提出方案与FLA-FSS方案的收敛性影响。从图3和图4中可以看出提出的TTP-FSS方案均能在20～35次达到收敛，在相同情况下提出方案的性能也始终优于FLA-FSS方案。这是因为在TTP-FSS方案中，无人机能够更自由调整高度来控制同信道干扰，特别是当无人机在飞行时更靠近WiFi用户和D2D用户时，无人机能够通过增大飞行高度来减小D2D用户和WiFi设备对自身的干扰以及无人机对DCCUs的干扰。从图3和图4中还能够看出当干扰阈值T_r降低时，提出方案得到的总吞吐量明显降低。当K为4， $\beta$ 为3.76， ${T_{\text{r}}}$ 为–73 dBm， ${E_{{\text{max}}}}$ 为13 kJ以及 ${H_{\min }}$ 为90 m时，不同T和R_min下提出方案的收敛性对比和不同的S,V,P(T=47 s，R_min=1.5 bit/(slot·Hz))下提出方案的收敛性对比分别如图5和图6所示。从图5中可以看出伴随着T的增加或者R_min的降低，蜂窝用户的总吞吐量会逐渐增大。在图6中，S,V,P的数量越多，蜂窝用户的总吞吐量就越低。这是因为S和V的增多会导致无人机受到更多干扰，而P越多，需要降低无人机的发射功率来减小对DCCUs的干扰。

图 3 干扰阈值对收敛性影响

下载: 全尺寸图片幻灯片

图 4 不同K, V, P下不同方案收敛性

下载: 全尺寸图片幻灯片

图 5 不同参数下提出方案收敛性

下载: 全尺寸图片幻灯片

图 6 不同S, V, P下提出方案收敛性

下载: 全尺寸图片幻灯片

4.2.2 不同优化方案的对比

显示出当K为4， $\beta$ 为3.76，T为47 s， ${H_{\min }}$ 为90 m以及 ${R_{\min }}$ 为1.5 bit/(slot·Hz)时D2D对数对不同方案性能的影响。从中可以看出当D2D对数增加时，所有方案总吞吐量逐渐降低。因为D2D对数越多，无人机会受到更多来自D2D用户的干扰。显示出当K为3， $\beta$ 为3.76，T为48 s， ${H_{\min }}$ 为100 m以及 ${R_{\min }}$ 为1.9 bit/(slot·Hz)时 ${T_{\text{r}}}$ 对不同方案性能的影响。在中，当 ${T_{\text{r}}}$ 增加时，不同方案总吞吐量逐渐增加。因为当T_r增大时，UCCUs和无人机的发射功率增大。在4种方案中，TTP-FSS和FLA-FSS方案明显优于另外两种方案，因为它们同时优化了飞行轨迹和用户传输功率。而TTP-FSS方案联合优化了3维轨迹和用户传输功率，考虑到高度变化对吞吐量性能的增益，因此提出方案性能会高于2维优化方案。

图 7 D2D的对数对总吞吐量的影响

下载: 全尺寸图片幻灯片

图 8 T_r对总吞吐量的影响

下载: 全尺寸图片幻灯片

图9、和分别显示出T, ${R_{\min }}$ 和 ${E_{{\text{max}}}}$ 对总吞吐量的影响以及不同方案下的吞吐量性能对比。中K为4， $\beta$ 为3.76， ${T_{\text{r}}}$ 为–73 dBm， ${E_{{\text{max}}}}$ 为13 kJ， ${H_{\min }}$ 为90 m和 ${R_{\min }}$ 为1.8 bit/(slot·Hz)。在图9中，当T增加时，所有方案的吞吐量都得到提高，因为无人机能够有更长时间停留在待服务蜂窝用户附近。在中，K, $\beta$ , ${T_{\text{r}}}$ , ${E_{{\text{max}}}}$ 和 ${H_{\min }}$ 与相同，T为47 s。当 ${R_{\min }}$ 增大时，所有方案的总吞吐量逐渐降低。这是因为 ${R_{\min }}$ 增大，D2D用户的功率提高，此时无人机受到的干扰会增大，总吞吐量会降低。在中，K为3， $\beta$ 为3.76， ${T_{\text{r}}}$ 为–77 dBm， ${H_{\min }}$ 为100 m，T为48 s和 ${R_{\min }}$ 为1.6 bit/(slot·Hz)，当 ${E_{{\text{max}}}}$ 较低时， ${E_{{\text{max}}}}$ 限制了飞行高度的变化，提出方案的总吞吐量与FLA-FSS方案是一致的。但是随着 ${E_{{\text{max}}}}$ 增加，提出方案性能始终更优于其他方案。FST-FSS方案的性能随 ${E_{{\text{max}}}}$ 的增大保持不变，因为无人机飞行能耗远高于通信能耗。在TTP-FSS, FLA-FSS和FFP-FSS方案中，总吞吐量伴随 ${E_{{\text{max}}}}$ 增大而增大，并最终保持恒定。因为 ${E_{{\text{max}}}}$ 限制了无人机飞行，当增大 ${E_{{\text{max}}}}$ 时，无人机能自适应调整高度靠近蜂窝用户提供服务。

图 9 T对总吞吐量的影响

下载: 全尺寸图片幻灯片

图 10 R_min对总吞吐量的影响

下载: 全尺寸图片幻灯片

图 11 E_max对总吞吐量的影响

下载: 全尺寸图片幻灯片

4.2.3 不同参数对无人机飞行轨迹的影响

图12和图13与图11中参数相同，对无人机水平轨迹和飞行高度的影响如和所示。从中可以看出当 ${E_{{\text{max}}}}$ 增大时，无人机将有更多时隙飞向并靠近UCCUs，此时无人机有更多电池能量能支持无人机飞行。在中， ${E_{{\text{max}}}}$ 越大，无人机最大飞行高度越大，因为无人机垂直飞行需要更多能量。

图 12 E_max对无人机水平轨迹的影响

下载: 全尺寸图片幻灯片

图 13 E_max对无人机飞行高度的影响

下载: 全尺寸图片幻灯片

、和显示 ${R_{\min }}$ 对无人机2维以及3维轨迹的影响。在TTP-FFS方案中，当无人机更靠近D2D用户、DCCUs和WiFi设备时，无人机会增加飞行高度，在靠近UCCUs或者DUHR时，无人机会降低飞行高度。当 ${R_{\min }}$ 增大时，无人机会靠近UCCU并提供上行通信服务。这是因为DT的发射功率伴随 ${R_{\min }}$ 的增大而增大，无人机受到的干扰也会增加，若要满足UCCU的最低速率，需要靠近每个UCCU。然而在中，对于同一 ${R_{\min }}$ 下的2维飞行轨迹而言无人机会更靠近UCCUs。在中，无人机开始飞行时最大飞行高度伴随 ${R_{\min }}$ 减小而增大，这是因为 ${R_{\min }}$ 减小时，DT功率降低，对无人机干扰减小，无人机不需要靠近UCCU去满足速率需求，此时无人机有更多时间靠近蜂窝用户密集区域。而 ${R_{\min }}$ 减小时，无人机在刚开始飞行时有更多时间靠近干扰区域，距离蜂窝用户较远，因此无人机会用更多时间飞离干扰区域。

图 14 R_min对2维飞行轨迹的影响

下载: 全尺寸图片幻灯片

图 15 R_min对水平飞行轨迹的影响(3维轨迹)

下载: 全尺寸图片幻灯片

图 16 R_min对飞行高度的影响(3维轨迹)

下载: 全尺寸图片幻灯片

当参数取值与图9的参数一致时，图17和图18显示T对3维轨迹以及2维轨迹的影响。在图17和图18中，当T增大时，无人机有更多时隙靠近UCCUs和DUHR。在FLA-FSS方案中，因为无法自适应调整无人机的高度，因此在飞离UCCUs和DUHR时会显著偏离干扰区域，同时需要更靠近UCCUs满足它们的通信需求。从图18中可以看出无人机在前十几个时隙会增加飞行高度，因为此时距离UCCUs以及DUHR较远，无人机受到的干扰以及对DCCUs的干扰随着无人机高度增加而减小。当周期T增加时，无人机在每个时隙中上升或下降的高度将会增加，而无人机靠近蜂窝用户时会降低高度提供更好的信道增益。

图 17 T对无人机水平与2维飞行轨迹的影响

下载: 全尺寸图片幻灯片

图 18 T对无人机飞行高度的影响

下载: 全尺寸图片幻灯片

5. 结论

本文首先针对DAUAV系统提出了一种全频谱共享方法，然后基于提出的共享方法，在无人机电池能量约束下通过联合优化3维轨迹和3类用户的发射功率从而最大化蜂窝用户的总吞吐量。为了求解系统产生的非凸优化问题，本文采用块坐标下降和连续凸逼近方法将原问题转化为3维轨迹优化和功率控制两个凸优化问题并迭代求解。仿真部分针对不同参数对蜂窝用户总吞吐量的影响进行了分析，验证了在基于频谱共享的 DAUAV系统中，无人机能够通过自适应调整高度来提高总吞吐量，并进一步验证了提出的共享方法显著提高了频谱利用效率，证明了本文算法的收敛性和有效性。

参考文献(1)

C. Silva. Shilnikov theorem-a tutorial, IEEE Trans. on Circuits Syst.-I, 1993, 40(10), 675 682.[2]Shui-Sheng Qiu.[J].A cell model of chaotic attractor, IEEE ISCAS97, Hong Kong.1997,:-[3]丘水生,混沌吸引子细胞模型研究的新结果,中国第十四届电路与系统学术年会论文集,福州,1998,37-41.[4]Shui-Sheng Qiu, Calculation of steady-state oscillations in non-linear circuits, Int. J. Electronics,1989, 67(3), 403-414.[5]Shui-Sheng Qiu, On verification of limit cycle stability in autonomous nonlinear systems, IEEE Trans. on Circuits Syst., 1988, 35(8), 1062-1064.[6]钟晓旭,三阶混沌系统的仿真与分析,[硕士论文],广州,华南理工大学,1998.

施引文献