物联网技术的飞速发展推动许多新兴应用的出现，如智能放牧、环境监测、自动驾驶和互动游戏等，但移动用户终端通常由于有限的资源难以满足新兴应用在通信和计算方面低时延的要求[1]．移动边缘计算(mobile edge computing，MEC)技术通过位于网络边缘的接入点提供较低传输延迟的密集计算服务，已经成为解决上述问题的一种有效的解决方案[2]．然而，传统的MEC在一些复杂的地形或恶劣天气中存在无法接入的情况，密集部署MEC服务器又是昂贵且不合实际的[3]．无人机(unmanned aerial vehicle，UAV)由于其突出的低成本性、机动性强和工作效率高等，被认为是安装MEC服务器为移动用户提供计算支持的重要手段，受到学术界和工业界的广泛关注[4]．无人机辅助的MEC系统可以在地面基础设施不足的偏远地区或抗险救灾等紧急情况下实现蜂窝网络覆盖，并通过利用良好的视距链路为地面用户提供计算服务[5-6]．然而，由于无人机和地面用户的尺寸及电池容量均严重受限，频繁进行充电或更换电池会带来巨大的人工和环境成本，能耗问题直接关系着无人机被调回充电之前地面用户和无人机所能完成的计算任务量，如何充分利用无人机的机动性等优势辅助MEC实现能耗最小化成为学术界研究的热点．文献[7-8]研究了地面用户的能耗最小化问题，但并未考虑无人机能耗的优化，而无人机能耗直接关系到无人机任务执行时间的长短．在文献[9-10]分别建立的固定翼和旋翼无人机推动力能耗模型的基础上，文献[11-13]研究了无人机能耗最小化问题，但并未联合考虑终端用户能耗的优化，而用户能耗往往关系到用户对提供服务的体验感．文献[14-16]研究了用户能耗和无人机总能耗最小问题．但上述研究没有考虑到无人机和地面用户的能耗大多不在同一水平的问题，文献[17]考虑了基于资源约束和轨迹约束下的MEC通信系统的加权总能耗最小化问题．文献[18]通过构建一个加权和函数，研究了无人机MEC系统中无人机和地面用户的总能耗最小化问题．上述研究虽然考虑了系统的加权能耗，却未研究权重因子会给无人机能耗和地面用户能耗带来怎样的影响，其主要关注的只是如何最小化加权能耗，没有关注权重因子对系统能耗影响．同时以上研究工作也没有考虑无人机作为中继对地通信时存在的回程链路功率约束问题，导致数据传输的时间和成本增加．此外，如果地面用户只有在无人机靠近时才卸载数据，地面用户可以降低任务卸载能耗，但此时，无人机却因更长的飞行轨迹，增大推动力能耗．因此，地面用户能耗的降低通常以增加无人机的飞行能耗为代价，这导致在无人机辅助MEC通信系统中存在一种新的空地能耗折中关系．因此不同于以上工作，这里研究了受限回程链路下无人机辅助MEC的多用户通信系统的空地能耗折中关系．且考虑到无人机能耗和地面能耗相差3~6个数量级，引入放大因子保证两者能耗的公平性．建立了回程链路约束下无人机辅助MEC的空地能耗折中框架．在回程链路速率约束和地面用户计算任务的约束下对卸载时间、计算任务分配、地面用户发射功率和无人机轨迹进行联合优化来最小化系统加权总能耗．此外，在给定放大因子条件下，通过改变权衡因子，验证空地能耗折中关系的同时兼顾两者能耗的公平性．为求解高度非凸的加权总能耗最小化的优化问题，提出一种基于次梯度算法和连续凸近似技术的两阶段交替迭代优化算法．1 系统模型及优化问题考虑一架单天线蜂窝连接的旋翼无人机为地面K个单天线地面用户(ground users，GUs)提供计算服务，且存在一个地面基站BS为无人机提供受限的回程链路以卸载从GUs接收的计算任务．GUs和无人机均配备机载通信电路和计算处理器，而地面基站具有超高性能的处理服务器，能够在很小的时延下提供较大速率的计算服务．为减少GUs之间的竞争以提高卸载效率，在该回程链路速率受限的通信网络中采用时分多址协议进行通信．不失一般性，采用三维欧几里得坐标系．假设BS和所有GUs的位置都固定在高度为零的地面上，则BS的水平位置表示为qa=(xa,ya)，GU k∈𝒦={1,2,⋯,K}的水平位置表示为qk=(xk,yk)，这里假设BS和GUs的位置均已知．无人机在任务完成时间T内保持固定高度(H0)飞行，T被等间隔分为N个足够小的时隙οn，即οn=T/N．令n∈𝒩={1,2,⋯,N}表示所有时隙的集合，则在时隙n无人机的水平位置坐标可表示为qu[n]=(xu[n],yu[n])，并满足以下移动性约束qu[n+1]-qu[n]≤οnVmax  (∀n);qu[N]=qF,qu[1]=qI;v[n]=qu[n+1]-qu[n]/οn, (1)式中：Vmax为无人机的最大飞行速度；qI和qF为无人机的初始和最终水平位置．1.1　信道模型与文献[8-10]类似，考虑无人机对地卸载的通信信道是视距链路，因此在时隙n，无人机到地面Gus、基站的信道功率增益分别为hui[n]=β0 /H2+qu[n]-qi2      (i∈{k,a}) ，(2)式中β0为参考距离d0=1 m处的信道功率增益．令puk[n]表示GU k的发射功率，pva[n]表示无人机的发射功率，则GUs到无人机之间的卸载速率ruk[n]和无人机到基站之间的卸载速率rua[n]分别为ruk[n]=Blog21+puk[n]huk[n]σ02;rua[n]=Blog21+pva[n]hua[n]σ12, (3)式中：B为通信带宽；σ02和σ12分别为GUs和基站处的接收噪声功率．1.2　能耗模型在建立的MEC场景中，忽略无人机发送给GUs的时延和能耗，系统的总能耗由卸载能耗、本地计算能耗和无人机飞行能耗3个部分组成．a．卸载能耗模型对于每个GU k在时隙n∈𝒩处被划分为两个子时隙，包括GU k卸载到无人机的时间和无人机卸载到基站的时间，分别由τk[n]和ya,k[n]表示，并满足约束条件：∑k=1K(τk[n]+ya,k[n])≤οn  (∀k,n), (4)则系统中无人机和GUs的卸载能耗可以分别写为Eu,koff=ya,k[n]pva[n],Eg,koff=τk[n]puk[n]．b．本地计算能耗模型为了有效提高计算效率，减小计算能耗，GUs和无人机采用动态电压频率缩放技术．令lg,k[n]和lu,k[n]分别表示GU k和无人机在时隙n的本地计算数据量，则GUs和无人机在每个时隙n内的计算能耗可以分别表示为[14]：Ekcom[n]=κk(cklg,k[n])3/οn2;Eu,kcom[n]=κu(culu,k[n])3/οn2,式中符号含义及取值见表1．10.13245/j.hust.240073.T001表1仿真参数设置参数名称符号及取值参数名称符号及取值信号带宽B=10 MHz无人机飞行高度H=20 m无人机CPU有效电容系数κu=1×10-27无人机计算单位比特CPU循环次数Cu=1 000 cycles/bit空气密度ρ=1.225 0 kg/m3GUs计算单位比特CPU循环次数Ck=1 000 cycles/bitGUs CPU有效电容系数κk=1×10-27无人机叶片扫过的面积A=0.503 0 m2无人机叶尖角速度Utip=120 m/s无人机悬停时平均旋翼诱导速度v0=4.03 m/s机身阻力比d0=0.301 0总叶片面积与叶片扫过面积之比s=0.049 9GUs任务量Lk=4 Mbit无人机悬停状态下的诱导功率Pi=88.63 W无人机最大CPU频率fumax=6 GHzGUs最大CPU频率fgmax=3 GHz参考信道增益β0=-50 dB无人机悬停状态下型阻功率P0=158.76 Wc．无人机飞行能耗模型旋翼无人机的推动力功耗可表示为[11]Pv[n]=P01+3v[n]2/Utip2+Pi1+v[n]4/4v04-v[n]2/(2v02)1/2+d0ρsAv[n]3/2,式中关于无人机悬停状态下的功率参数及飞行状态下的空气动力学有关参数的含义见表1，因此在时隙n无人机的飞行能耗简写为Efly[n]=οnPv[n]．综上所述，GUs和无人机的总能耗可以表示为：Eg,k[n]=Eg,koff[n]+Ekcom[n];Eu[n]=∑k=1KEu,koff[n]+∑k=1KEu,kcom[n]+Efly[n];Eg=∑k=1K∑n=1NEg,k[n];Eu=∑n=1`NEu[n].1.3　优化问题定义MEC系统中卸载时间Ya={τk[n],ya,k[n],k∈𝒦,n∈𝒩}，计算任务比特划分La={lg,k[n],lu,k[n],la,k[n],k∈𝒦，n∈𝒩，GUs发射功率Pa={puk[n],k∈𝒦，n∈𝒩，无人机轨迹Qa={qu[n],v[n],n∈𝒩}，则系统加权总能耗ET最小化的优化问题可以写为P1:minYa,La,Pa,Qa(1-θ)wEg+θEu,s.t.    (1),(4);         C1:∑i=1n(lu,k[i]+ya,k[i]rua[i])≤∑i=1n-1τk[i]ruk[i]  (∀k,n),C2:la,k(n)≤yk[n]rua[n]  (∀k,n),C3:lu,k[n]+lg,k[n]+la,k[n]≥Lkmin  (∀k,n),C4:ruk[n]≤rua[n]  (∀k,n),C5:culu,k[n]≤οnfumax,cklg,k[n]≤οnfgmax  (∀k,n),C6:0≤puk[n]≤Pmax  (∀k,n),C7:lu,k[n]≥0,lg,k[n]≥0,la,k[n]≥0  (∀k,n),式中：考虑到GUs与无人机能耗之间的公平性，引入w作为GUs总能耗的放大因子；0≤θ≤1为无人机能耗与GUs能耗的权衡因子；C1和C2是两个信息因果约束；C3是GUs的计算任务约束，以确保所有GUs的计算任务在规定时间内处理完毕；C4为回程链路速率约束；C5和C6对无人机和GUs的计算能力及GUs的发射功率进行限制；C7确保优化变量的非负性．显然，问题P1是一个复杂的非凸优化问题，既然非凸目标函数之间存在非线性耦合变量Ya，La，Pa和Qa．为了解决这个问题，提出一种两阶段交替优化算法，此算法可以保证收敛到局部最优解．2 优化问题求解为了求解问题P2，在块坐标下降法和连续凸近似技术的基础上，提出两阶段交替优化算法进行求解．2.1　优化计算任务分配La和地面用户发射功率Pa首先，对于任意给定的卸载时间Ya和无人机轨迹Qa，原问题P1可以重写为：P2:minLa,Pa  (1-θ)w∑n=1N∑k=1KEg,k[n]+θ∑n=1N∑k=1KEu,kcom[n]; s.t.       C1~C7.容易证明，问题P2是凸的并且符合Slater条件[15]，可以采用拉格朗日对偶方法求得计算任务分配和GUs发射功率的闭式最优解Laopt和Paopt分别为：lg,kopt[n]=οnμk,n+ckω¯k,n3w(1-θ)κkck31/2(0≤lg,kopt[n]≤fgmax/ck);lu,kopt[n]=οnρk,n+μk,n+cuλk,n3θκucu31/2(0≤lu,kopt[n]≤fumax/cu);pukoptn=Bαk,n-∑i=1nτk[i]ρk,iln2wθ-1τk[n]-σ02huk[n](0≤pukoptn≤Pmax),式中ρk,n≥0，μk,n≥0，αk,n≥0，λk,n≥0，ω¯k,n≥0，分别表示与约束条件相关的对偶变量．2.2　优化卸载时间Ya和无人机轨迹Qa给定第一阶段获得的GUs的发射功率最优解Paopt和计算任务分配最优解Laopt，则卸载时间Ya和无人机轨迹Qa的优化子问题可以写为：P3:minYa,Qa(1-θ)w∑n=1N∑k=1KEg,koff[n]+θ∑n=1N∑k=1KEu,koff[n]+Efly[n],s.t.        (1),(4);C8:∑i=1nlu,kopt[i]≤∑i=1nτk[i]φk,n[i]  (∀k,n),C9:la,k(n)≤Bya,k[n]∙log21+pva[n]γ1H2+qu[n]-qa2  (∀k,n),C10:φk,n[n]≤B∙log21+pva[n]γ1H2+qu[n]-qa2  (∀k,n),C11:∑i=1nlu,kopt[i]+ya,k[i]B∙log21+pva[n]γ1H2+qu[i]-qa2≤∑i=1n-1τk[i]φk,n[i]  (∀k,n).显然P3是难以求解的非凸优化问题，因此引入SCA技术将其转化为凸优化问题．具体而言，为了解决约束条件中辅助变量φk,n[n]的非凸性，在第j次迭代中，对于给定的局部点quj[n]，利用一阶泰勒展开获得quj[n]的下界αuk,j[n]-γlog2(e)pukopt[n]χuk,j[n]βuk,j[n](χuk,j[n]+pukopt[n]γ)，式中：αuk,j[n]=log21+pukopt[n]γH2+quj[n]-qk2；βuk,j[n]=qu[n]-qk2-quj[n]-qk2；χuk,j[n]=quj[n]-qk2+H2．令：αua,j[n]=log21+pva[n]γ1H2+quj[n]-qa2；βua,j[n]=qu[n]-qa2-quj[n]-qa2；χua,j[n]=quj[n]-qa2+H2．引入松弛变量ϑkn2，约束C8可以转化为以下凸约束：ϑkn2τk[n]≤Bαuk,j[n]-Bγlog2(e)pukopt[n]χuk,j[n]βuk,j[n](χuk,j[n]+pukopt[n]γ);(ϑknj)2+2ϑknj(ϑkn-ϑknj)≥lu,kopt[n].(5)同理，约束C9通过引入松弛变量ϕkn2可转化为以下凸约束形式：ϕkn2ya,k[n]≤Bαua,j[n]-Bγ1log2(e)pva[n]χua,j[n]βua,j[n](χua,j[n]+pva[n]γ1);(ϕknj)2+2ϕknj(ϕkn-ϕknj)≥la,k[n].(6)引入松弛变量υkn，约束C10可以分解为以下两个子约束，从而转化为凸约束：υkn≤αuk,j[n]-γlog2(e)pukopt[n]χuk,j[n]βuk,j[n]χuk,j[n]+pukopt[n]γ;αua,j[n]-γ1log2(e)pva[n]χua,j[n]βua,j[n]χua,j[n]+pva[n]γ1≥υkn. (7)同理，非凸约束C11的凸函数形式可以表示为∑i=1nlu,kopt[i]+ya,k[i]υkn≤∑i=1n-1τk[i]υkn  (∀k,n)．(8)引入辅助变量𝓁n2=1+v[n]4/4v041/2-v[n]2/2v02  (∀n)，应用一阶泰勒级数将目标函数转化为凸形式，即在第j次迭代中，给定𝓁nj和vj[n]，无人机能耗可以重写为Eapprofly[n]=οnP01+3v[n]2Utip2+Pi𝓁n+d0ρsAv[n]3/2.综上所述，问题P3的近似凸问题可写为：P4:minYa,Qa,𝓁n,ϑkn,ϕkn,υkn(1-θ)w∑n=1N∑k=1KEg,koff[n]+θ∑n=1NEuoff[n],                     s.t.     (1),(4),(5)~(8);C12:1𝓁n2≤𝓁n(j)2+2𝓁nj(𝓁n-𝓁nj)-vj[n]2v02+2v02∙(vj[n])Tv[n],C13:𝓁n≥0,ϑkn≥0,ϕkn≥0,υkn≥0  (∀k,n),其中Euoffn=∑k=1KEu,koff[n]+Eapprofly[n]，且约束C13保证松弛变量大于零．显然，问题P4是一个标准的凸优化问题，可使用CVX等工具求解[19].3 仿真结果及分析通过大量仿真结果验证所提出的多目标优化算法的有效性，为便于性能对比，基准方案w=1，w=100和w=1 000中，优化变量为Ya，La，Pa和Qa；“固定功率”方案中固定GUs发射功率Pa，此时优化变量为Ya，La和Qa．具体仿真设置参数如表1所示．图1显示了不同放大因子下无人机辅助MEC系统总能耗随加权因子θ的变化．可以看出：当w=1 000，θ=0.7时，系统的总能耗最小，因此以w=1 000作为后续仿真中放大因子设置的依据，验证了通过设置放大因子来考虑空地能耗公平性在提升系统性能方面的必要性．图2给出了GUs能耗与无人机能耗随加权因子θ变化的关系．不难看出：随着θ增加，所有优化方案下GUs的能耗增加，无人机的能耗减小．这是因为当θ增加时无人机能耗占总能耗的比例变大，此时优化无人机的能耗更为重要．当放大因子w增加，GUs的能耗减小，无人机能耗增加，此时优化GUs的能耗越来越重要．10.13245/j.hust.240073.F001图1系统总能耗随θ变化的关系10.13245/j.hust.240073.F002图2GUs无人机能耗随θ变化的关系图3所示为当w=1 000时无人机在不同加权因子θ下的最优飞行轨迹．可以看出随着θ的增加，无人机从整体上更靠近GUs进行数据卸载，并且在较长的时间内保持相对离BS距离更近的飞行轨迹．同时这样的轨迹保证了在回程链路中无人机对基站的卸载速率大于GUs的发射速率．结合图3和图4可知：旋翼无人机处于加速飞行和悬停状态时将会产生较大的能耗．10.13245/j.hust.240073.F003图3无人机轨迹随θ变化的关系10.13245/j.hust.240073.F004图4不同任务量下的计算任务分配图4所示为在GUs不同任务需求下的计算任务分配数据量．可以看出：当任务需求较小时，GUs选择的本地计算模式在减小能耗方面更为有效．当任务量增加，无人机选择将一部分计算数据卸载到BS上，在此过程中会增加无人机的卸载能耗，但是会降低GUs和无人机的计算能耗，从而达到降低MEC系统总能耗的目的．图5给出了w=1 000时无人机能耗和加权地面用户能耗的Pareto曲线．可以看出降低无人机的能耗是以增加GUs能耗为代价的，这个结论和前面的理论分析结果一致．10.13245/j.hust.240073.F005图5无人机能耗与GUs能耗的关系面向回程链路速率受限的无人机辅助的MEC系统，以最小化系统加权总能耗为目标，分别引入权重因子和放大因子，考虑系统的空地能耗折中问题的同时兼顾两者的公平性．基于块坐标下降法和连续凸优化技术提出了一种卸载时间、GUs发射功率、计算任务分配和无人机轨迹联合优化的两阶段交替迭代算法．仿真结果进一步展示了不同放大因子和加权因子对系统总能耗的影响和保证公平性条件下空地能耗的折中关系，以及所提出方案对比其他基准方案在降低总能耗方面具有显著的优越性．下一步继续考虑由于GUs的高移动性造成的动态信道环境下无人机辅助MEC的资源分配策略．