-
本文设计了由四聚长方体组成的全介质超表面, 其中每个长方体刻蚀两个椭圆柱并填装空气. 当分别为超表面单独引入面内对称破缺、位移扰动和周期扰动时, 可在近红外波段产生稳健的准连续域束缚态模式(quasi-bound states in the continuum). 通过测量准BIC (quasi-BIC)模式的谐振波长, 计算准BIC模式的Q因子(quality factor)与不对称参数的关系, 可进一步证实不对称参数对准BIC共振频率和Q因子的可调谐性. 在此基础上, 当同时引入面内对称破缺、位移扰动和周期扰动时, 可获得5个高Q因子的准BIC模式. 共振峰的数量、位置以及Q因子都可通过调整面内破缺、位移扰动和周期扰动的程度进行调控. 该超表面的设计可为传感器的多参数传感以及灵敏度等性能的提升提供一种全新思路.High-quality factor (high-Q) resonance has broad prospects in applications such as in narrow-band filtering, slow-light devices, and nonlinear optics interaction enhanced to highly sensitive sensing. Previous methods of designing high-Q resonance suffered intrinsic drawbacks such as high-volume cavities or large-scale bending radii. However, recently, a new approach to designing high-Q resonances has begun to attract public attention on the basis of asymmetric metasurfaces that are related to the bound states in the continuum (BIC) phenomenon. Constructing BIC resonance in electromagnetic metasurface can generate sharp resonant transmission peak. Therefore, there is growing interest in utilizing BIC to achieve metasurface with high-Q. However, most of existing studies are based on single BIC, and few studies focusing on multiband BICs and multiple forms of symmetry breaking. In this work, we propose an all-dielectric metasurface composed of tetrameric cuboids. By etching two elliptical cylinders in each cuboid, the metasurface can simultaneously support in-plane symmetry breaking, displacement perturbations and periodic perturbations. We first use multipole calculations to analyze the physical mechanism by which the metasurface generates quasi-BIC under these three conditions. It is confirmed that the Q factor and resonant peak position of quasi-BIC can be controlled by adjusting the asymmetry parameters. Subsequently, we introduce the in-plane symmetry breaking, displacement perturbations and periodic perturbations into the metasurface simultaneously and generate five quasi-BIC modes, whose numbers and positions can be flexibly adjusted, and the largest Q factor is 58039. In summary, this work provides a new practical design concept for realizing high-Q all-dielectric metasurfaces, which can be used to improve the sensitivity of multi-parameter sensors.
-
Keywords:
- metasurface /
- quasi-bound states in the continuum /
- symmetry breaking
引 言
空腔外形广泛应用于航空飞行器中, 如飞机起落架、武器舱等. 空腔流动包含剪切层、旋涡、复杂波系等典型流场结构, 在超声速流动中还包含激波. 这些复杂波系之间以及与空腔壁面的相互作用会产生自持振荡而引起强烈的噪声, 噪声声压级可达到170 dB以上[1], 强噪声不仅会导致环境污染, 还会引起飞行器结构疲劳, 影响飞行安全[2]. 准确认识空腔流动中复杂波系相互作用产生噪声的物理机理可为研究工程上可行的降噪方法提供理论基础. 为了掌握空腔内部流动结构和振荡频谱特性, 大量学者通过实验和数值模拟研究了不同类型的空腔流动, 相关研究工作可参考相关综述文献[3, 4].
空腔流动的分类与来流条件和空腔长宽比直接相关, 如超声速流动长深比在10以下为开式流动, 长深比大于13为闭式流动, 长深比介于二者之间为过渡式流动[4]. 开式与闭式流动的主要区别在于, 对于开式流动, 空腔上方的剪切层会流过腔口而撞击到空腔后壁; 而对于闭式流动, 剪切层会撞击空腔底部, 类似于后台阶流动[4]. 在工程实际中, 空腔的长深比大约在4—10之间, 表现为开式流动[5]. 开式空腔会因流动的自持振荡产生强噪声, 主要过程为[6]: 上游边界层脱离空腔前缘, 由于流体的Kelvin-Helmholtz不稳定性, 剪切层逐渐失稳并与空腔后壁面发生碰撞形成反射声波, 反射声波向上游传播到达空腔前缘进一步激发剪切层失稳, 形成反馈回路. 这种声反馈回路机制最早可见于Powell[7]关于边缘音的研究中. 通过一系列实验, Krishnamurty[8]发现亚声速和超声速空腔流动会产生强烈的具有离散频率的噪声. Rossiter[6]进一步拓展了实验参数范围并给出了预测空腔振荡离散频率的半经验公式, 该公式推导基于声速为常数, 因而在马赫数较高时(M > 2)与实验测量偏差明显. Heller等[9]考虑空腔中的声速变化, 给出了一个改进的半经验预测公式. 部分学者试图给出无经验参数的空腔振荡离散频率预测公式[10-13], 其中Kerschen等[12,13]通过将剪切层、空腔中向上下游传播的波系以及空腔两端的散射情况考虑在内提出了一个理论模型, 该模型不包含经验参数且可以预测共振模态的增长与衰减情况.
上述半经验公式和模型基于二维波系结构建立, 然而在实际空腔流动中, 空腔还存在展向流动, 此外空腔侧壁也会对腔内振荡产生影响[14,15]. Kegerise等[16]在长深比为2、马赫数为0.4的空腔流动实验中, 观察到腔内压力信号频谱中存在一个未知的低频成分, 通过小波分析发现, 该低频成分会对Rossiter主模态产生调制作用. Brès和Colonius [17]基于周期性的展向边界条件数值模拟了三维空腔流动, 他们在速度脉动中同样发现了低频成分, 并认为这种低频成分为三维模态; 同时发现三维模态与Rossiter振荡模态可同时存在, 二者间的相互作用可能会影响空腔流动自持振荡主模态的阶数; 三维模态的波长与空腔深度接近, 比Rossiter振荡模态在幅值上小一个量级. Brès和Colonius [17]认为, 三维模态的产生机制与空腔下游壁面附近再循环涡流中的闭合涡流相关的离心不稳定性有关. 因此, 理论模型还须考虑三维模态的影响.
空腔流动产生的噪声在频谱上可能包含多个具有离散频率的声模态, 设计高效空腔噪声控制方案需要厘清这些模态的演化规律. 对于二维空腔流动, Rowley等[18]采用直接数值模拟研究了亚声速情形, 结果表明一阶、二阶Rossiter模态是共存的. 我们前期针对二维空腔流动的直接数值结果表明[19—21], 在亚声速情况下, 低阶模态会随着流动稳定逐渐消失; 在超声速情况下, 主模态会产生谐振模态, 主模态与其是共存的. 上述模态演化规律与二维层流这种特定条件有关, 在三维流动中, 声模态的演化特征与上述情形有明显区别. 1998年, Cattafesta等[22]开展的马赫数0.6的空腔流动实验结果表明, 主要的Rossiter模态之间存在模态切换现象. Gloerfelt等[23]采用LES数值模拟了展向为周期边界、马赫数0.8的三维空腔流动, 在流动结构中发现了模态切换的证据. Thangamani[24]实验研究了马赫数1.58的三维超声速空腔流动, 同样发现一阶和二阶主模态之间存在模态切换现象. 我们前期针对三维空腔流动的数值模拟结果同样表明[25], 空腔振荡的主导模态之间存在切换现象. 但上述研究对模态切换性质的认识仍有不足, 对于空腔流动中各阶模态在流动振荡过程中的非线性相互作用仍需进一步研究.
受Kerschen等[12,13]工作的启发, 本文通过研究空腔中包含展向波系在内的主要波系结构, 建立相关波系模型. 各部分内容如下: 在第1节中, 区别于以往的二维模型[6,9], 通过在空腔两端引入展向波系, 同时考虑亚声速与超声速的区别, 分析空腔两端波系的散射过程, 对亚声速与超声速两种情况下的三维空腔波系分别进行建模, 同时分析波系之间的耦合作用; 在第2节中, 介绍马赫数0.9和马赫数1.5空腔流动的实验情况; 在第3节中, 基于马赫数0.9和马赫数1.5空腔流动实验数据, 采用信号分析工具, 结合所建模型分析空腔振荡模态的非线性耦合作用, 同时研究空腔流动中各阶模态的演化规律; 第4节为结论部分.
1. 散射波系模型
1.1 波系模型
1.1.1 亚声速情形
在亚声速情况下, 三维空腔流动中的波系结构示意图如图1所示. 空腔的长、宽、深分别设为L, W, D. 空腔内部包含向上游传播的声波
wb , 其主要是由剪切层撞击后壁导致; 空腔内部向下游传播的声波wf , 其主要由wb 传播至空腔前端时由前壁面反射形成; 剪切层中的不稳定波ws , 这是由于流体介质的Kelvin-Helmholtz不稳定性导致的; 空腔口沿剪切层向下游传播的声波wd 以及向上游传播的声波wu ,wd 的产生与向上游传播的波与空腔前缘相互作用有关,wu 与剪切层撞击后壁过程直接关联;wt 表示流动在展向上的波动. 下文中, 在空腔前端的上述波系及相关参数用^[⋅] 表示, 以区别于在空腔后端的波系及相关参数.参考Kerschen等[12,13]的做法, 接下来考虑上述波系在空腔两端的散射情况. 在空腔前端, 考虑波系之间的相互关系有:
{ˆws=ˆCsbˆwb+ˆCsuˆwu+ˆCstˆwt,ˆwf=ˆCfbˆwb+ˆCfuˆwu+ˆCftˆwt,ˆwd=ˆCdbˆwb+ˆCduˆwu+ˆCdtˆwt. (1) 其中,
ˆC[⋅] 表示波系之间的散射系数. 类似地, 在空腔后端同样有:{wu=Cusws+Cudwd+Cufwf+Cutwt,wb=Cbsws+Cbdwd+Cbfwf+Cbtwt. (2) 对于空腔前后同一个波, 根据波的增长关系, 可以将波的发展写成复数形式, 对于空腔内的波和剪切层不稳定波有:
ˆwb=wbeiαbL,ws=ˆwseiαsL,wf=ˆwfeiαfL. (3) 对于三维展向波有:
ˆwt=ˆσteiˆαtW,wt=σteiαtW. (4) 对于空腔外部的波, 其沿空腔流向以
l−3/2 衰减[12, 13], l表示离空腔前缘的距离, 因此其形式为ˆwu=wueiαuL/L3/2,wd=ˆwdeiαdL/L3/2. (5) 定义变量
X 为X=(ˆws,ˆwf,ˆwd,wb,wu)T. (6) 将(3)式—(5)式代入(1)式和(2)式并整理, 即可得到矩阵形式的方程组:
AX=b , 其中,A=(−100ˆCsbeiαbLˆCsueiαuL/L3/20−10ˆCfbeiαbLˆCfueiαuL/L3/200−1ˆCdbeiαbLˆCdueiαuL/L3/2CbseiαsLCbfeiαfLCbdeiαdL/L3/2−10CuseiαsLCufeiαfLCudeiαdL/L3/20−1), (7) b=−(ˆCstˆσteiˆαtW,ˆCftˆσteiˆαtW,ˆCdtˆσteiˆαtW,CbtσteiαtW,CutσteiαtW)T. (8) 于是,
|A|X=A∗b , 其中|A| 表示矩阵A 的行列式,A∗ 是矩阵A 的伴随矩阵. 当忽略展向波时, 上述波系为二维情形, 此时b=0 , 此时方程退化为Kerschen等[12]描述的情形.显然, 上述方程组解的属性与
|A| 的性质密切相关, 其表达式为|A|=−1+CbsˆCsbei(αb + αs)L+CbfˆCfbei(αb + αf)L+L−3/2CbdˆCdbei(αb + αd)L+L−3/2CusˆCsuei(αu + αs)L−L−3/2CufˆCfuei(αf + αu)L+L−3CudˆCduei(αu + αd)L+Γ4, (9) 其中,
Γ4 表示行列式值的剩余的四阶组合项:Γ4=L−3ei(αb+αd+αf+αu)L (−CbfCudˆCduˆCfb+CbdCufˆCduˆCfb+ CbfCudˆCdbˆCfu−CbdCufˆCdbˆCfu)+L−3ei(αb+αd+αs+αu)L (−CbsCudˆCduˆCsb+CbdCusˆCduˆCsb+ CbsCudˆCusˆCsu−CbdCusˆCdbˆCsu)+L−3/2ei(αb+αf+αs+αu)L (−CbsCufˆCfuˆCsb+CbfCusˆCfuˆCsb+ CbsCufˆCfbˆCsu−CbfCusˆCfbˆCsu). (10) 1.1.2 超声速情形
对于超声速空腔流动, 不考虑激波结构, 其波系结构示意图如图2所示. 在超声速情况下, 空腔中的波系与亚声速情况略有不同, 在空腔上方, 由于流动是超声速的, 声波不能向上游传播, 但向下游传播的声波包含慢声波和快声波, 这里分别记为
wds 和wdf . 整体波系结构中, 只有空腔中声波wb 向上游传播.类似于亚声速情况, 分别考虑空腔前后两端的波系散射情况. 在空腔前端有:
{ˆws=ˆCsbˆwb+ˆCstˆwt,ˆwf=ˆCfbˆwb+ˆCftˆwt,ˆwds=ˆCdsbˆwb+ˆCdstˆwt,ˆwdf=ˆCdfbˆwb+ˆCdftˆwt. (11) 在空腔后端有:
wb=Cbsws+Cbfwf+Cbdswds+Cbdfwdf+Cbtwt, (12) 其中, 慢、快声波的演化衰减关系为[12]
wds=ˆwdseiαdsL/L3/2,wdf=ˆwdfeiαdfL/L3/2. (13) 定义变量
X 为X=(ˆws,ˆwf,ˆwds,ˆwdf,wb)T, (14) 将(3)式、(4)式和(13)式代入(11)式和(12)式并整理, 即可得到矩阵形式的方程组:
AX=b , 其中,A=(−1000ˆCsbeiαbL0−100ˆCfbeiαbL00−10ˆCdsbeiαbL000−1ˆCdfbeiαbLCbseiαsLCbfeiαfLCbdseiαdsLL3/2CbdfeiαdfLL3/2−1), (15) b=−(ˆCstˆσteiˆαtW,ˆCftˆσteiˆαtW,ˆCdstˆσteiˆαtW,ˆCdftˆσteiˆαtW,CbtσteiˆαtW)T, (16) 行列式值
|A| 为|A|=−1+CbsˆCsbei(αb + αs)L+CbfˆCfbei(αb + αf)L+L−3/2CbdsˆCdsbei(αb + αds)L+L−3/2CbdfˆCdfbei(αb + αdf)L. (17) 不难发现, 由于超声速情形在空腔两端的散射特征较亚声速情形简单, 因此, 方程组的系数矩阵和对应的行列式值也更为简洁.
1.2 模型分析
在Rossiter[6]得到的半经验公式中, 图1和图2所示的波系只考虑了
ws 和wb , 方程组变为AX=0 . 其中, 矩阵A 简化为A=(−1ˆCsbeiαbLCbseiα sL−1). (18) 因为方程组有非零解, 于是必有
|A|=0 , 即|A|=−1+CbsˆCsbei(αb + αs)L=0. (19) 利用复数欧拉公式, 可以得到:
Arg[Cbs]+Arg[ˆCsb]+Re[αb+αs]L=2nπ, (20) 其中, Arg[▪]表示复数的幅角, Re[▪]表示复数的实部. (20)式可以进一步表示为
L=2πm−Arg[Cbs]+Arg[ˆCsb]2πRe[αb+αs], (21) 其中
αs 是剪切层不稳定波的波数, 设剪切层的对流速度为Uc , 即有αs=ω/Uc .αb 是空腔后壁面向上游传播的声波的波数, 其传播速度近似为声速, 因此,αb=ω/a=ω⋅M/U . 将αs 和αb 代入(21)式并化简可以得到:Stn=ωL2πU=fnLU=n−ˆγM+1/ˆk,n=1,2,⋯, (22) 其中,
{ˆγ=(Arg[Cbs]+Arg[ˆCsb])/(2π),ˆk=Uc/U. (23) (22)式和(23)式是Kerschen等[12]通过空腔两端的散射模型获得的空腔振荡频率预测公式, 可以看出该公式与Rossiter半经验公式[6]形式上完全一致. 在(23)式中, 函数
ˆγ 和ˆk 与具体的空腔尺寸和来流条件有关, 在实际中一般变化较小, 因此可以近似取为常数. 在Rossiter半经验公式中, 一般取ˆγ=0.25 ,ˆk 为剪切层相对对流速度, 一般取ˆk=0.57 . 在我们的二维层流空腔流动高精度数值模拟中[19-21], 获得的ˆk=0.55 , 与常用值较为接近. 在实际三维空腔流动中, 由于湍流、侧壁等的影响, 其值可能会有所区别, 可以根据实际测量值进行估计, 如在Ahuja和Mendoza[14]以及Kegerise[26]的实验中, 得到的估计值为ˆk=0.66 .当进一步考虑三维展向波系时, 右端项
b 简化为:b=−(ˆCstˆσteiˆαtW,CbtσteiαtW)T. (24) 即有:
X=(ˆCbtˆCsbσtei(Lαb+Wαt)+ˆCstˆσteiWˆσt−1+CbsˆCsbeiL(αb+αs)ei(Lαs+Wˆσt)Cbsxstˆσt+eiWαtˆCbtσt−1+eiL(αb+αs)CbsˆCsb). (25) 此时方程组
AX=b 决定了三维展向波系将与Rossiter模态产生调制. 尽管三维展向模态频率和幅值较低, 但在特定情况下, 三维模态仍然会与主模态产生非线性作用. 在Neary和Stephanoff[27]基于水介质的层流空腔实验中, 他们发现了三维模态与主模态调制产生的谐频. 在Kegerise等[16]马赫数为0.4的空腔流动实验中, 三维模态对Rossiter主模态产生了调制作用. 同时, 空腔宽度参数W 也会影响不同离散频率组分的幅值[3].当考虑空腔中的其余波系时, 方程组
AX=b 的形式变得较为复杂, 矩阵A 的行列式值如(9)式和(17)式所示. 此时, 不同波系之间将产生非线性的相互作用, 尽管部分散射过程波系相互作用的系数值较小, 但仍将产生相互调制作用, 这种调制过程可能产生其余频率成分. Kegerise等[16]通过对空腔中测得的压力信号采用高阶谱分析发现, 在特定的流动条件下, 不同Rossiter模态之间会产生非线性作用, 这种非线性作用将导致产生不同于Rossiter模态的谐频成分, 频率大小与两个Rossiter模态频率的差值相关.2. 三维空腔实验
2.1 实验设置
空腔模型C201由中国空气动力研究与发展中心高速空气动力研究所设计[28], 模型示意图如图3所示. 空腔的长(L)、宽(W)、深(D)分别为200, 66.7, 33.3 mm, 即满足L∶W∶D = 6∶2∶1, 模型的其余尺寸见图3标注. 在实验中, 亚声速和超声速的前缘设计略有区别, 具体可参考引用文献中的说明[28].
实验在中国空气动力研究与发展中心高速空气动力研究所亚跨超声速风洞中完成, 实验包含亚声速和超声速两种来流条件, 马赫数分别为0.9和1.5, 其余参数条件如表1所列.
表 1 不同马赫数下的实验参数Table 1. The experimental parameters for different casesCase Mach number Pressure/Pa Temperature/K Reynold number(ReL) 1 0.9 67938.3 247.85 3.41 × 106 2 1.5 38080.2 198.62 4.28 × 106 2.2 实验验证
实验中的脉动压力测量点设置在空腔底部中截面上, 下面通过重复性实验, 根据壁面压力信号总声压级验证实验测量的可靠性. 总声压级(overall sound pressure level, OASPL)计算公式如下:
OASPL=20log10(pePref), (26) 其中,
pe=√1NN∑i=1(p−ˉp)2, (27) ˉp 为压力平均值, 空气中参考声压一般取pref=2.0×10−5Pa , 对应正常人耳对 1 kHz声音的可听阈声压. 对于马赫数1.5, 两次实验的对比结果如图4所示, 两次测压最大差别为1.1 dB, 表明实验重复性较好, 实验系统可靠.3. 模态演化特性分析
3.1 频谱调制特性
声压级的频域分布通过压力脉动的声功率谱密度(pressure spectral density, PSD)进行计算, 公式如下,
SPL=10log10(PSDP2ref). (28) 对于当前两个马赫数下, 空腔前壁测点的压力信号声压级频谱计算结果如图5所示. 两种马赫数的主导频率如表2所列. 在(23)式中, 参数
ˆγ 和ˆk 与具体的空腔尺寸和流动条件相关, 这里可以通过得到的表2中的频谱数据拟合得到. 对于马赫数0.9, 通过线性拟合, 可以得到ˆr=0.45 ,ˆk=0.68 ; 类似地, 对于马赫数1.5, 可以得到ˆr=0.40 ,ˆk=0.86 . 将获得的参数代入(22)式即可得到离散频率预测公式, 公式预测值见图5中灰色竖线所示. 可以看到, 两种马赫数下, 空腔流动的振荡主导模态都为Rossiter二阶模态, 同时预测的各阶模态频率位置与压力信号频谱峰值符合良好.表 2 不同马赫数下的主导峰值频率Table 2. The frequencies of the most energetic peaks for different cases.Case f1 f2 f3 f4 f5 f6 f7 1 330 940 1500 2112 2735 3344 3925 2 480 1264 2042 2871 3673 4468 — 从图5可以观察到, 两个马赫数下的频谱中都存在一个低频成分, 对马赫数0.9, 这个低频频率为
f=134.3 Hz , 对马赫数1.5, 这个低频频率为f=122.2 Hz . 根据Brès和Colonius[17]的研究, 这种低频成分与三维模态有关, 其波长与空腔深度接近, 因此建议采用空腔深度作为其特征长度. 对于马赫数0.9, 低频成分的无量纲频率为StD=flD/U=0.016 , 对于马赫数1.5, 低频成分的无量纲频率为StD=flD/U=0.01 . 在Kegerise等[16]马赫数0.4的空腔流动实验中, 观察到的低频成分的无量纲频率为StD=flD/U=0.011 , 与当前两个马赫数下观测到的低频成分无量纲频率较为接近.在第2节的分析中, 空腔中出现的主要波系的主导频率可以通过方程组确定, 不同波系之间可能产生相互调制, 这种调制过程可能产生区别于主导Rossiter模态之外的其余频率成分. 这种现象可以在本文两个实验的频谱结果中观察到. 在图5(a)所示的马赫数0.9的频谱结果中, 在圆圈圈出的位置可以观察到一个显著的区别于主频的频率:
fs=3613 Hz . 根据表2中马赫数0.9对应的结果可以发现fs≈f3+f4=3612 Hz . 即频率fs 是由Rossiter三阶模态与Rossiter四阶模态调制产生的. 在图5(b)所示的马赫数1.5的频谱结果中, 在圆圈圈出的位置可以观察到有两个显著的区别于主频的频率fs1=2526 Hz 和fs2=3314 Hz . 根据表2中马赫数1.5对应的结果可以发现fs1≈f2+f2=2528 Hz ,fs2≈f2+f3=3306 Hz . 即频率fs1 是由Rossiter二阶模态与自身调制产生的, 频率fs2 是由Rossiter二阶模态与Rossiter三阶模态调制产生的. 这里仅列出了个别调制产生的幅值显著的谐频, 接下来采用高阶谱分析进一步分析这种调制作用.从上述分析可以看出, 对于马赫数1.5, 其Rossiter二阶主模态的幅值显著高于其余振荡模态, 并可与多个其余模态调制, 产生区别于主模态的谐频. 这里采用双谱分析(bispectral analysis)来分析调制过程. 对于一个离散时间信号
x(t) , 其对应的频域设为Xk , 即满足:x(t)=∞∑k=−∞Xke−iωkt (29) 和
Xk=limk→∞1T∫T/2−T/2eiωktx(t)dt, (30) 其中,
ωk=2πk/T ,T 是信号x(t) 的长度. 相应的三阶累计谱为B(k,l)=E[XkXlX∗k+l], (31) 其中, E[▪]表示期望值. 双谱分析相关系数bic定义为[29]
bic2(k,l)=|B(k,l)|2E[|XkXl|2]E[|Xk+l|2]. (32) 利用Cauchy-Schwarz不等式容易证明上述相关系数满足
0⩽bic⩽1 . 上述定义中的双谱相关系数具有对称性质. 特别是, 当频率ωk 和频率ωl 存在非线性耦合时,bic→1 ; 而当频率ωk 和频率ωl 相互独立时,bic→0 ; 当bic为介于0—1之间的其余数值时, 其值的大小表征频率ωk 和频率ωl 非线性作用产生频率ωk+l 占实际频率ωk+l 幅值的比例. 于是可以利用上述性质来分析空腔振荡中不同模态的非线性耦合情况.马赫数1.5空腔前面测量信号的双谱分析系数计算结果如图6所示, 图6(a)为整体结果, 图6(b)为局部区域放大情况, 结果关于斜中心线是对称的. 由于二阶Rossiter主模态较强, 这里主要关注与其相关的非线性作用情况. 从图6(a)可以看出, 与主频
f2 相关的非线性作用十字区域相较于其余区域十分明显. 在坐标(f2,f2) 处, bic = 0.94, 即主模态自相关性很强, 由此会产生相关频率2f2 , 这与上文中结论一致, 从该频率成分在图5(b)中的频谱结果可以看到, 幅值十分显著. 从图6(b)还可以看到, 在坐标(f2,f=1361 Hz) 处双谱相关系数较大: bic = 0.85, 意味着将产生频率为f2+f=2625 Hz 的成分, 该频率成分见图5(b)中的标注.3.2 模态演化规律
在3.1节中, 使用快速傅里叶变换(FFT)等获得了主要振荡模态的频谱特征, 但FFT缺乏时频定位功能, 且无法获得非平稳信号的频率演化信息. 从前文论述中可知, 在三维空腔流动中, 流动振荡的各阶模态之间存在模态切换现象, 因此在接下来的分析中, 采用时频分析方法来分析其压力脉动特征. 采用连续小波变换(continuous wavelet transform, CWT)来分析压力脉动信号, 小波变换相较于短时傅里叶变换等时频分析方法, 可以实现随分辨率变化而自动调节分析带宽, 且计算过程存储量小, 适合分析长时间尺度的信号, 其具体变换公式如下[30]:
Wx(a,b)=1√|a|∫Rx(t)ψ∗(t−ba)dt=∫Rx(t)ψ∗a,b(t)dt=⟨x(t),ψa,b(t)⟩, (33) 其中,
Ψ 为小波函数. 本文选取Morlet小波, 形式为Ψ(t)=e−t2/2eiΩ0t. (34) 小波变换计算结果如图7所示, 图中反映了压力脉动中不同频率成分随时间的演化情况. 在前面关于二维层流空腔流动高精度数值模拟的腔内压力信号的时频分析结果中显示, 主导模态随时间演化过程中幅值稳定, 当前的时频结果与二维简单情形有明显的不同. 根据表2中的数值, 从图7可以看出, 在两个马赫数下, 主要振荡模态的幅值随时间有明显的改变, 存在模态切换现象, 但Rossiter 2模态占主导时幅值要更高, 与前面的频谱分析结果一致. 为了进一步分析主要振荡模态随时间的演化特征, 以马赫数0.9为例, 将一阶模态
f1=330 Hz 和二阶模态f2=940 Hz 从小波分析结果中单独截取出来, 结果如图8所示. 需要说明的是, 从前面的理论模型分析中可知, 由于不同波系之间存在非线性作用, 主要振荡模态的频率在随时间演化过程中可能会发生微小的改变, 由于变化很小, 这里将其忽略. 从图8可以看出, 模态f1 和f2 的幅值随时间变化剧烈, 对截取结果进行FFT分析结果如图9所示. 从图9可以看出, 模态成分f1 和f2 的幅值随时间变化的波动频率较低, 但没有显著的主频. 进一步地, 使用概率密度估计来计算模态f1 和f2 随时间演化信息的概率密度函数分布, 结果如图10所示. 从图10可以看出, 两个模态随时间演化的概率密度分布形态与正态分布形态较为接近, 即模态演化整体表现出随机性.图 7 两个马赫数的空腔前壁测点压力信号的连续小波变换结果 (a) M = 0.9 (注: 等值线范围(200, 3700)); (b) M = 1.5 (注: 等值线范围(300, 6000))Fig. 7. The CWT results of the pressure perturbation signals at the front wall of the cavity for two Mach numbers: (a) M = 0.9 (contour levels between 200 to 3700); (b) M = 1.5 (contour levels between 300 to 6000).图 9 一阶、二阶模态小波变换系数的FFT频谱 (a) 一阶模态系数FFT; (b) 二阶模态系数FFTFig. 9. The FFT results of the CWT coefficients of the dominant modes in Fig. 8: (a) FFT of the coefficient of the first mode; (b) FFT of the coefficient of the second mode.图 10 一阶、二阶模态小波变换系数的PDFFig. 10. The probability density estimate of the CWT coefficients of the dominant modes in Fig. 8.4. 结 论
本文通过研究空腔中包含展向波系在内的主要波系结构, 分析不同波系在空腔两端的散射情况, 同时考虑亚声速和超声速流动下波系传播的差异, 分别建立了针对亚声速和超声速空腔流动的三维波系模型. 三维波系模型包含了空腔中不同波系之间的非线性相互作用, 这种非线性作用可导致产生不同于Rossiter模态的其余频率成分. 当仅保留剪切层不稳定波及其撞击后壁引起的反射声波时, 模型退化为二维简单情况, 此时方程表征空腔流动主要振荡模态—Rossiter模态.
基于三维空腔流动实验测量的马赫数0.9和马赫数1.5的腔内压力信号数据, 通过计算其声压级频谱, 获得了两种典型流动条件下, 空腔自持振荡的主要频谱特征. 根据获得的主要振荡频率, 对所建模型中的参数进行了线性估计, 模型预测的各阶模态频率位置与压力信号频谱峰值符合良好. 对压力信号采用双谱分析, 结果表明, 振荡主模态之间会产生非线性作用, 这种非线性作用产生了幅值较高的谐频, 与频谱分析结果一致. 采用连续小波变换方法对压力信号进行了时频分析, 结果表明, 主要的振荡模态之间存在模态切换现象, 且模态切换呈低频特征, 整体表现出随机性特性.
需要注意的是本文的空腔波系结构模型未考虑有舱门的情况, 舱门的存在会对空腔口的波系产生反射作用, 且舱门在不同角度时, 这种反射作用会有所差别. 针对上述情形, 仍需结合仿真和实验开展进一步研究.
[1] von Neumann J, Wigner E P 1929 Phys. Z. 30 467
[2] Kaelberer T, Fedotov V A, Papasimakis N, Tsai D P, Zheludev N I 2010 Science 330 1510
Google Scholar
[3] Plotnik Y, Peleg O, Dreisow F, Heinrich M, Nolte S, Szameit A, Segev M 2011 Phys. Rev. Lett. 107 183901
Google Scholar
[4] Hsu C W, Zhen B, Lee J, Chua S L, Johnson S G, Joannopoulos J D, Soljacic M 2013 Nature 499 188
Google Scholar
[5] Urseel F 1997 Phys. Eng. Sci. 435 575
[6] Porter R, Evans D V 2005 Wave Motion 43 29
Google Scholar
[7] Xiao Y X, Ma G, Zhang Z Q, Chan C T 2017 Phys. Rev. Lett. 118 166803
Google Scholar
[8] Retzler C H 2001 Appl. Ocean Res. 23 249
Google Scholar
[9] Cobelli P J, Pagneux V, Maurel A, Petitjeans P 2009 Europhys. Lett. 88 20006
Google Scholar
[10] Cobelli P J, Pagneux V, Maurel A, Petitjeans P 2011 J. Fluid Mech. 666 445
Google Scholar
[11] EVERY A G 2008 Phys. Rev. B 78 174104.
Google Scholar
[12] Tittl A, Leitis A, Liu M, Yesilkoy F, Choi D Y, Neshev D N, Kivshar Y S, Altug H 2018 Science 360 1105
Google Scholar
[13] Foley J M, Young S M, Phillips J D 2014 Phys. Rev. B 89 165111
Google Scholar
[14] Chen S S, Zhang W X, Yang B, Wu T, Zhang X D 2019 Sci. Rep. 9 5551
Google Scholar
[15] Zhang Y B, Liu W W, Li Z C, Li Z, Cheng H, Chen S, Tian J G 2018 Opt. Lett. 43 1842
Google Scholar
[16] Huo Y Y, Zhang X, Yan M, Sun K, Jiang S Z, Ning T Y, Zhao L N 2022 Opt. Express 30 19030
Google Scholar
[17] Yu J B, Ma B Z, OuYang A, Ghosh P, Luo H, Pattanayak A, Kaur S, Qiu M, Belov P, Li Q 2021 Optica 8 1290
Google Scholar
[18] Huang L J, Jin R, Zhou C B, et al. 2023 Nat. Commun. 14 3433
Google Scholar
[19] Zhou C B, Huang L J, Jin R, et al. 2023 Laser Photonics Rev. 17 2200564
Google Scholar
[20] Koshelev K, Lepeshov S, Liu M K, Bogdanov A, Kivshar Y 2018 Phys. Rev. Lett. 121 193903
Google Scholar
[21] Mikheeva E, Koshelev K, Choi D Y, Kruk S, Lumeau J, Abdeddaim R, Voznyuk I, Enoch S, Kivshar Y 2019 Opt. Express 27 33847
Google Scholar
[22] Chen X, Fan W H 2020 Nanomaterials 10 623
Google Scholar
[23] Cai Y P, Huang Y, Zhu K Y, Wu H H 2021 Opt. Lett. 46 4049
Google Scholar
[24] 杨磊 2022 硕士学位论文 (北京: 北京邮电大学)
Yang L 2022 M. S. Thesis (Beijing: Beijing University of Posts and Telecommunications
[25] Tang C J, Yan B, Wang Q G, Chen J, Yan Z D, Liu F X, Chen N B, Sui C H 2017 Sci. Rep. 7 582
Google Scholar
[26] Jin R, Huang L J, Zhou C B, et al. 2023 Nano Lett. 23 9105
Google Scholar
[27] Yu J B, Qin R, Ying Y B, Qiu M, Li Q 2023 Adv. Mat. 35 2302478
Google Scholar
[28] Chai R H, Liu W W, Li Z C, Cheng H, Tian J G, Chen S Q 2021 Phys. Rev. B 104 075149
Google Scholar
-
图 2 (a)面内破缺后超表面单元平面俯视图; (b)旋转角θ不同时, x偏振光入射对应的透射率曲线; (c)准BIC的Q因子随1/α2的变化曲线; (d)厚度h不同时, x偏振光入射对应的透射率曲线
Fig. 2. (a) Top view of the metasurface unit after in-plane symmetry perturbation; (b) transmittance spectra for x-polarized incident light and different rotation angle θ; (c) Q-factors of QBIC plotted as functions of 1/α2; (d) transmittance spectra for x-polarized incident light and different thickness h.
图 3 (a) x偏振光入射时的多极分解结果, ED, MD, TD, EQ, MQ分别表示电偶极子、磁偶极子、环偶极子、电四极子、磁四极子; (b)磁偶极子的散射功率在x, y, z方向上的分量; (c) x-y平面的电磁场分布; (d) y-z平面的电磁场分布
Fig. 3. (a) Multipole decomposition results during x-polarized normal incidence light, ED, MD, TD, EQ, MQ represent electric dipole, magnetic dipole, toroidal dipole electric, quadrupole and magnetic quadrupole, respectively; (b) x, y and z components of the MD; (c) electromagnetic field distribution in the x-y plane; (d) electromagnetic field distribution in the y-z plane.
图 4 (a)位移扰动后超表面单元平面俯视图; (b)参数Δg不同时, x偏振光入射对应的透射率曲线; (c)准BIC模式1的Q因子随1/β2的变化曲线; (d)准BIC模式2的Q因子随1/β2的变化曲线; (e)厚度h不同时, x偏振光入射对应的透射率曲线
Fig. 4. (a) Top view of the metasurface unit after displacement perturbation; (b) transmittance spectra for x-polarized incident light and different parameters Δg; (c) Q-factors of QBIC mode 1 plotted as functions of 1/β2; (d) Q-factors of QBIC mode 2 plotted as functions of 1/β2; (e) transmittance spectra for x-polarized incident light and different thickness h.
图 5 (a) x偏振光入射时产生的准BIC模式1的多极分解结果; (b)准BIC模式1的磁偶极子的散射功率在x, y, z方向上的分量; (c)准BIC模式1在x-z平面的电磁场分布; (d)准BIC模式1在y-z平面的电磁场分布; (e) x偏振光入射时产生的准BIC模式2的多极分解结果; (f)准BIC模式2的环偶极子的散射功率在x, y, z方向上的分量; (g)准BIC模式2在x-y平面的电磁场分布; (h)准BIC模式2在y-z平面的电磁场分布
Fig. 5. (a) Multipole decomposition results of QBIC mode 1 during x-polarized normal incidence light; (b) components of the MD in the x, y, z directions of QBIC mode 1; (c) electromagnetic field distribution in the x-z plane of QBIC mode 1; (d) electromagnetic field distribution in the y-z plane of QBIC mode 1; (e) multipole decomposition results of QBIC mode 2 during x-polarized normal incidence light; (f) components of the TD in the x, y, z directions of QBIC mode 2; (g) electromagnetic field distribution in the x-y plane of QBIC mode 2; (h) electromagnetic field distribution in the y-z plane of QBIC mode 2.
图 6 (a)当Px > Py时, 超表面单元平面俯视图; (b) Px不同时, x偏振光入射对应的透射率曲线; (c)准BIC模式的Q因子随1/γ22的变化曲线; (d)厚度h不同时, x偏振光入射对应的透射率曲线
Fig. 6. (a) When Px > Py, the top view of the metasurface unit; (b) transmittance spectra for x-polarized incident light and different Px; (c) Q-factors of QBIC plotted as functions of 1/γ22; (d) transmittance spectra for x-polarized incident light and different thickness h.
图 7 (a) x偏振光入射时的多极分解结果; (b)环偶极子的散射功率在x, y, z方向上的分量; (c) x-y平面的电磁场分布; (d) y-z平面的电磁场分布
Fig. 7. (a) Multipole decomposition results during x-polarized normal incidence light; (b) x, y and z components of the TD; (c) electromagnetic field distribution in the x-y plane; (d) electromagnetic field distribution in the y-z plane.
图 8 (a)当Py > Px时, 超表面单元平面俯视图; (b) Py不同时, x偏振光入射对应的透射率曲线; (c)准BIC模式的Q因子随1/γ22的变化曲线; (d)厚度h不同时, x偏振光入射对应的透射率曲线
Fig. 8. (a) When Py > Px, the top view of the metasurface unit; (b) transmittance spectra for x-polarized incident light and different Py; (c) Q-factors of QBIC plotted as functions of 1/γ22; (d) transmittance spectra for x-polarized incident light and different thickness h.
图 9 (a) x偏振光入射时的多极分解结果; (b)环偶极子的散射功率在x, y, z方向上的分量; (c) x-y平面的电磁场分布; (d) y-z平面的电磁场分布
Fig. 9. (a) Multipole decomposition results during x-polarized normal incidence light; (b) x, y and z components of the TD; (c) electromagnetic field distribution in the x-y plane; (d) electromagnetic field distribution in the y-z plane.
图 11 (a) x偏振光入射时产生的准BIC模式1的多极分解结果; (b)准BIC模式1的磁偶极子的散射功率在x, y, z方向上的分量; (c)准BIC模式1在x-y平面的电磁场分布; (d)准BIC模式1在y-z平面的电磁场分布; (e) x偏振光入射时产生的准BIC模式2的多极分解结果; (f)准BIC模式2的磁偶极子的散射功率在x, y, z方向上的分量; (g)准BIC模式2在x-z平面的电磁场分布; (h)准BIC模式2在y-z平面的电磁场分布; (i) x偏振光入射时产生的准BIC模式3的多极分解结果; (j)准BIC模式3的环偶极子的散射功率在x, y, z方向上的分量; (k)准BIC模式3在x-y平面的电磁场分布; (l)准BIC模式3在y-z平面的电磁场分布; (m) x偏振光入射时产生的准BIC模式4的多极分解结果; (n)准BIC模式4的磁偶极子的散射功率在x, y, z方向上的分量; (o)准BIC模式4在x-y平面的电磁场分布; (p)准BIC模式4在y-z平面的电磁场分布; (q) x偏振光入射时产生的准BIC模式5的多极分解结果; (r)准BIC模式5的磁偶极子的散射功率在x, y, z方向上的分量; (s)准BIC模式5在x-y平面的电磁场分布; (t)准BIC模式5在y-z平面的电磁场分布
Fig. 11. (a) Multipole decomposition results of QBIC mode 1 during x-polarized normal incidence light; (b) Components of the MD in the x, y, z directions of QBIC mode 1; (c) electromagnetic field distribution in the x-y plane of QBIC mode 1; (d) electromagnetic field distribution in the y-z plane of QBIC mode 1; (e) multipole decomposition results of QBIC mode 2 during x-polarized normal incidence light; (f) components of the MD in the x, y, z directions of QBIC mode 2; (g) electromagnetic field distribution in the x-z plane of QBIC mode 2; (h) electromagnetic field distribution in the y-z plane of QBIC mode 2; (i) multipole decomposition results of QBIC mode 3 during x-polarized normal incidence light; (j) components of the TD in the x, y, z directions of QBIC mode 3; (k) electromagnetic field distribution in the x-y plane of QBIC mode 3; (l) Electromagnetic field distribution in the y-z plane of QBIC mode 3; (m) multipole decomposition results of QBIC mode 4 during x-polarized normal incidence light; (n) components of the MD in the x, y, z directions of QBIC mode 4; (o) electromagnetic field distribution in the x-y plane of QBIC mode 4; (p) electromagnetic field distribution in the y-z plane of QBIC mode 4; (q) Multipole decomposition results of QBIC mode 5 during x-polarized normal incidence light; (r) components of the MD in the x, y, z directions of QBIC mode 5; (s) electromagnetic field distribution in the x-y plane of QBIC mode 5; (t) electromagnetic field distribution in the y-z plane of QBIC mode 5.
-
[1] von Neumann J, Wigner E P 1929 Phys. Z. 30 467
[2] Kaelberer T, Fedotov V A, Papasimakis N, Tsai D P, Zheludev N I 2010 Science 330 1510
Google Scholar
[3] Plotnik Y, Peleg O, Dreisow F, Heinrich M, Nolte S, Szameit A, Segev M 2011 Phys. Rev. Lett. 107 183901
Google Scholar
[4] Hsu C W, Zhen B, Lee J, Chua S L, Johnson S G, Joannopoulos J D, Soljacic M 2013 Nature 499 188
Google Scholar
[5] Urseel F 1997 Phys. Eng. Sci. 435 575
[6] Porter R, Evans D V 2005 Wave Motion 43 29
Google Scholar
[7] Xiao Y X, Ma G, Zhang Z Q, Chan C T 2017 Phys. Rev. Lett. 118 166803
Google Scholar
[8] Retzler C H 2001 Appl. Ocean Res. 23 249
Google Scholar
[9] Cobelli P J, Pagneux V, Maurel A, Petitjeans P 2009 Europhys. Lett. 88 20006
Google Scholar
[10] Cobelli P J, Pagneux V, Maurel A, Petitjeans P 2011 J. Fluid Mech. 666 445
Google Scholar
[11] EVERY A G 2008 Phys. Rev. B 78 174104.
Google Scholar
[12] Tittl A, Leitis A, Liu M, Yesilkoy F, Choi D Y, Neshev D N, Kivshar Y S, Altug H 2018 Science 360 1105
Google Scholar
[13] Foley J M, Young S M, Phillips J D 2014 Phys. Rev. B 89 165111
Google Scholar
[14] Chen S S, Zhang W X, Yang B, Wu T, Zhang X D 2019 Sci. Rep. 9 5551
Google Scholar
[15] Zhang Y B, Liu W W, Li Z C, Li Z, Cheng H, Chen S, Tian J G 2018 Opt. Lett. 43 1842
Google Scholar
[16] Huo Y Y, Zhang X, Yan M, Sun K, Jiang S Z, Ning T Y, Zhao L N 2022 Opt. Express 30 19030
Google Scholar
[17] Yu J B, Ma B Z, OuYang A, Ghosh P, Luo H, Pattanayak A, Kaur S, Qiu M, Belov P, Li Q 2021 Optica 8 1290
Google Scholar
[18] Huang L J, Jin R, Zhou C B, et al. 2023 Nat. Commun. 14 3433
Google Scholar
[19] Zhou C B, Huang L J, Jin R, et al. 2023 Laser Photonics Rev. 17 2200564
Google Scholar
[20] Koshelev K, Lepeshov S, Liu M K, Bogdanov A, Kivshar Y 2018 Phys. Rev. Lett. 121 193903
Google Scholar
[21] Mikheeva E, Koshelev K, Choi D Y, Kruk S, Lumeau J, Abdeddaim R, Voznyuk I, Enoch S, Kivshar Y 2019 Opt. Express 27 33847
Google Scholar
[22] Chen X, Fan W H 2020 Nanomaterials 10 623
Google Scholar
[23] Cai Y P, Huang Y, Zhu K Y, Wu H H 2021 Opt. Lett. 46 4049
Google Scholar
[24] 杨磊 2022 硕士学位论文 (北京: 北京邮电大学)
Yang L 2022 M. S. Thesis (Beijing: Beijing University of Posts and Telecommunications
[25] Tang C J, Yan B, Wang Q G, Chen J, Yan Z D, Liu F X, Chen N B, Sui C H 2017 Sci. Rep. 7 582
Google Scholar
[26] Jin R, Huang L J, Zhou C B, et al. 2023 Nano Lett. 23 9105
Google Scholar
[27] Yu J B, Qin R, Ying Y B, Qiu M, Li Q 2023 Adv. Mat. 35 2302478
Google Scholar
[28] Chai R H, Liu W W, Li Z C, Cheng H, Tian J G, Chen S Q 2021 Phys. Rev. B 104 075149
Google Scholar
计量
- 文章访问数: 4380
- PDF下载量: 177