刘紫玉 亓丽梅2)† 道日娜 戴林林 武利勤

1) (北京邮电大学电子工程学院,北京 100876)

2) (山东师范大学,光场调控及应用协同创新中心,济南 250358)

3) (北京市科学技术研究院,北京市辐射技术研究所,北京 100875)

基于VO2 的相变特性,提出了仅用两种混合结构实现2-bit (四种状态)编码的太赫兹编码超表面.结构单元中贯通的金属线可用作电压引线对单行阵列进行控制,使得固定阵列结构能够呈现不同的状态编码,实现对波束的动态调控.此外,采用MATLAB 软件对编码超表面阵列天线进行了可视化设计,通过对工作频率、波束偏转角度等参数的设置,实现了对状态序列与辐射结果的预测.该可视化系统不限于具体的结构单元,对一切满足编码条件的阵列均具有普适性.最后,采用深度神经网络进行了逆向天线设计,通过与模拟对比验证了其在波束偏转角度和单元排布的有效性.本文为主动灵活调控太赫兹波提供了新途径,在太赫兹成像、相控雷达、通信等领域具有潜在的应用价值.

1 引言

太赫兹波是指频率范围在0.1—10.0 THz 的电磁波.近年来,太赫兹科学技术的应用已趋于成熟,在生物医学、安全检测、高速通信、隐身技术等方面的潜在应用受到广泛关注[1-3].由于太赫兹波频率高、器件结构尺寸小,有源器件难以集成,器件的局部操控变得非常困难.东南大学崔铁军团队[4,5]提出编码超表面的概念,为实现多样化的太赫兹波束调控提供了新的思路.传统超表面采用尺寸诱导或Pancharatnam-Berry (PB) 相位诱导实现波前控制与编码,通常n种编码状态需要n个结构单元,一旦加工完成,其功能也固定[6-8].为实现一个单元的状态调控,研究人员提出将石墨烯、液晶、Ge2Sb2Te5(GST)、VO2等可调谐材料引入到超表面设计中[9-12].考虑到实际加工的可操作性和调控的便利性,本文选择VO2作为可调谐材料进行编码单元设计.VO2是一种特殊的相变材料,其相变温度约为68 ℃,且在相变前后伴随着显著的电磁和光学特性突变.相变前,VO2呈现绝缘态,相变后,VO2呈现金属态[13-15],这一性质引起国内外研究者的广泛讨论,并将其引入太赫兹编码超表面的设计中.已有研究者在圆环[16,17]、十字双箭头[18]、方形环[19]等金属贴片中加入VO2,结合VO2的相变状态和单元旋转实现太赫兹波段的数字编码[20];也有研究者仅使用VO2进行贴片图案设计,此类设计多需要不同形状或尺寸的VO2贴片,通过调节各贴片的性质进行组合,实现阵列编码[21,22].以上提到的编码单元各自分立,VO2结构部分无法分别激励,实验难以精准控制超表面单元的编码状态,实际可操作性较差.天津大学李佳辉等[20]在VO2薄膜上设计金属图形,能够实现整体调控,但是VO2薄膜的设计使电磁波经过该超表面时仅反射两种波束形态,限制了超表面的可编码特性.Kim 等[23]在条状金属纳米线下方设计一定厚度相同尺寸的VO2纳米线形成栅状结构,通过不同数字组合在理论上验证了仅使用一个结构单元就可以实现的1-bit 数字光学超表面.

现阶段报道的基于VO2的太赫兹编码超表面大多难以应用到实际中,其设计多以正向设计为导向.近年来,机器学习算法已成功应用于多个波段的电磁学研究与设计中,通过输入几个特征值或谱线能够对波导天线[24]、微带天线[25]、超表面[26]等器件进行结构优化和逆向设计.对于工作频率为25.5 GHz 的反射阵列天线,Prado 等[27]用支持向量机法加快其设计效率,该方法能够有效预测波束赋形阵列且与仿真结果高度一致;Hou 等[28]基于深度学习的目标驱动方法定向设计了GHz 超材料吸波结构;Ma 等[29]通过两个双向神经网络的部分堆叠,可以自动设计和优化具有强手性响应的三维太赫兹手性超材料.此类机器学习方法以结果为导向,我们也将深度神经网络算法应用到反射阵列天线设计中,该算法能够简化设计过程,缩短计算时间,提高电磁天线的设计效率.

为此,本文结合VO2提出了一种仅用两个单元结构就可实现2-bit (四种状态) 编码的方法,彼此连接的金属条使单行控制状态成为可能;结合单元的可调谐特性,预设波束偏转角得出阵列状态序列,并且通过MATLAB 软件进行了可视化系统设计;结合机器学习算法实现了阵列单元的逆向设计,提出的波束可调编码超表面阵列以结果为导向,增强波束调节的主动性与灵活性,为可调太赫兹天线的设计提供新的思路和方法.

2 理论计算与仿真

2.1 编码单元结构设计

在石英基底上设计的两种嵌入VO2贴片的天线单元结构见图1(a)和图1(d),为方便后续描述,图1(a)展示的“G”字结构标记为单元A,图1(d)展示的“艹”字型结构标记为单元B.它们都由三层结构组成,最上层是电导率为 4.56×107S/m 的金属金(黄色贴片)和VO2(红色贴片)组成的图案层,为了方便成行或成列地控制材料状态,每个单元都有一条贯穿的金属引线.中间层是厚度为500 µm的石英片,其介电常数为3.75;石英的另一侧,即结构最底层覆盖有一层金属金,具体尺寸参数如表1 和表2 所列.本次研究中涉及到太赫兹波段VO2的电磁参数,在没有外部激励的情况下,VO2是一个对太赫兹波透明的绝缘体,外界温度或激光激励使VO2发生相变后,电导率急剧增加,数量级可增加到105,完全相变为金属态的VO2介电常数εm(ω)满足Drude 模型[30],

表1 单元A 结构尺寸参数Table 1.Geometric parameters of unit A.

表2 单元B 结构尺寸参数Table 2.Geometric parameters of unit B.

其中,ωp是等离子体频率且满足=σ/(ε0τ),松弛时间τ为2.27 fs,εi指无限频率处的材料介电常数,ε0则是真空介质的介电常数,σ为电导率.在仿真中,仅使用VO2的完全绝缘态和完全金属态来实现编码单元的设计,所以在仿真软件中设置通用材料的属性为Normal 类,此时材料的介电损耗和电导率可以在页面中直接定义,根据实验验证[31],设置电导率分别为200 S/m 和200000 S/m[16,17].

图1(b)和图1(c)分别是单元A 的反射幅度和相位曲线,红色实线和蓝色虚线分别表示VO2是金属态和绝缘态时单元结构的电磁响应.在0.22 THz 处,可以观察到单元A 在VO2呈现绝缘态和金属态时反射相位差约180°,其反射幅度均大于0.7.图1(e)和图1(f)分别展示了以相同参考平面模拟观测单元B 的反射幅度和相位曲线,在0.22 THz 处单元B 在VO2呈现绝缘态和金属态时反射相位差约180°,其反射幅度均高于0.8.

图1 编码单元A (a)和B (d)示意图(黄色为金属,红色为VO2);VO2 为金属态和绝缘态时编码单元A 和B 的(b),(e)反射幅度图及(c),(f)反射相位图;在工作频率0.22 THz 处,编码单元A 在(g)金属态和(h)绝缘态时电场分布;编码单元B 在(i)绝缘态和(j)金属态时电场分布Fig.1.Model of unit A (a) and B (d) (The yellow part is the metal and the red part is VO2).Magnitude (b),(e) and phase (c),(f)of reflection for unit A and B when VO2 is in metallic state and insulation state.At operating frequency of 0.22 THz,the electric field distribution of unit A in the (g) metallic and (h) insulating state;the electric field distribution of unit B in the (i) metallic and(j) insulating state.

从电场分布的角度阐明其物理机理.对于单元A,中间的VO2贴片呈现金属态时,与引线相连的结构等效为一个方形纯金属贴片,该类结构在反射谱图中仅产生一个谐振,入射电磁波的激励集中在贴片的顶点与引线上,如图1(g)所示.图1(h)给出了VO2贴片呈现绝缘态时的电场分布,该结构相当于一个挂在引线下的开口环,由于等效开口环的宽度不一致,即图中尺寸g≠w,反射曲线会出现多个谐振点,在0.22 THz 时电磁波被束缚在开口处.对于单元B,VO2连接了中间引线与外围金属框,当它呈绝缘态时,外围金属框与引线之间出现缝隙,反射响应则受中间两个条形谐振器和缝隙影响,如图1(i)所示,电场能量也集中在条形谐振器两头和绝缘态VO2产生的缝隙处.当VO2贴片被电压加热呈现金属态时,顶层图案相当于金属引线、方框与条形谐振器的组合,可观察到工作频率处于两个谐振点之间,此时反射幅度较高,图1(j)中的电场分布则说明此时电磁波与垂直于电场方向的边框发生耦合.因此,当VO2贴片呈现不同相态时,表面的结构变化引起电磁场分布的变化,从而产生不同幅度与不同相位.

由图1 的反射曲线可观察到,在固定参考平面的前提下,单元A 与单元B 的两种状态彼此可以构成梯度90°的相位排列,如图1(c)和图1(f)所示,将单元A 的金属态及绝缘态分别标记为状态“1”(图1(g))和状态“3” (图1(h)),单元B 的绝缘态及金属态分别标记为状态“2”(图1(i))和状态“4”(图1(j)).此时,这4 个状态单元的反射幅度都很高,相位形成90°梯度差,可以覆盖一个完整周期内的360°相位,满足2-bit 编码的前提.

此项设计的突出亮点在于仅用两个单元结构就可实现2-bit 编码,其中单元A 可以在状态“1”和“3”之间转换,单元B 可以在状态“2”和“4”之间转换,对于固定的结构排列,通过控制不同行或列的状态可以在一维空间内实现多种编码序列,具有可调谐特性.从实验的角度来说,整行或整列的金属线贯穿设计为电压加热调控VO2特性提供了便利性和可行性.

2.2 固定阵列单元的状态序列可视化设计

超表面阵列天线主要由馈源和反射面组成,当超表面阵列被简化表示为在一维方向沿直线等间距排列的N个天线单元时,电磁波以角度θin入射到平面阵列上,经过阵列单元的相位调控以出射角θout反射到空间中,阵列在x方向上第n个单元的连续相位可以表示为

其中-π＜＜π,d表示单元周期,λ为工作波长.对于反射型天线,当入射波束从天线阵正面沿法线方向入射,即θin=0,在第n个天线单元处,反射波产生的辐射场强表示为

其中,Eout(rn,t)和Ein(rn,t)分别表示第n个单元的出射和入射电磁波的辐射场,G(θin)和G(θout) 分别是入射激励源和反射单元的增益,Rn和是第n个单元的反射和相移系数.根据上述理论分析,考虑到该反射阵列的性能参数与每一个阵列单元有关,在阵列中用给出第n个单元的反射幅度,当馈源发射的入射波经过超表面单元后,继续向均匀媒质传播,此时的电磁场方程是线性方程,满足叠加定理的条件.值得注意的是,构成反射阵列的单元仅能实现4 个相位,(2)式得到的相位值无法完全满足,因此在数值计算过程中,相位值可按照如下规则进行量化并用表示:

在远场观察的总场强可以看成N个单元在远场叠加之和,因此阵列的远场辐射功率P(θ) 为所有单元在任意方向上的辐射场功率之和:

综上所述,当改变线性阵列每一个超表面单元的幅度和相位时,远场波束的功率分布可以被改变,这也是用超表面结构实现波束控制的基础[32,33].

图2 给出了该类超表面天线阵列的示意图,在后文中,由A 和B 构成的是单元排布序列,由1,2,3 和4 构成的是状态排布序列.图2 中预先设定阵列在一维方向上按照“ABAB······”间隔周期排列,为了说明阵列的大小,定义一组阵列单元为超级单元(supercell),超级单元重复排列构成完整阵列.图2 中示例一个超级单元包含6 个阵列单元,4 个超级单元依次排列则可构成完整阵列,此时在一维方向上阵列共有24 个单元.

图2 超表面天线阵列示意图Fig.2.Schematic diagram of the metasurface antenna array.

根据(6)式将阵列天线的计算程序封装为一个可视化界面,输入部分包含工作频率(working frequency)、辐射图输出形式 (far-field)、期望偏转角(deflection angle)、单元周期 (unit length of each elements)、一个超级单元里包含的单元个数(array elements of supercell)和超级单元的个数(number of supercell) 6 个参数,输出部分包括实际辐射图和状态编码序列(coding sequence).首先在程序代码中设定阵列的单元排布序列,用户在文本框内输入上述6 个参数,程序通过量化计算可以给出该种单元排布下,实现预设偏转角的状态序列,并且给出这种状态序列排布下编码超表面的实际辐射图,该系统不限定单元结构,对一切满足编码条件的阵列具有普适性.

下面对编码超表面阵列天线进行计算和验证,在程序代码中设定单元按照“ABAB······”间隔周期排列,计算中预设天线阵列的工作频率为0.22 THz,每个单元周期为300 µm,一个超级单元由6 个单元构成,阵列中超级单元的个数为4.在统一工作频率和单元阵列的前提下,图3 给出了预设不同偏转角得到的状态序列和波束辐射图,其中状态序列编号与图1 一致.当预设偏转角为15°时,根据量化程序可得超表面状态序列为“121232323434341 212121232”,图3(a)和图3(b)分别展示了此种状态序列的极坐标辐射图和直角坐标辐射图.根据程序计算,实际主瓣波束偏转角度指向15.3°,几乎没有旁瓣干扰.当阵列预设偏转角为30°时,量化程序得到的阵列状态序列为“1232341212123434121 23234”,图3(c)和图3(d)分别给出了此种状态序列的极坐标辐射图和直角坐标辐射图.此时实际主瓣波束偏转角度指向30.3°,辐射图中有一束很小的旁瓣,因此主瓣的能量稍有减弱.另外,根据模拟分析,单元序列按照“ABAB…”构成反射阵列时,反射波束可在—65.3°—66.8°之间实现角度可调.

另一方面,图4 采用电磁仿真软件CST Studio Suite 中的快捷阵列计算模块对上述两组序列进行仿真.用户在模拟中仅对少数单元进行辐射图计算,并设置三个方向上的空间位移和相位差形成阵列对辐射情况进行模拟,这样大幅缩短了仿真时间.当阵列的反射幅度和相位按照序列“121232323 434341212121232”依次排布时,得到的辐射结果如图4(a)所示,其主瓣波束偏转角为15°;同理,当阵列的相位和反射幅度按照序列“1232341212123 43412123234”进行排列时,阵列主瓣偏转角为30°,这与图3 反映的情况基本一致.

图3 单元阵列按照“ABAB······”周期排列,预设角度不同时对应的状态序列的(a),(c)极坐标辐射图和(b),(d)直角坐标辐射图(a),(b) 预设角度为15°;(c),(d)预设角度为30°.图中标注出了实际主瓣波束偏转角度Fig.3.Under the condition of different deflection angle,(a),(c) polar radiation map and (b),(d) cartesian radiation map of corresponding coding sequence when the structure arrays are arranged alternately by “ABAB······” : (a),(b) The deflection angle is 15°;(c),(d) the deflection angle is 30°.The actual deflection angle of main lobe beam is marked in the figures.

图4 CST 仿真辐射结果 (a) 状态序列为“121232323434341 212121232”,主瓣波束偏转角为15°;(b) 状态序列为“123234 121212343412123234”,主瓣波束偏转角为30°Fig.4.CST simulation radiation: (a) State sequence of“121232323434341212121232”,the deflection angle of main lobe beam is 15°;(b) state sequence of “12323412121234341 2123234”,the deflection angle of main lobe beam is 30°.

本文提出的2-bit 编码超表面在固定单元排布的情况下能够实现反射波束可调,通过量化选择结构状态可以预测阵列的状态序列,同时程序也能对实际辐射图作出预判断.此外,该界面化的程序更方便用户操作.为了进一步拓展超表面阵列的研究维度,突破固定结构的限制,我们采用机器学习的方法对单元排布序列进行反向设计和研究.

2.3 基于机器学习的单元排布序列反向设计

通过以上的模拟对比可以发现,量化后序列的实际辐射角与预设角度不完全一致,这是因为程序计算过程中对理想相位排列进行了量化,相位排布不再连续,使实际计算得到的辐射图与预设值有一定误差.为更好比较不同状态序列的辐射效果,定义能量占比Eratio和角度误差Rθ如下:

其中,Emain和Eall分别表示实际和理想状态辐射主瓣的半高全宽功率,θreal和θpreset分别是实际和预设的偏转角度.

神经网络在过去的二十年中已经在很多领域实现应用,一般来说,神经网络由输入层、隐藏层和输出层组成.每个隐藏层中层和神经元的数量越大,模型的复杂性就越大.当隐藏层的数量和神经元数量增加时,该神经网络构成一个深度神经网络(deep neural network,DNN)[34-36].将DNN 引入天线阵列的设计,以获得更准确的超表面单元排布序列.图5(a)给出了设计中所用的深度神经网络示意图,输入预期的波束偏转角、能量占比和角度误差,经过算法计算,得到对应单元排布序列.在训练过程中,首先对单元结构周期排布的多种阵列进行预设偏转角为30°,40°,50°,60°和70°的模拟,得到共560 组能量占比和角度误差值,并将其作为训练数据集,把预设偏转角为35°,45°,55°和65°的数据组作为测试数据集.本次训练中,网络结构为3 层隐藏层,每层64 个节点,并且用误差指标均方误差(MSE)来评估该网络,计算公式为

训练过程中的损耗曲线如图5(b)所示,可见,随着训练次数的增加,误差得到明显改善,经过50 次训练后收敛到一个常数值,训练完成时,MSE 为0.3181.

图5 (a) 深度神经网络示意图;(b) DNN 预测误差曲线Fig.5.(a) Structure of proposed DNN;(b) the training loss of DNN.

选择测试集中的3 组数据[0.87229362,35,0],[0.72561196,45,0]和[0.83603583,55,0]进行MATLAB 计算和CST 阵列仿真,每组数据中的数字分别表示能量占比、偏转角度和角度误差.数据组[0.87229362,35,0]的预测单元序列为AABBAAABAAABAABBAAABAAABAAAA,图6(a)和图6(b)是对应MATLAB 和CST 阵列的模拟结果,可以看出此时波束指向35°,与预设结果完全一致,但主瓣能量占比只有63.17%.数据组[0.72561196,45,0]的算法预测单元序列为ABBAAABBABAAABBBABAAAABB.从图6(c)和图6(d)可以看出MATLAB 计算和CST 模拟波束分别指向45.7°和45.0°,与预设结果高度一致,主瓣能量仅为57.7%.图6(e)和图6(f)则是数据组[0.83603583,55,0]预测单元结构序列AAAAAABBAAAAAABBAAAAAAAB 的仿真结果,MATLAB 和CST 计算实际偏转角度为55.8°和52.0°,与预设角度存在一定误差,其主瓣能量占比为60.19%.

图6 MATLAB (a),(c),(e)和CST (b),(d),(f)对不同偏转角的辐射验证 (a),(b) 35°;(c),(d) 45°;(e),(f) 55°.图中标注出了实际偏转角度Fig.6.Verification with different deflection angles by MATLAB (a),(c),(e) and CST (b),(d),(f): (a),(b) 35°;(c),(d) 45°;(e),(f) 55°.The actual deflection angle is marked in the figures.

通过对单元阵列的预测结果进行验证,发现该算法得到的实际偏转角度与预测基本一致,然而辐射波束的能量占比与预计结果有一些差距,这是因为此次实验中为缩短计算时间只进行了有限次训练,但是上述结果在一定程度上说明深度学习算法在阵列反向设计的有效性.综上所述,将DNN 算法运用在太赫兹阵列反向设计中能够有效预测单元阵列的排布形式,使反射波束以期望角度出射.

3 结论

本文提出了一种仅用两个单元结构就可实现2-bit 编码超表面的方法,通过在单元结构中引入可调谐的VO2贴片,在一维方向上实现了反射波束的角度可调.在此基础上,通过MATLAB 对编码超表面阵列天线计算进行可视化封装,使计算操作更为简易,该系统提前设定单元排布序列,通过预设波束偏转角反向计算状态序列,使固定阵列能够实现反射波束角度可调谐.最后,将超表面阵列排布与深度神经网络相结合,实现了阵列单元的逆向设计,MATLAB 和CST Studio Suite 仿真结果验证了所用算法在波束偏转角度和单元排布方面的有效性.