基于电磁超材料的两种等效参数提取算法的比较分析*

丁敏 薛晖 吴博 孙兵兵 刘政 黄志祥† 吴先良2)‡

1)(安徽大学,计算智能与信号处理教育部重点实验室,合肥 230039)

2)(合肥师范学院物理与电子工程系,合肥 230061)

(2012年6月29日收到;2012年9月17日收到修改稿)

1 引言

电磁超材料(metamaterials)是一种异于天然的人工复合结构或复合材料[1],具有诸如负折 射 效 应[2]、 反 常 Doppler效 应[3]、 反 常Cherenkov效应[4]、完美透镜[5]等超常物理特性.近年来,对超材料结构实验方面的研究,已经从最初微波频段的开口谐振环(SRR)[6,7]和电导线结构逐步深入到光波频段的渔网结构[8-10].很多学者选用传统的反演算法[11,12]或基于Kramers-Kronig(K-K)关系[13-15]的改进算法,提取微波频段超材料结构的等效参数进行研究分析.本文主要利用这两种算法提取光波频段超材料结构的等效参数,分析比较两种算法的准确度与适用性.

为了确定两种算法的准确度和适用性,我们考虑了两类具有不同特性的渔网结构.将若干个超材料渔网结构小单元以固定的间隔d沿传播方向叠加,根据d的取值不同,将结构块分为弱耦合和强耦合的情况,分别利用两种算法对不同结构提取等效参数.计算表明,传统的反演算法适用性很强,对两类结构的超材料都可以准确地提取等效参数,但计算过程比较复杂;而基于K-K关系的改进算法巧妙地避免了复杂的计算,并能简单准确地提取材料的等效参数,但其使用范围具有一定的局限性.

2 算法理论

2.1 传统反演算法的基本理论

在研究超材料时,一般将其描述为等效媒质.在波长远大于结构单元的尺寸时,电磁波将不能识别出其内部结构,从而使得超材料可以用有效介电常数和有效磁导率来描述.计算中将超材料结构描述为一种有一定厚度的平板均匀结构,如果得到其传输、反射系数的强度大小和相位[8-10],则通过反演计算可以求得复数折射率n和波阻抗z,进一步可以得到媒质的等效介电常数和磁导率[16,17].

本文利用上述原理,对复数折射率n的实部进行研究分析.当平面波垂直入射到一定厚度的超材料结构上时,通过计算或测量得到其传输系数T,反射系数R,再利用(1),(2a)和(2b)式,可以得到其波阻抗z和折射率n：

其中,k0是自由空间的波数,deff是超材料结构的有效厚度,m是整数,代表不同分支.

利用超材料自身的特性条件可以确定波阻抗z和折射率n的正负符号.一般非增益超材料要求Re(z)＞0,因此可以得出(1)式的符号.同样地,Im(n)＞0,就惟一地确定了(2a)及(2b)式的符号,至于分支m的确定需要结合具体的物理概念,比如可依据一般超材料折射率具有连续性等限制条件决定m的具体取值.若是有效厚度deff取值较大,这些分支将十分靠近,很难准确地选择正确的分支数.传统反演算法的计算复杂度也就体现在分支数的选择上.尤其是多层渔网结构[18]的模型,分支选择的计算量就更大更复杂.

2.2 基于K-K关系的改进算法的基本理论

K-K关系的原理与传统反演算法基本一致.这种算法的新颖之处在于根据折射率的虚部利用KK关系得到折射率的实部,确保了等效参数的惟一性,避免了分支选择的复杂性,从而降低了计算的复杂度.K-K关系的本质是利用解析函数的实部和虚部的内在联系,其物理本质则在于物理系统因果率的体现.

基于K-K关系的算法步骤如下：波阻抗z可由T及R惟一地确定.折射率n的计算涉及一个复对数计算(5)式,这个复对数函数是一个多值函数,但多值性只针对n的实部.这时可用K-K关系由n的虚部来联系其实部,那么n的实部值就是取最接近K-K关系预测的对数函数对应的分支.但是,此算法要求折射率和频率是连续的.

该算法需要两个输入数据,一是有效厚度deff,即电磁波的传播方向上所包含的单元格的长度的总和;二是参数T及R.平面波垂直入射到均匀介质板上,波阻抗和折射率与T及R的关系如下：

其中,R01=(Zeff-1)/(Zeff+1),Zeff(ω)是复波阻抗,Neff(ω)=neff(ω)+ikeff(ω),是复折射率,neff(ω)是折射率,keff(ω)是消耗系数,k0是自由空间波数,ω是角频率.

其中,P.V.是柯西主值.

对(9)式可以利用梯形积分规则进行数值计算.但为了避免不当的K-K积分的奇异性,我们可以将积分分成两个部分,从而可以得到近似的K-K积分：

这里的Round()函数是取最接近的整数.因此,折射率的计算也就是选择最接近K-K关系预测值的那个分支.

3 数值算例

3.1 开口谐振环组成的超材料模型

首先,选用开口谐振环[6,7]组成的负磁导率超材料模型,分别用传统的反演算法和基于K-K关系的改进算法提取其等效参数,用来验证改进算法的适用性.如图1所示,边长为2.5 mm的空气立方体单元模型,中间嵌入厚度为0.25 mm的FR4基板,在基板的一侧是厚度为0.017 mm的铜材质开口谐振环.开口环外侧边长2.2 mm,开口间隙0.3 mm,内外环间距0.15 mm,环宽0.2 mm.图2是开口谐振环模型计算得到的传输系数、反射系数和吸收系数的结果.

图1 开口谐振环模型

如图3所示,传统的反演算法和基于K-K关系的改进算法结果基本符合,从而验证了改进算法的有效性.图3中曲线neff是用传统的反演算法计算出的折射率的实部,而曲线nK-K是用基于K-K关系的改进算法计算出的结果,曲线keff是折射率的虚部,根据结算结果可以看出针对此结构传统方法与基于K-K方法的结果符合得较好,这样从侧面证实了K-K方法的有效性.图4也给出了基于K-K关系反演得到的磁导率实部随频率变化曲线.

图2 开口谐振环模型的传输、反射和吸收系数

图3 开口谐振环模型的折射率

图4 开口谐振环模型的磁导率

3.2 渔网结构模型

3.2.1 强、弱耦合模型的确定

图5(a)是超材料渔网结构的一个单元模型[8,18],结构由金属银和电介质氟化镁(MgF2)组成,金属银分布在氟化镁(其介电常数为1.9)两侧.金属银是用Drude模型来模拟的,其等离子频率为1.37×1016rad/s,电碰撞系数为8.19×1013rad/s.传播方向是沿正z轴方向,结构在传播方向上是几何对称的.其电场方向(E)、磁场方向(H)和波矢方向(k),如图5(a)所示.其中,ax=ay=860 nm,wx=565 nm,wy=265 nm,s=50 nm,t=30 nm.

图5 (a)一个单元模型;(b)7层单元模型叠加的侧视图

图5(b)是两个单元模型的叠加,中间由空气层隔开,间距为d.根据d的取值不同,将叠加的模型分为弱耦合模型和强耦合模型.当d取某个阈值时,每个单元模型之间的相互影响很小,以这个特定的d而叠加的多个单元模型,结构的等效参数将保持稳定,不会随着结构的叠加而改变,这样的情况看作是弱耦合模型;当d取值超过阈值时,每个单元模型之间的相互影响很大,近似于光子晶体,这样的情况看作是强耦合模型.

图6是该渔网结构的谐振频率随d变化的曲线图.随着d的增加,结构的谐振频率也不断变化.在d取值大于90 nm时,谐振频率趋于固定值141.6 THz,不再随着d的增加而变化,即在90 nm处,谐振频率达到饱和.因此,界定当多个单元模型以d=90 nm为阈值间距沿传播方向增加时,构成的多层模型称为弱耦合模型.当d取值小于90 nm的情况均可视为强耦合模型.

图6 渔网结构的谐振频率随间距d的变化曲线图

3.2.2 多层弱耦合模型

确定强、弱耦合的阈值后,首先研究的是多层弱耦合模型,由7层单元模型沿传播方向叠加,每个单元模型之间的间距d=90 nm.图7是计算得到的传输系数、反射系数和吸收系数的结果.

图7 7层弱耦合模型的传输、反射和吸收系数

根据传输和反射系数,就可以利用传统的反演算法和基于K-K关系的改进算法计算出渔网结构的折射率.由于两种算法计算折射率的虚部的原理基本一致,所以其虚部曲线是完全一致的,也就是图8中的曲线keff.在图8中,曲线neff是用传统的反演算法计算出的折射率的实部,而曲线nK-K是用基于K-K关系的改进算法计算出的结果.可以看出频率在大概136—150 THz的区域内是负值,这体现了超材料负折射率的奇异特性.由于K-K积分中存在积分近似的过程,折射率的实部曲线nK-K与曲线neff大体趋势相同,计算结果表明两种方法在处理弱耦合时的一致有效性.而传统的反演算法在计算折射率的实部时,需要进行分支选择,计算过程比较复杂,而基于K-K关系的算法却巧妙地利用K-K关系得到折射率的实部,具有一定的计算优势.

图8 7层弱耦合模型的折射率

3.2.3 多层强耦合模型

多个单元模型叠加时,其间距d取小于90 nm的情况称为强耦合模型.此处以d=0 nm为进行数值模拟的强耦合模型.组成强耦合模型的单个模型之间的相互作用较强,所以它的折射率等效参数就会呈现出一些特别的情况,比如不连续的现象.通过计算分析,强耦合模型的折射率确实出现了不连续的情况,所以这时基于K-K关系的算法就无法适用了,因为它要求计算的参数是解析、连续的.

传统的反演算法可以用于计算不连续的模型,但是比较复杂,因为它涉及具体分支选择的问题.这里,计算的是19层、23层和27层强耦合模型,通过人工选择分支确定折射率的值.

图9是19层、23层和27层强耦合模型的传输系数曲线.当频率f＜136 THz时,传输系数T≈0,这是一种金属特性,也是f＜136 THz时Re(n)≈0的原因.

图10是折射率实部的不同分支曲线图,其中：m=1(青色),m=0(蓝色),m=-1(红色),m=-2(绿色),m=-3(品红色),m=-4(黑色).19层、23层和27层强耦合模型分别是用交叉符号,菱形符号和圆圈表示.

从图中可以看出,在频率136 THz附近,折射率的实部出现了不连续现象.频率小于136 THz时,折射率的实部趋于0.基于传统的反演算法,根据(2b)式中m的不同选择及不同层数重叠情形,可以确定出最终折射率实部的曲线如图10中阴影部分所示.

图9 19层、23层和27层强耦合模型的传输系数

图10 19层、23层和27层强耦合模型的折射率实部反演结果

4 结论

本文主要利用两种反演算法分析光波频段电磁超材料结构的等效媒质理论.给出了多层弱耦合和强耦合结构模型的阈值,利用多层弱耦合和多层强耦合渔网结构模型探讨了传统的反演算法和基于K-K关系的改进算法的适用性及有效范围.结果表明传统的反演算法适用性很强,对两种情况的模型都可以计算其等效参数,但是比较复杂,需人工选择分支来确定等效参数;而基于K-K关系的改进算法巧妙地利用K-K积分来避免了分支的选择,但是这种算法又存在一定的局限性,即计算模型的等效参数必须是解析的和连续的.我们的结果为电磁超材料的等效媒质理论提供了新的技术支持.

[1]Veselago V G 1968Sov.Hys.Usp.10 509

[2]Houck A A,Brock J B,Chuang I L 2003Phys.Rev.Lett.90 137401

[3]Seddon N,Bearpark T 2003Science302 1537

[4]Parazzoli C G,Greegor R B,Li K,Koltenbah B E C,Tanielian M 2003Phys.Rev.Lett.90 107401

[5]Pendry J B 2000Phys.Rev.Lett.85 3966

[6]Wang J F,Zhang J Q,Ma H,Yang Y M,Gu C,Qu S B,Xu Z 2009Acta Phys.Sin.58 3224(in Chinese)[王甲富,张介秋,马华,杨一鸣,顾超,屈绍波,徐卓2009物理学报58 3224]

[7]Li J C,Guo L X,Liu S H 2012Acta Phys.Sin.61 124102(in Chinese)[李俊成,郭立新,刘松华2012物理学报61 124102]

[8]Dolling G,Enkrich C,Wegener M,Soukoulis C M,Linden S 2006Opt.Lett.31 1800

[9]Smith D R 2010Phys.Rev.E 81 036605

[10]Minowa Y,Nagai M,Tao H,Fan K,Strikwerda A C,Zhang X,Averitt R D,Tananka K 2011IEEE Trans.Terahertz Sci.Tech.1 441

[11]Koschny T,Markos P,Economou E N,Smith D R,Vier D C,Soukoulis C M 2005Phys.Rev.B 71 245105

[12]Smith D R,Vier D C,Koschny T,Soukoulis C M 2005Phys.Rev.E 71 036617

[13]Lucarini V,Saarinen J J,Peiponen K E,Vartiainen E M 2005Kramers-Kronig Relations in Optical Materials Research(Berlin：Springer-Verlag)p28

[14]Peiponen K E,Saarinen J J 2009Rep.Prog.Phys.72 056401

[15]Akyurtlu A,Kussow A G 2010Phys.Rev.A 82 055802

[16]Tatartschuk E,Radkovskaya A,Shamonina E,Solymar L 2010Phys.Rev.B 81 115110

[17]Powell D A,Lapine M,Gorkunov M V,Shadrivov I V,Kivshar Y S 2010Phys.Rev.B 82 155128

[18]Valentine J,Zhang S,Zentgraf T,Ulin-Avila E,Genov D A,Bartal G,Zhang X 2008Nature455 376