1.   引言

随着电子系统向信息化和智能化发展，新一代多功能一体化电子系统通过共用信号、信道、天线、数据处理等硬件和软件资源以实现多功能共用，成为国内外学术界和工业界的研究热点^[1-4]。其中，多输入多输出(Multi-Input Multi-Output, MIMO)技术由于具有空间分集与复用增益等优势，在探通一体化系统中可增强系统性能，从而引发了广泛关注^[5-8]。

MIMO探通一体化系统根据载频数量可分为两种：(1)多载频MIMO一体化系统：以跳频(Frequency Hopping, FH)和正交频分复用(Orthogonal Frequency Division Multiplexing, OFDM)信号为代表^[9]，通信信息可调制于频率、频率增量、相位、方向图等参数中；(2)单载频MIMO一体化系统：通过设计发射波形矩阵或正交波形加权矩阵^[10]，用户端通信信息可采用幅移键控(Amplitude Shift Keying, ASK)、相移键控(Phase Shift Keying, PSK)、正交调幅(Quadrature Amplitude Modulation, QAM)等调制。

对于多载频MIMO一体化系统，Hassanien等人^[11]提出将通信信息PSK调制在MIMO雷达正交FH信号中。文献[12]放松了不同FH信号的正交性以增加通信数据量和优化探通一体化性能，并在文献[13]中通过空域调制进一步提升数据传递速率。文献[14]提出了MIMO系统多载波PM-FMCW架构。通过选择阵列天线、调制相位、载波和天线发射信号置换选择调制通信信息，增加了通信信息量，且雷达分辨率接近于等孔径宽带雷达。文献[15,16]采用多子脉冲OFDM信号实现雷达通信一体化，并提出了基于通信信息补偿和解相干处理的目标距离和速度超分辨估计算法。对于单载频MIMO一体化系统，部分MIMO探通一体化系统利用阵列空间几何结构，通过多正交信号的线性加权将通信信息嵌入其旁瓣方向图幅度或(和)相位中。其中，文献[17]可同时实现主、旁瓣PSK下行链路通信，且有更高的精度和数据率，文献[18]进一步介绍了如何利用接收波束形成器将上行通信信号从接收混合信号中分离。文献[19]提出采用置换矩阵混排天线间发射信号顺序以调制通信信息，设计正交信号加权矩阵以实现理想雷达波束图，该调制方法对雷达方向图性能无影响。文献[20]针对分离和共享式天线阵列，提出了基于方向图模板匹配的一体化波束成型方法，保证了下行用户的信干噪比(Signal to Interference plus Noise Ratio, SINR)以及系统的发射总功率预算。文献[21]提出了基于空-频域信息调制的MIMO一体化波形设计方法，通过通信方向接收信号多频谱约束下最小化发射方向图积分旁瓣电平(Integrated Sidelobe Level, ISL)，保证雷达探测功能的同时嵌入通信信息。

本文针对基于预编码矩阵调制的MIMO探通一体化系统，提出了基于交替方向乘子(Alternation Direction Method of Multipliers, ADMM)的一体化信号矩阵设计方法。具体而言，在用户和窃听用户参考密码本约束下通过最大化MIMO一体化系统波束峰值主瓣旁瓣电平比(Peak Mainlobe to Sidelobe level Ratio, PMSR)，保证探测方向图性能的同时防止通信信息被窃听。针对预编码矩阵解调问题，提出了基于交替方向惩罚(Alternating Direction Penalty Method, ADPM)的排序学习优化解调算法，解决了穷举法导致的维数灾难问题，可快速有效恢复置换矩阵。

本文后续使用以下符号定义：粗体小写字母(例如：a)用于表示向量；粗体大写字母(例如：A)用于表示矩阵；${(\cdot)}^{\text{T}}$, ${(\cdot)}^{\text{*}}$和${(\cdot)}^{\text{H}}$分别表示转置、共轭和共轭转置。${\mathbb{C}^{N \times M}},{\mathbb{R}^N},{\mathbb{C}^N}$分别表示$N \times M$维复矩阵、N维实向量和复向量。$\left\| {} \right\|$表示向量2范数。${{\boldsymbol{I}}_N}$和${{{ {\textit{1}}}}_N}$分别表示$N \times N$维单位阵和N维全1向量。$\Re \{\cdot\}$和$\Im \{\cdot\}$分别表示实部和虚部。$\mathrm{arg}(\cdot)$和$|\cdot|$分别表示复数或复向量的相位和模。$\lfloor \cdot\rfloor$表示向下取整。

2.   信号模型

图1描绘了基于预编码矩阵调制的MIMO探通一体化系统框架^[19]。一体化系统通过置换矩阵P混排正交波形${\boldsymbol{\varPhi }}$位置调制通信信息，通过优化设计权重矩阵W得到发射信号矩阵${{\boldsymbol{S}}_i}$，限制密码本约束从而便于通信方向解调、防止窃听方向截获信息，同时抑制发射方向图旁瓣电平以保证探测性能。通信和窃听用户为解调通信信息，令接收信号与正交波形组${\boldsymbol{\varPhi }}$脉压。由于通信方向密码本各元素不同(图1中以PSK调制为例)，可通过恢复置换位置进行解调。而窃听用户方向脉压后信号元素全部相同，无法解调获取通信信息。后续首先介绍窄带MIMO体制探通一体化系统的发射信号模型和通信信息嵌入方法，进而介绍其发射波束方向图；然后针对通信用户端介绍了对应的接收信号模型。

图  1  基于预编码矩阵调制的MIMO探通一体化系统框架示意图

Figure  1.  The framework diagram of the MIMO DFRC system based on permutation matrix modulation

下载: 全尺寸图片 幻灯片

2.1   一体化系统发射信号模型

考虑探通一体化系统发射阵列为一个具有M个正交共置天线，间距为d的均匀线阵。令${\boldsymbol{\varPhi}} \in {\mathbb{C}^{K \times L}}$为正交离散波形矩阵，并满足所有时滞和多普勒位移下的正交条件，其中K和L分别为正交波形个数和每个脉冲的快拍数。则${\boldsymbol{\varPhi}}$满足^[22,23]

为了在MIMO一体化系统雷达探测任务中嵌入通信信息，采用预编码矩阵以置乱预设正交波形集${\boldsymbol{\varPhi}}$。则第i个脉冲的重构正交波形集可表示为

其中，${{\boldsymbol{P}}_i} \in {\mathbb{C}^{K \times K}}$为任意置换矩阵并满足

因此，单个脉冲嵌入信息量为

为了控制发射方向图和通信性能，采用权重矩阵${\boldsymbol{W}} = \left[ {{{\boldsymbol{w}}_1}\;{{\boldsymbol{w}}_2}\; \cdots \;{{\boldsymbol{w}}_K}} \right] \in {\mathbb{C}^{M \times K}}$，则基带发射信号可表示为

假设传播是非分散的，在第i个脉冲目标方向$\theta$上的合成信号为

其中，${\boldsymbol{a}}\left( \theta \right)$是空域导向矢量，表达式为

$\lambda$为波长。则探测信号在$\theta$方向上的能量，即发射方向图，可写作

其中，${\boldsymbol{w}} = {{\rm{vec}}} ({\boldsymbol{W}}),{\boldsymbol{A}}(\theta ) = {{\boldsymbol{I}}_K} \otimes {\boldsymbol{a}}(\theta )$。注意到$P(\theta )$与${{\boldsymbol{P}}_i}$无关，即通信信息的传输对于方向图性能无影响。

为了提高MIMO探通一体化系统探测目标回波SINR、增强系统的探测性能，MIMO探通一体化系统发射波束赋形旨在通过控制发射波形能量尽可能集中于主瓣区域，从而获得更大能量感兴趣方位目标信息，降低旁瓣区域辐射能量以减少信号相关干扰回波返回，是重要的探测性能指标之一。

2.2   通信接收信号模型

假设${{\boldsymbol{S}}_i}$是第i个脉冲发射波形矩阵，C个单天线通信用户分别位于$\theta _c^{{\text{com}}}$, $c = 1,2, \cdots ,C$。则第c个通信接收机的输入基带信号可表示为

其中，$\alpha _c^i,c = 1,2, \cdots ,C$是第i个脉冲MIMO探通一体化系统发射阵列与第c个通信接收机之间传播环境的信道系数；${\boldsymbol{n}}_c^i$为加性噪声，并建模为一个零均值圆对称高斯随机复向量。

然后，接收信号通过匹配滤波器${\boldsymbol{\varPhi}}$，可得

其中，${\boldsymbol{\hat n}}_c^i$是协方差矩阵为$\sigma _n^2{{\boldsymbol{I}}_K}$的加性噪声；且

表示通信密码本并满足$P\left( \theta \right) = {\left\| {{\boldsymbol{s}}\left( \theta \right)} \right\|^2}$。

可得${\boldsymbol{x}}_c^i$是${\boldsymbol{s}}\left( {\theta _c^{{\text{com}}}} \right)$的放缩、置乱和加噪。因此，用户端调制方式取决于${\boldsymbol{s}}\left( {\theta _c^{{\text{com}}}} \right)$，若${\boldsymbol{s}}\left( {\theta _c^{{\text{com}}}} \right)$为PSK密码本，则无噪下接收信号${\boldsymbol{x}}_c^i$可视为特殊的PSK调制信号，其中${\boldsymbol{s}}\left( {\theta _c^{{\text{com}}}} \right)$各元素均出现且只出现一次。解调的过程即为通过比较${\boldsymbol{x}}_c^i$和${\boldsymbol{s}}\left( {\theta _c^{{\text{com}}}} \right)$的各元素，确定置换顺序，从而恢复${{\boldsymbol{P}}_i}$。因此在一体化波形优化设计过程中须令通信方向${\boldsymbol{s}}\left( {\theta _c^{{\text{com}}}} \right)$固定以保证通信性能。此外，假设窃听用户方位角为$\theta _h^{{\text{eav}}},h = 1,2, \cdots ,E$，可在设计过程中令${\boldsymbol{s}}\left( {\theta _e^{{\text{eav}}}} \right)$各元素完全相同，从而令其无法在该位置进行解调。

3.   基于ADMM的一体化信号矩阵设计方法

根据上述探通一体框架和信号模型，本节提出了基于ADMM的一体化加权矩阵设计算法，通过密码本约束下最大化PMSR，以实现方向图赋形的同时传输通信信息。

3.1   问题模型

所有可能的方位角被分为以下几种离散格点集：主瓣${\mathcal{R}_{{\text{main}}}} = \left\{ {{\theta _i}} \right\}_{i = 1}^I$、旁瓣${\mathcal{R}_{{\text{side}}}} = \{ {\vartheta _s}\} _{s = 1}^S$、通信方向${\mathcal{R}_{{\text{com}}}} = \{ \theta _c^{{\text{com}}}\} _{c = 1}^C$、窃听用户方向${\mathcal{R}_{{\text{eav}}}} = \{ \theta _h^{{\text{eav}}}\} _{h = 1}^E$和过渡带${\mathcal{R}_{{\text{tran}}}}$。通过发射波束PMSR评估发射方向图性能，PMSR定义为^[24,25]

此外，在通信和窃听用户方向，加权矩阵需满足密码本约束，即

由式(8)可知，式(13)，式(14)满足

为调制通信信息，令${\boldsymbol{s}}\left( {\theta _c^{{\text{com}}}} \right) = {b_c}{{\boldsymbol{\bar a}}_1}$，其中${{\boldsymbol{\bar a}}_1}$可根据需求设定为PSK或QAM调制密码本并满足单位平均功率，即${\left\| {{{{\boldsymbol{\bar a}}}_1}} \right\|^2}{\text{ = }}K$(例如：令${{\boldsymbol{\bar a}}_1} = {[ {1\;{{{\rm{e}}}^{{{\rm{j}}}{\textstyle\frac{{{2\pi }}}{K}}}}\; \cdots \;{{{\rm{e}}}^{{{\rm{j}}}{\textstyle\frac{{{2\pi }(K - 1)}}{K}}}}} ]^{{\rm{T}}}}$为PSK调制端点分布向量)，且${b_c} \in \mathbb{R}$表示其幅度。通信用户接收信噪比(Signal-to-Noise Ratio, SNR)越小，其误符率(Symbol Error Rate, SER)越大。然而，若通信方向波束能量过高，用户设备则会出现接收信号幅度过饱和。因此，为满足不同距离通信用户接收信号功率需求，限制${b_c}$幅度范围为${l_c} \le {b_c} \le {u_c},c = 1,2, \cdots ,C$。此外，为了阻止位于$\theta _h^{{\text{eav}}},h = 1,2, \cdots ,E$的窃听用户截取通信信息的同时不影响发射方向图，令${\boldsymbol{s}}\left( {\theta _h^{{\text{eav}}}} \right)$中各元素相同，即${\boldsymbol{s}}\left( {\theta _h^{{\text{eav}}}} \right) = {r_h}{{{ {\textit{1}}}}_K},h = 1,2, \cdots ,E$。其中，${r_h} \in \mathbb{C}$为${\boldsymbol{s}}\left( {\theta _h^{{\text{eav}}}} \right)$的元素值。式(13)，式(14)可重写为

最后，该MIMO探通一体化系统的方向图设计问题可表示为

其中，$\log$函数是为了便于对分式目标函数处理^[24]，$\varDelta$为发射能量。问题(19)旨在抑制方向图峰值旁瓣，同时控制密码本相位和发射能量。接下来将引入基于ADMM的迭代算法求解该问题。

3.2   问题求解

首先，通过引入辅助变量$\eta$, $\epsilon$, ${{\boldsymbol{y}}_i}$, ${{\boldsymbol{z}}_s}$, ${{\boldsymbol{x}}_c}$, ${{\boldsymbol{v}}_h}$和n，问题(19)可等价转化为^[24]

因此，增广拉格朗日函数定义为

其中，${\boldsymbol{\rho }} = {[{\rho _1}\;{\rho _2}\;{\rho _3}\;{\rho _4}]^{\text{T}}}$为惩罚因子向量；${{\boldsymbol{\mu}} _i},{\text{ }}{{\boldsymbol{\iota}} _s}, {{\boldsymbol{\xi}} _c},{{\boldsymbol{\lambda}} _h},{\boldsymbol{\zeta }}$为对偶变量。

采用ADMM算法迭代最小化${L}_{{\boldsymbol{\rho}} }\left({\boldsymbol{w}},{{\boldsymbol{y}}}_{i},{{\epsilon}},{{\boldsymbol{z}}}_{s}, \eta ,{{\boldsymbol{x}}}_{c},{{\boldsymbol{v}}}_{h},{{\boldsymbol{\mu}} }_{i},{{\boldsymbol{\iota}} }_{s},{{\boldsymbol{\xi }}}_{c},{{\boldsymbol{\lambda}} }_{h},{\boldsymbol{n}},{\boldsymbol{\zeta}} \right)$以迭代更新${\boldsymbol{w}},{{\boldsymbol{y}}}_{i},{{\epsilon}},\text{ }{{\boldsymbol{z}}}_{s}, {{\eta}} ,{{\boldsymbol{x}}}_{c}, {{\boldsymbol{v}} }_{h},{{\boldsymbol{\mu}} }_{i}, \text{ }{{\boldsymbol{\iota}} }_{s},{{\boldsymbol{\xi}} }_{c},{{\boldsymbol{\lambda}} }_{h}$。令${{\boldsymbol{w}}}^{\left(t\right)},{{\boldsymbol{y}}}_{i}^{\left(t\right)},{{{\epsilon}}}^{\left(t\right)},\text{ }{{\boldsymbol{z}}}_{s}^{\left(t\right)},{{{\eta}} }^{\left(t\right)}, {{\boldsymbol{x}}}_{c}^{\left(t\right)}, {{\boldsymbol{v}}}_{h}^{\left(t\right)},{{\boldsymbol{\mu}} }_{i}^{\left(t\right)}, \text{ }{{\boldsymbol{\iota}} }_{s}^{\left(t\right)},{{\boldsymbol{\xi}} }_{c}^{\left(t\right)},{{\boldsymbol{\lambda}} }_{h}^{\left(t\right)},{{\boldsymbol{n}}}^{\left(t\right)},{{\boldsymbol{\zeta}} }^{\left(t\right)}$表示${\boldsymbol{w}},{{\boldsymbol{y}}}_{i},{{\epsilon}},\text{ }{{\boldsymbol{z}}}_{s}, {{\eta}} , {{\boldsymbol{x}}}_{c}, {{\boldsymbol{v}}}_{h},{{\boldsymbol{\mu}} }_{i},\text{ }{{\boldsymbol{\iota}} }_{s},{{\boldsymbol{\xi}} }_{c},{{\boldsymbol{\lambda}} }_{h},{\boldsymbol{n}},{\boldsymbol{\zeta}}$的第t次迭代值，具体迭代流程如算法1所示。注意到如果无窃听用户，即无约束(14)，仍可通过上述思路，令${\rho _3} = 0$并忽略式(26)和式(31)求解过程，以解决简化的${\mathcal{P}_0}$问题。通过算法1得到加权向量${{\boldsymbol{w}}^ \star }$后，随即可得MIMO一体化系统调制通信信息${{\boldsymbol{P}}_i}$时的发射信号矩阵${{\boldsymbol{S}}_i}$。算法1的收敛性参见文献[24]定理1。

算法1 基于ADMM的一体化信号矩阵设计方法
Alg. 1 DFRC waveform matrix design method based on ADMM
输入：$\left\{ {\boldsymbol{y} }_{i}^{(0)}\right\},{ {\epsilon} }^{(0)},\left\{ {\boldsymbol{z} }_{s}^{(0)}\right\},{ {\eta} }^{(0)},\left\{ {\boldsymbol{x} }_{c}^{(0)}\right\},\left\{ {\boldsymbol{v} }_{h}^{(0)}\right\},\left\{ {\boldsymbol{\mu} }_{i}^{(0)}\right\},\left\{ {\boldsymbol{\iota} }_{s}^{(0)}\right\},\left\{ {\boldsymbol{\xi} }_{c}^{(0)}\right\},\left\{ { {\boldsymbol{\lambda} } }_{h}^{(0)}\right\},{\boldsymbol{n} }^{\left(0\right)},{\boldsymbol{\zeta} }^{\left(0\right)},{\boldsymbol{\rho} },\varDelta ,{\delta }_{1}$；
输出：MIMO一体化系统加权向量${{\boldsymbol{\boldsymbol{w}}}^ \star }$；
步骤1. $t = 0$；
步骤2. 通过求解以下问题更新${\boldsymbol{w}}^{\left(t\text{+}1\right)},\left\{{\boldsymbol{y}}_{i}^{\left(t\text{+}1\right)}\right\},{{\epsilon}}^{\left(t\text{+}1\right)},\left\{{\boldsymbol{z}}_{s}^{\left(t\text{+}1\right)}\right\},{{\eta}}^{\left(t\text{+}1\right)},\left\{{\boldsymbol{x}}_{c}^{\left(t\text{+}1\right)}\right\},\left\{{\boldsymbol{v}}_{h}^{\left(t\text{+}1\right)}\right\},$${{\boldsymbol{\boldsymbol{n}}}^{\left( {t{\text{ + }}1} \right)}}$：
${\boldsymbol{w} }^{(t+1)}:=\mathrm{arg}\;\underset{\boldsymbol{w} }{\mathrm{min} }\;{L}_{ { {\boldsymbol \rho} } }\left({\boldsymbol{w}},\left\{ {\boldsymbol{y} }_{i}^{(t)}\right\},{ {\epsilon} }^{(t)},\left\{ {\boldsymbol{z} }_{s}^{(t)}\right\},{ {\eta} }^{(t)},\left\{ {\boldsymbol{x} }_{c}^{(t)}\right\},\left\{ {\boldsymbol{v} }_{h}^{(t)}\right\},\left\{ {\boldsymbol{\mu} }_{i}^{(t)}\right\},\left\{ {\boldsymbol{\iota} }_{s}^{(t)}\right\},\left\{ {\boldsymbol{\xi} }_{c}^{(t)}\right\},{ {\boldsymbol{\lambda} } }_{h}^{(t)},{\boldsymbol{n} }^{\left(t\right)},{\boldsymbol{\zeta} }^{\left(t\right)}\right)$ (22)
$\begin{aligned}& \left\{ {\boldsymbol{y} }_{i}^{(t+1)},{ {\epsilon} }^{(t+1)}\right\}:=\mathrm{arg}\;\underset{ {\boldsymbol{y} }_{i},\epsilon}{\mathrm{min} }\;{L}_{ { {\boldsymbol \rho} } }\left({\boldsymbol{w} }^{(t+1)},\left\{ {\boldsymbol{y} }_{i}\right\},\epsilon,\left\{ {\boldsymbol{z} }_{s}^{(t)}\right\},{ {\eta} }^{(t)}, \left\{ {\boldsymbol{x} }_{c}^{(t)}\right\},\left\{ {\boldsymbol{v} }_{h}^{(t)}\right\},\left\{ {\boldsymbol{\mu} }_{i}^{(t)}\right\},\left\{ {\boldsymbol{\iota} }_{s}^{(t)}\right\},\left\{ {\boldsymbol{\xi} }_{c}^{(t)}\right\},{ {\boldsymbol{\lambda} } }_{h}^{(t)},{\boldsymbol{n} }^{\left(t\right)},{\boldsymbol{\zeta} }^{\left(t\right)}\right)\text{ }\\ & \quad\text{s}\text{.t}\text{. }{\Vert {\boldsymbol{y} }_{i}\Vert }^{2}\ge \epsilon,i=1,2,\cdots ,I; \qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\quad\;\, (23)\end{aligned}$
$\begin{array}{l}\left\{ {\boldsymbol{z} }_{s}^{(t+1)},{ {\eta} }^{(t+1)}\right\}:=\mathrm{arg}\;\underset{ {\boldsymbol{z} }_{s},\eta }{\mathrm{min} }\;{L}_{ { {\boldsymbol \rho} } }\left({\boldsymbol{w} }^{(t+1)},\left\{ {\boldsymbol{y} }_{i}^{(t+1)}\right\},{ {\epsilon} }^{(t+1)},\left\{ {\boldsymbol{z} }_{s}\right\}, \eta ,\left\{ {\boldsymbol{x} }_{c}^{(t)}\right\},\left\{ {\boldsymbol{v} }_{h}^{(t)}\right\},\left\{ {\boldsymbol{\mu} }_{i}^{(t)}\right\},\left\{ {\boldsymbol{\iota} }_{s}^{(t)}\right\},\left\{ {\boldsymbol{\xi} }_{c}^{(t)}\right\},{ {\boldsymbol{\lambda} } }_{h}^{(t)},{\boldsymbol{n} }^{\left(t\right)},{\boldsymbol{\zeta} }^{\left(t\right)}\right)\text{ }\\ \quad\text{s}\text{.t}\text{. }{\Vert {\boldsymbol{z} }_{s}\Vert }^{2}\le \eta ,s=1,2,\cdots ,S; \qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\;\;\, (24)\end{array}$
$\begin{array}{l}\left\{ {\boldsymbol{x} }_{c}^{(t+1)},{b}_{c}^{(t+1)}\right\}:=\mathrm{arg}\underset{ {\boldsymbol{x} }_{c},{b}_{c} }{\mathrm{min} }{L}_{ { {\boldsymbol \rho} } }\left({\boldsymbol{w} }^{(t+1)},\left\{ {\boldsymbol{y} }_{i}^{(t+1)}\right\},{ {\epsilon} }^{(t+1)},\left\{ {\boldsymbol{z} }_{s}^{(t+1)}\right\}, { {\eta} }^{(t+1)},\left\{ {\boldsymbol{x} }_{c}\right\},\left\{ {\boldsymbol{v} }_{h}^{(t)}\right\},\left\{ {\boldsymbol{\mu} }_{i}^{(t)}\right\},\left\{ {\boldsymbol{\iota} }_{s}^{(t)}\right\},\left\{ {\boldsymbol{\xi} }_{c}^{(t)}\right\},{ {\boldsymbol{\lambda} } }_{h}^{(t)},{\boldsymbol{n} }^{\left(t\right)},{\boldsymbol{\zeta} }^{\left(t\right)}\right)\text{ }\\ \quad\text{s}\text{.t}\text{. }{\boldsymbol{x} }_{c}={b}_{c}{\overline{{\boldsymbol{a}}} }_{1},{l}_{c}\le {b}_{c}\le {u}_{c},c=1,2,\cdots ,C; \qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\, (25)\end{array}$
$\begin{aligned}& \left\{ {\boldsymbol{v} }_{h}^{(t+1)},{r}_{h}^{(t+1)} \right\} :=\mathrm{arg}\underset{ {\boldsymbol{v} }_{h},{r}_{h} }{\mathrm{min} }\;{L}_{ { {\boldsymbol \rho} } }\left({\boldsymbol{w} }^{(t+1)},\left\{ {\boldsymbol{y} }_{i}^{(t+1)}\right\},{ {\epsilon} }^{(t+1)},\left\{ {\boldsymbol{z} }_{s}^{(t+1)} \right\}, { {\eta} }^{(t+1)},\left\{ {\boldsymbol{x} }_{c}^{(t+1)} \right\},\left\{ {\boldsymbol{v} }_{h}\right\},\left\{ {\boldsymbol{\mu} }_{i}^{(t)}\right\},\left\{ {\boldsymbol{\iota} }_{s}^{(t)}\right\},\left\{ {\boldsymbol{\xi} }_{c}^{(t)} \right\},{ {\boldsymbol{\lambda} } }_{h}^{(t)},{\boldsymbol{n} }^{\left(t\right)},{\boldsymbol{\zeta} }^{\left(t\right)}\right)\\ & \quad\text{s}\text{.t}\text{. }{\boldsymbol{v} }_{h}={r}_{h} { {\textit{1}}}_{K},h=1,2,\cdots ,E \qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\;\, (26)\end{aligned}$
$\begin{aligned}& {\boldsymbol{n} }^{\left(t\text{+}1\right)}:=\mathrm{arg}\;\underset{\boldsymbol{n} }{\mathrm{min} }\;{L}_{ { {\boldsymbol \rho} } }\left({\boldsymbol{w} }^{(t+1)},\left\{ {\boldsymbol{y} }_{i}^{(t+1)}\right\},{ {\epsilon} }^{(t+1)},\left\{ {\boldsymbol{z} }_{s}^{(t+1)}\right\}, { {\eta} }^{(t+1)},\left\{ {\boldsymbol{x} }_{c}^{(t+1)}\right\},\left\{ {\boldsymbol{v} }_{h}^{(t+1)}\right\},\left\{ {\boldsymbol{\mu} }_{i}^{(t)}\right\},\left\{ {\boldsymbol{\iota} }_{s}^{(t)}\right\},\left\{ {\boldsymbol{\xi} }_{c}^{(t)}\right\},{ {\boldsymbol{\lambda} } }_{h}^{(t)},{\boldsymbol{n} },{\boldsymbol{\zeta} }^{\left(t\right)}\right)\\& \quad\text{s}\text{.t}\text{. }{\Vert {\boldsymbol{n} }\Vert }^{2}=\varDelta \qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad \quad\;\, (27)\end{aligned}$
步骤3. 通过下列公式更新$\left\{ {{\boldsymbol{\boldsymbol{\mu}}}_i^{(t + 1)}} \right\},\left\{ {{\boldsymbol{\boldsymbol{\iota}}}_s^{(t + 1)}} \right\},\left\{ {{\boldsymbol{\boldsymbol{\xi}}}_c^{(t + 1)}} \right\},\left\{ {{\boldsymbol{\lambda }}_h^{(t + 1)}} \right\},{{\boldsymbol{\boldsymbol{\zeta}}}^{\left( {t{\text{ + }}1} \right)}}$：
${\boldsymbol{\boldsymbol{\mu}}}_i^{(t + 1)}: = {\boldsymbol{\boldsymbol{\mu}}}_i^{(t)} + {\boldsymbol{\boldsymbol{y}}}_i^{(t + 1)} - {{\boldsymbol{A}}^{\text{H}}}\left( {{\theta _i}} \right){{\boldsymbol{\boldsymbol{w}}}^{(t + 1)}}$ (28)
${\boldsymbol{\boldsymbol{\iota}}}_s^{(t + 1)}: = {\boldsymbol{\boldsymbol{\iota}}}_s^{(t)} + {\boldsymbol{\boldsymbol{z}}}_s^{(t + 1)} - {{\boldsymbol{A}}^{\text{H}}}\left( {{\vartheta _s}} \right){{\boldsymbol{\boldsymbol{w}}}^{(t + 1)}}$ (29)
${\boldsymbol{\boldsymbol{\xi}}}_c^{(t + 1)}: = {\boldsymbol{\boldsymbol{\xi}}}_c^{(t)} + {\boldsymbol{\boldsymbol{x}}}_c^{(t + 1)} - {{\boldsymbol{A}}^{\text{H}}}\left( {\theta _c^{{\text{com }}}} \right){{\boldsymbol{\boldsymbol{w}}}^{(t + 1)}}$ (30)
${\boldsymbol{\lambda }}_h^{(t + 1)}: = {\boldsymbol{\lambda }}_h^{(t)} + {\boldsymbol{\boldsymbol{v}}}_h^{(t + 1)} - {{\boldsymbol{A}}^{\text{H}}}\left( {\theta _h^{{\text{eav }}}} \right){{\boldsymbol{\boldsymbol{w}}}^{(t + 1)}}$ (31)
${{\boldsymbol{\boldsymbol{\zeta}}}^{\left( {t{\text{ + }}1} \right)}} = {{\boldsymbol{\boldsymbol{\zeta}}}^{\left( t \right)}} + {{\boldsymbol{\boldsymbol{n}}}^{(t + 1)}} - {{\boldsymbol{\boldsymbol{w}}}^{(t + 1)}}$ (32)
步骤4. 如果原始可行性容差$V_1^{(t + 1)} = \displaystyle\sum\limits_{i = 1}^I {{{\left\| {{\boldsymbol{\boldsymbol{y}}}_i^{(t + 1)} - {{\boldsymbol{A}}^{\text{H}}}\left( {{\theta _i}} \right){{\boldsymbol{\boldsymbol{w}}}^{(t + 1)}}} \right\|}^2}} +$$\displaystyle\sum\limits_{s = 1}^S {{{\left\| {{\boldsymbol{\boldsymbol{z}}}_s^{(t + 1)} - {{\boldsymbol{A}}^{\text{H}}}\left( {{\vartheta _s}} \right){{\boldsymbol{\boldsymbol{w}}}^{(t + 1)}}} \right\|}^2}} +$$\displaystyle\sum\limits_{c = 1}^C { { {\left\| { {\boldsymbol{\boldsymbol{x} } }_c^{(t + 1)} - { {\boldsymbol{A} }^{\text{H} } }\left( {\theta _c^{ {\text{com } } } } \right){ {\boldsymbol{\boldsymbol{w} } }^{(t + 1)} } } \right\|}^2} }$ 　$+\displaystyle\sum\limits_{h = 1}^E { { {\left\| { {\boldsymbol{\boldsymbol{v} } }_h^{(t + 1)} - { {\boldsymbol{A} }^{\text{H} } }\left( {\theta _h^{ {\text{eav } } } } \right){ {\boldsymbol{\boldsymbol{w} } }^{(t + 1)} } } \right\|}^2} }$$+ {\left\| { { {\boldsymbol{\boldsymbol{n} } }^{(t + 1)} } - { {\boldsymbol{\boldsymbol{w} } }^{(t + 1)} } } \right\|^2} \le {\delta _1}$，则输出${{\boldsymbol{\boldsymbol{w}}}^ \star } = {{\boldsymbol{\boldsymbol{w}}}^{\left( {t + 1} \right)}}$；否则$t: = t + 1$，回到步骤2。

下载: 导出CSV 

| 显示表格

(1) 更新${{\boldsymbol{w}}^{(t + 1)}}$

问题(22)可等价变换为

其中

其中，$\mathcal{R} = {\mathcal{R}_{{\text{main}}}}\bigcup {{\mathcal{R}_{{\text{side}}}}} \bigcup {{\mathcal{R}_{{\text{com}}}}} \bigcup {{\mathcal{R}_{{\text{eav}}}}}$。因此，该问题的闭式解为

(2) 更新${{\boldsymbol{y}}}_{i}^{(t+1)},{\epsilon}^{(t+1)},{{\boldsymbol{z}}}_{s}^{(t+1)},{\eta }^{(t+1)}$

忽略无关常数项，问题(23)和问题(24)可分别等价转换为^[24]

和

其中，${\boldsymbol{\bar y}}_i^{(t)} = {{\boldsymbol{A}}^{\text{H}}}\left( {{\theta _i}} \right){{\boldsymbol{w}}^{(t + 1)}} - {\boldsymbol{\mu }}_i^{(t)}$, ${\boldsymbol{\bar z}}_s^{(t)} = {{\boldsymbol{A}}^{\text{H}}}\left( {{\vartheta _s}} \right) \cdot{{\boldsymbol{w}}^{(t + 1)}} - {\boldsymbol{\iota }}{\kern 1pt} _s^{(t)}$。

首先，针对问题(39)，一旦获得${\epsilon}^{(t+1)}$，$\{ {\boldsymbol{y}}_i^{(t + 1)}\}$即可通过${\boldsymbol{\bar y}}_i^{(t)}$放缩得到，即

然后，将式(41)反代回式(39)，即可得到关于$\epsilon$的优化问题

其中

相似的，给定${\eta ^{(t + 1)}}$，$\{ {\boldsymbol{z}}_s^{(t + 1)}\}$可通过${\boldsymbol{\bar z}}_s^{(t)}$放缩得到，即

将式(45)反代回式(40)可得

其中

注意到问题(42)和问题(46)具有相似结构，均为无约束优化问题，可通过文献[24]中方法求得其闭式解。

(3) 更新${\boldsymbol{x}}_c^{(t + 1)},b_c^{(t + 1)},{\boldsymbol{v}}_h^{(t + 1)},r_h^{(t + 1)}$

问题(25)可以重写为

其中，${\boldsymbol{\bar x}}_c^{(t)} = {{\boldsymbol{A}}^{\text{H}}}\left( {\theta _c^{{\text{com }}}} \right){{\boldsymbol{w}}^{(t + 1)}} - {\boldsymbol{\xi }}_c^{(t)}$。注意到目标函数是关于${b_c}$的二次函数，令

则该问题闭式解为

则有${\boldsymbol{x}}_c^{(t + 1)} = b_c^{(t + 1)}{{\boldsymbol{\bar a}}_1}$。

相似地，问题(26)可重写为

其中，${\boldsymbol{\bar v}}_h^{(t)} = {{\boldsymbol{A}}^{\text{H}}}\left( {\theta _h^{{\text{eav }}}} \right){{\boldsymbol{w}}^{(t + 1)}} - {\boldsymbol{\lambda }}_h^{(t)}$。该问题的闭式解为

则有${\boldsymbol{v}}_h^{(t + 1)} = r_h^{(t + 1)}{{{ {\textit{1}}}}_K}$。

(4) 更新${{\boldsymbol{n}}^{\left( {t{\text{ + }}1} \right)}}$

问题(27)可等价转化为

其中，${{\boldsymbol{\bar n}}^{\left( t \right)}} = {{\boldsymbol{w}}^{\left( {t + 1} \right)}} - {{\boldsymbol{\zeta}} ^{\left( t \right)}}$。其闭式解为

3.3   计算复杂度和收敛性

算法1的计算复杂度分析如下：${{\boldsymbol{R}}^{ - 1}}$和d可在算法开始前计算并保存，其计算复杂度分别为$O(Q{M^2} + {(KM)^{2.373}})$^[24]和$O(QMK)$，其中$Q = I + S + C + E$为约束个数。每次迭代中，式(22)—式(31)可用简单的闭式解来解决。因此，只需要基本的矩阵到向量的乘法。具体的，式(22)、式(23)—式(26)、式(28)—式(31)的计算复杂度分别为$O({K^2}{M^2})$, $O(KQ)$和$O(KMQ)$。因此，算法1的总计算复杂度为$O(Q{M^2} + {(KM)^{2.373}} + {T_0}({K^2}{M^2} + KMQ))$，其中${T_0}$为迭代次数。

4.   基于ADPM的排序学习解调方法

针对通信接收信号解调问题，本节提出了基于ADPM的排序学习解调方法以优化求解置换矩阵。假设每个通信接收机都完全已知正交波形矩阵${\boldsymbol{\varPhi }}$和${{\boldsymbol{\bar a}}_1}$，解调过程可建模为预编码矩阵P的优化问题^[19]

其中，$\varLambda = \{ {\boldsymbol{P}}|{\boldsymbol{P}}(i,j) \in \{ 0,1\} , {\boldsymbol{P}}{{{ {\textit{1}}}}_K} = {{{ {\textit{1}}}}_K},{{\boldsymbol{P}}^{\text{T}}} {{{ {\textit{1}}}}_K} = {{{ {\textit{1}}}}_K}\}$定义了置换矩阵集合。

注意到问题(56)是一个混合-布尔优化问题，文献[19]中通过穷举法解决该问题，然而当K很大时穷举法无法应对维度爆炸。为实现快速有效信息解调，首先将目标函数重新参数化为

其中，$b = {\left\| {{\boldsymbol{x}}_c^i} \right\|^2} + P\left( {\theta _c^{{\text{com}}}} \right)$, ${\boldsymbol{g}} = \Re \left\{ {{\boldsymbol{x}}_c^i \otimes {\boldsymbol{\bar a}}_1^*} \right\} \in {\mathbb{R}^{{K^2}}}$, ${\boldsymbol{p}} = {{\rm{vec}}} ({\boldsymbol{P}}) \in {\mathbb{R}^{{K^2}}}$。进一步地，引入辅助变量${{\boldsymbol{u}}_1}, {{\boldsymbol{u}}_2} \in {\mathbb{R}^{{K^2}}}$，问题(56)可等价转化为

其中，${\boldsymbol{B}} = {{ {\textit{1}}}}_K^{\text{T}} \otimes {{\boldsymbol{I}}_K},{\boldsymbol{C}} = {{\boldsymbol{I}}_K} \otimes {{ {\textit{1}}}}_K^{\text{T}}$。因此，增广拉格朗日函数可定义为

其中，${{\boldsymbol{\kappa}} _1}$和${{\boldsymbol{\kappa}} _2}$为对偶变量。

令${{\boldsymbol{p}}^l},{\boldsymbol{u}}_1^l,{\boldsymbol{u}}_2^l,{\boldsymbol{\kappa}} _1^l,{\boldsymbol{\kappa}} _2^l,{\rho ^l}$表示第l次迭代${\boldsymbol{p}},{\boldsymbol{u}}_1^{},{\boldsymbol{u}}_2^{}, {\boldsymbol{\kappa}} _1^{},{\boldsymbol{\kappa}} _2^{},\rho$的值，基于ADPM的排序学习优化解调方法流程如算法2所示。ADPM和ADMM算法的主要区别在于惩罚因子$\rho$是否随迭代改变。由于问题(60)是混合布尔优化问题，${{\boldsymbol{p}}^l}$易陷入局部最优解，若采用ADMM算法则此后各元素在迭代过程中无更新。因此采用了ADPM算法框架求解解调优化问题，令${\text{ }}\left\| {{{\boldsymbol{p}}^{l + 1}} - {{\boldsymbol{p}}^l}} \right\| = 0$时，增大惩罚因子${\rho ^l}$，从而增强原变量与辅助变量的约束并加速收敛。

算法2 基于ADPM的排序学习优化解调方法
Alg. 2 The permutation learning demodulation method based on ADPM
输入：${{\boldsymbol{p}}^0},{\boldsymbol{u}}_1^0,{\boldsymbol{u}}_2^0,{\boldsymbol{\kappa }}_1^0,{\boldsymbol{\kappa }}_2^0,{\rho ^0},{\boldsymbol{g}},{\boldsymbol{B}},{\boldsymbol{C}},{\delta _2},{\delta _3}$；
输出：问题(56)的最优解${{\boldsymbol{p}}^ \star }$；
步骤1. $l = 0$；
步骤2. 通过求解以下问题更新${{\boldsymbol{p}}^{l{\text{ + }}1}},{\boldsymbol{u}}_1^{l{\text{ + }}1},{\boldsymbol{u}}_2^{l{\text{ + }}1}$： 　$\begin{aligned} { {\boldsymbol{p} }^{l + 1} }: =& \arg \mathop {\min }\limits_{\boldsymbol{p} } {\mathcal{L}_{ {\rho ^l} } }\left( { {\boldsymbol{p} },{\boldsymbol{u} }_1^l,{\boldsymbol{u} }_2^l,{\boldsymbol{\kappa } }_1^l,{\boldsymbol{\kappa } }_2^l} \right){\text{ } } \\& {\text{s} }{\text{.t} }{\text{. } }{p_i} \in \{ 0,1\} ,i = 1,2, \cdots ,{K^2} \qquad\qquad\qquad\qquad\; (60)\end{aligned}$
$\begin{aligned} {\boldsymbol{u} }_1^{l + 1}: =& \arg \mathop {\min }\limits_{ {\boldsymbol{u} }_1^{} } {\mathcal{L}_{ {\rho ^l} } }\left( { { {\boldsymbol{p} }^{l + 1} },{\boldsymbol{u} }_1^{},{\boldsymbol{u} }_2^l,{\boldsymbol{\kappa } }_1^l,{\boldsymbol{\kappa } }_2^l} \right){\text{ } } \\ & {\text{s} }{\text{.t} }{\text{. } }{\boldsymbol{B} }{ {\boldsymbol{u} }_1} = { { { {\textit{1} } } }_K} \qquad\qquad\qquad\qquad\qquad\qquad\qquad\quad\;\, (61)\end{aligned}$
$\begin{aligned} {\boldsymbol{u} }_2^{l + 1}: =& \arg \mathop {\min }\limits_{ { {\boldsymbol{u} }_2} } {\mathcal{L}_{ {\rho ^l} } }\left( { { {\boldsymbol{p} }^{l + 1} },{\boldsymbol{u} }_1^{l + 1},{\boldsymbol{u} }_2^{},{\boldsymbol{\kappa } }_1^l,{\boldsymbol{\kappa } }_2^l} \right){\text{ } } \\ & {\text{s} }{\text{.t} }{\text{. } }{\boldsymbol{C} }{ {\boldsymbol{u} }_2} = { { { {\textit{1} } } }_K} \qquad\qquad\quad\qquad\qquad\qquad\qquad\quad\; (62)\end{aligned}$
步骤3. 通过下列公式更新${\boldsymbol{\kappa }}_1^{l + 1},{\boldsymbol{\kappa }}_2^{l + 1},{\rho ^{l + 1}}$：
${\rho }^{l+1}=\left\{\begin{aligned}& {\delta }_{2}{\rho }^{l},\text{ }\Vert {p}^{l+1}-{p}^{l}\Vert =0 \\ & {\rho }^{l},\text{ }其他\text{ } \end{aligned}\right.$ (63)
${\boldsymbol{\kappa }}_1^{l + 1}: = {\boldsymbol{\kappa }}_1^l + {\rho ^{l + 1}}\left( {{\boldsymbol{u}}_1^{l + 1} - {{\boldsymbol{p}}^{l + 1}}} \right)$ (64)
${\boldsymbol{\kappa }}_2^{l + 1}: = {\boldsymbol{\kappa }}_2^l + {\rho ^{l + 1}}\left( {{\boldsymbol{u}}_2^{l + 1} - {{\boldsymbol{p}}^{l + 1}}} \right)$ (65)
4. 如果原始可行性容差 　${V}_2^{l + 1} = {\left\| { {\boldsymbol{u} }_1^{l + 1} - { {\boldsymbol{p} }^{l + 1} } } \right\|^2} + {\left\| { {\boldsymbol{u} }_2^{l + 1} - { {\boldsymbol{p} }^{l + 1} } } \right\|^2} \le {\delta _3}$，则输出 　${{\boldsymbol{p}}^ \star } = {{\boldsymbol{p}}^{l + 1}}$；否则$l: = l + 1$，回到步骤2。

下载: 导出CSV 

| 显示表格

(1) 更新${{\boldsymbol{p}}^{l{\text{ + }}1}}$

忽略常数项，问题(60)可重写为

其中，${f_1}({\boldsymbol{p}}) = {\rho ^l}{{\boldsymbol{p}}^{\text{T}}}{\boldsymbol{p}} - {({\boldsymbol{g}}_2^l)^{{\text{T}}}}{\boldsymbol{p}} = \displaystyle\sum\limits_{i = 1}^{{K^2}} {\left( {{\rho ^l}p_i^2 - g_{2,i}^l{p_i}} \right)}$, ${\boldsymbol{g}}_2^l = {\boldsymbol{g}} + {\rho ^l}\left( {{\boldsymbol{u}}_1^l + {\boldsymbol{u}}_2^l} \right) + {\boldsymbol{\kappa}} _1^l + {\boldsymbol{\kappa}} _2^l$。

观察问题(66)可知，目标函数和约束对于p中各元素相互独立，则式(66)的最优解为

其中，${\boldsymbol{\tilde p}} = {\boldsymbol{g}}_2^l/(2{\rho ^l})$。

(2) 更新${\boldsymbol{u}}_1^{l{\text{ + }}1},{\boldsymbol{u}}_2^{l{\text{ + }}1}$

问题(61)可重写为

将等式约束并入增广拉格朗日函数中，可得

则由KKT条件可得优化问题(68)的闭式解

另外，问题(62)可以重写为

类似可得问题(62)的最优解为

(3) 计算复杂度

算法2的计算复杂度与迭代次数以及P的大小有关。每次迭代中，${\boldsymbol{\tilde p}}$和更新${\boldsymbol{\kappa }}_1^{l + 1},{\boldsymbol{\kappa }}_2^{l + 1}$所需计算复杂度为$O({K^2})$，更新${\boldsymbol{u}}_1^{l{\text{ + }}1}$和${\boldsymbol{u}}_2^{l{\text{ + }}1}$所需计算复杂度为$O({K^3})$。因此，每次迭代所需总计算复杂度为$O({K^3})$。

5.   仿真实验

本节主要从发射波束方向图、星座图、SER以及收敛性等方面对所提出的算法框架进行性能评估。考虑具有$M = 10$个发射天线且间隔为半波长的均匀线阵MIMO探通一体化系统，空间角间隔1°均匀划分离散格点，${\mathcal{R}_{{\text{main}}}} = [ - {10^\circ },{10^\circ }]$。令算法1退出条件为${\delta _1} = K \times {10^{ - 8}}$；惩罚因子为${\rho _1} = {\rho _2} = 1$；当$E \ge 1$时，${\rho _3} = 1$，否则${\rho _3} = 0$；${\rho _4} = 10$。w的能量设为$\varDelta = M$。

5.1   单通信用户探通一体化性能

当单通信用户采用PSK密码本且不考虑通信接收信号幅度约束时，可由发射方向图不变和选择(Transmit Radiation Pattern Invariance and Selection, TRPIS)方法^[19]实现。其中，TRPIS的波束形成权向量由文献[26]生成，其母权向量采用最小最大准则优化，旁瓣保持在20 dB以下。

因此，首先考虑无窃听用户($E = 0$)时，一个通信接收机(即$C = 1$)位于$\theta _1^{{\text{com}}} = - {60^\circ }$, ${\mathcal{R}_{{\text{tran}}}} = \left[ { - {{68}^\circ }, - {{52}^\circ }} \right] \cup \left[ { - {{18}^\circ }, - {{11}^\circ }} \right] \cup \left[ {{{11}^\circ },{{18}^\circ }} \right]$, ${\mathcal{R}_{{\text{side }}}} = \left[ { - {{90}^\circ }, - {{69}^\circ }} \right] \cup \left[ { - {{51}^\circ }, - {{19}^\circ }} \right] \cup \left[ {{{19}^\circ },{{90}^\circ }} \right]$。假设一体化系统探测主瓣方向能量较大，通信用户方向所需满足能量约束参数为${l}_{1}\text{=0}\text{.05},\;{u}_{1}\text{=}0.1$。图2(a)—图2(c)描述了$K = 4,8,16$时PSK和QAM调制下所提ADMM算法(简记为PSK-ADMM,QAM-ADMM)，以及PSK调制下TRPIS算法(简记为PSK-TRPIS)所得发射方向图。其中，通信方向的能量上限和下限分别采用“+”和“×”标记，并将不同方法所得PMSR值和运行时间分别记录于表1和表2中。TRPIS算法从${2^{M - 1}} - 1$个复数根中选择K个，与母权向量相乘，形成相同的方向图和通信方向不同的相位角。因此，TRPIS在不同K下所得方向图相同，且时间相近。仿真结果表明，不同K下PSK-ADMM和QAM-ADMM均可控制通信方向发射能量且所得方向图相似，并可在不设定方向图模板下耗费与TRPIS算法相近时间获得更高的PMSR。

图  2  假设$\theta _1^{{\text{com}}} = - {60^\circ }$时发射方向图

Figure  2.  The transmit beampattern with $\theta _1^{{\text{com}}} = - {60^\circ }$

下载: 全尺寸图片 幻灯片

表  1  不同算法所得PMSR

Table  1.  PMSR derived by different methods

数值	PSK-ADMM	QAM-ADMM	PSK-TRPIS
$K = 4$	18.7693	18.7715	17.8392
$K = 8$	18.7693	18.7735
$K = 16$	18.7535	18.8135

下载: 导出CSV 

| 显示表格

表  2  不同算法所需时间

Table  2.  Computing time required by different methods

数值	PSK-ADMM	QAM-ADMM	PSK-TRPIS
$K = 4$	8.354	8.434	8.276
$K = 8$	5.502	5.182	8.281
$K = 16$	9.221	7.398	8.282

下载: 导出CSV 

| 显示表格

$K = 4,8,16$时PSK-ADMM,QAM-ADMM和PSK-TRPIS所得星座图(即极坐标下密码本幅度$\left| {{\boldsymbol{s}}\left( {\theta _1^{{\text{com}}}} \right)} \right|$随相位$\arg \left( {{\boldsymbol{s}}\left( {\theta _1^{{\text{com}}}} \right)} \right)$的分布)分别如图3(a)—图3(c)所示。其中，${\left\| {{\boldsymbol{s}}\left( {\theta _1^{{\text{com}}}} \right)} \right\|^2}$与图2中$P(\theta _1^{{\text{com}}})$的大小一致。由于PSK-ADMM和QAM-ADMM的设计中均考虑了约束(17)，因此所得的参考字典为PSK和QAM端点分布的放缩，且满足通信方向能量约束${l_1} \le {b_1} \le {u_1}$。而PSK-TRPIS算法由于采用了复根选择方案，与PSK-ADMM相比，${\boldsymbol{s}}\left( {\theta _1^{{\text{com}}}} \right)$端点分布略有偏差，且无法令${\boldsymbol{s}}\left( {\theta _1^{{\text{com}}}} \right)$各元素实现不同幅度，因此不能采用8QAM和16QAM调制。

图  3  通信星座图

Figure  3.  Communication constellation diagram

下载: 全尺寸图片 幻灯片

接下来，随机发射${10^4}$个符号对应一体化信号并评估所提算法通信性能。$K = 4,8,16$时总的符号个数分别为4!=24, 8!=40320和16!>2⁴⁴，文献[19]中所提穷举法以求解解调问题(56)在K较大时无法应对维度灾难问题，因此无法实现有效快速解调。定义第c个通信用户接收端SNR为${\text{SN}}{{\text{R}}_c} = {\left| {{\alpha _c}{b_c}} \right|^2}/ \sigma _n^2$，本文所提算法1(PSK-ADMM和QAM-ADMM)和TRPIS算法均采用本文所提算法2以恢复置换顺序，相关参数设为${\rho ^0} = 0.01$, ${\delta _2} = 1.01$和${\delta _3} = {10^{ - 4}}$。不同一体化信号通信用户端符号SER随SNR的变化曲线如图4所示。其中，“KPSK”，K=4, 8, 16标记的黑色实线表示在接收端采用PSK直接解调式(10)并恢复${{\boldsymbol{\bar a}}_1}$置换顺序所得SER结果。正如预期，SNR的增加使得SER性能改善。算法2利用了置换矩阵特性，通过最大化字典${\boldsymbol{\bar a}}_1^{}$与${\boldsymbol{Px}}_c^i$的互相关实现了SNR增益，因此相比PSK直接解调，在相同SNR下可得更小SER，验证了所提算法2在排序学习优化问题上的有效性。由于PSK信号星座图和QAM信号矩形星座图点间最小距离特性，当K=4, 8时，PSK-ADMM, QAM-ADMM和PSK-TRPIS算法所得SER随SNR变化曲线相近；当K=16时，QAM-ADMM性能优于PSK-ADMM，由于TRPIS所得参考密码本相位具有一定偏差，因此SER高于PSK-ADMM。

图  4  不同算法SER随SNR变化曲线

Figure  4.  SER versus SNR for different algorithms

下载: 全尺寸图片 幻灯片

以4QAM-ADMM为例，图5描绘了SNR=4 dB时不同角度解调SER。由于通信信息嵌入在$- {60^\circ }$中，因此在该方向具有最低的误码率。如图6所示，${\boldsymbol{s}}\left( {{{65}^\circ }} \right)$与${\boldsymbol{s}}\left( { - {{60}^\circ }} \right)$元素分布(即$\left| {{\boldsymbol{s}}\left( \theta \right)} \right|$随$\arg \left( {{\boldsymbol{s}}\left( \theta \right)} \right)$变化)相似，因此该方向解调也可得较小误码率。因此，若该方向存在窃听用户，则可在较大信噪比下解调调制信息。

图  5  SER随角度变化曲线

Figure  5.  SER versus angle

下载: 全尺寸图片 幻灯片

图  6  ${\boldsymbol{s}}\left( {{{65}^\circ }} \right)$与${\boldsymbol{s}}\left( { - {{60}^\circ }} \right)$元素分布

Figure  6.  The element distributions of ${\boldsymbol{s}}\left( {{{65}^\circ }} \right)$ and ${\boldsymbol{s}}\left( { - {{60}^\circ }} \right)$

下载: 全尺寸图片 幻灯片

图7(a)展示了不同K下算法1分别采用PSK和QAM调制所得原始可行性容差$V_1^{(t)}$随运行时间变化情况。图7(b)进一步描绘了当SNR=4 dB时算法2解调各发射信号所得$V_2^l$随运行时间变化曲线。仿真结果说明本文所提算法均逐渐收敛，因此满足约束条件。其中，算法2在K=16时运行时间小于0.25 s，当K=4, 8时，运行时间小于0.05 s，因此该算法适用于通信系统实时性需求。

图  7  算法1与算法2计算复杂度

Figure  7.  The computational complexities of Alg. 1 and Alg. 2

下载: 全尺寸图片 幻灯片

5.2   窃听下多通信用户探通一体化性能分析

考虑窃听用户位于$\theta _1^{{\text{eav}}} = {65^\circ }$，通信接收机分别位于$\theta _1^{{\text{com}}} = - {60^\circ }$和$\theta _2^{{\text{com}}} = {40^\circ }$, ${\mathcal{R}_{{\text{tran}}}} = \left[ { - {{68}^\circ }, - {{52}^\circ }} \right] \cup \left[ { - {{18}^\circ }, - {{11}^\circ }} \right] \cup \left[ {{{11}^\circ },{{18}^\circ }} \right] \cup \left[ {{{32}^\circ },{{48}^\circ }} \right]$, ${\mathcal{R}_{{\text{side }}}} = \left[ - {{90}^\circ }, - {{69}^\circ } \right] \cup \left[ { - {{51}^\circ }, - {{19}^\circ }} \right]$$\cup \left[ {{{19}^\circ },{{31}^\circ }} \right] \cup \left[ {{{49}^\circ },{{90}^\circ }} \right]$。假设两个通信用户方向与一体化系统距离不同，所需满足能量约束参数分别为${l_1}{\text{ = 0}}{\text{.05,}}\;{u_1}{\text{ = }}0.1,{l_2}{\text{ = 1}}, {u_2}{\text{ = }}2$。此时，TRPIS算法无法同时优化发射波束性能和控制通信与窃听用户参考密码本。本节以K=4下PSK调制为例，验证了所提算法的有效性。

图8(a)，图8(b)分别描绘了采用算法1所得发射方向图和${\boldsymbol{s}}\left( \theta \right),\theta \in \left\{ {\theta _c^{{\text{com}}},\theta _e^{{\text{eav}}}} \right\}$分布。由仿真结果可得，图8(a)中PMSR为16.54 dB，设计的一体化加权矩阵可实现不同方向通信用户发射不同能量通信信号的同时抑制方向图旁瓣电平，从而满足了雷达探测和不同距离的通信用户SNR需求。由图8(b)可得，算法1可精确控制密码本${\boldsymbol{s}}\left( \theta \right),\theta \in \left\{ {\theta _c^{{\text{com}}},\theta _e^{{\text{eav}}}} \right\}$分布，以阻止窃听并在通信方向实现解调。

图  8  探通一体化性能

Figure  8.  DFRC performance

下载: 全尺寸图片 幻灯片

图9展示了$\theta _c^{{\text{com}}},c = 1,2,\theta _1^{{\text{eav}}}$方向采用算法2的解调SER随SNR变化曲线。正如预期，在相同的调制分布${{\boldsymbol{\bar a}}_1}$和解调方法下，不同方向用户SER性能几乎相同，而由于窃听方向密码本各元素完全相同，窃听方向无法解调通信信息。

图  9  SER随SNR变化曲线

Figure  9.  SER versus SNR

下载: 全尺寸图片 幻灯片

6.   结语

针对基于预编码矩阵调制的MIMO探通一体化系统，本文提出了基于ADMM的一体化信号矩阵优化设计方法和基于ADPM的排序学习优化解调方法。首先，建立了通信和窃听用户方向密码本约束下最大化方向图PMSR优化问题；然后，引入了辅助变量将二次分式耦合问题转换为多个有闭式解的二次优化子问题，提出了基于ADMM的一体化加权矩阵优化设计方法并分析了该算法的计算复杂度和收敛性。最后，基于通信发射调制解调机理，建立了排序学习混合布尔优化问题，提出了基于ADPM的排序学习优化解调方法提升了解调效率。仿真结果评估了所提一体化信号矩阵设计和解调方法的有效性。后续研究过程将考虑宽带MIMO^[27,28]和信道估计误差下的稳健MIMO一体化系统波形设计^[29]，以及非完全正交波形集带来的性能损失^[30]。