用于熱光伏系統的近場輻射光譜控制表面結構?

于海童 劉東 楊震 段遠源

1)(清華大學,熱科學與動力工程教育部重點實驗室,二氧化碳資源化利用與減排技術北京市重點實驗室,北京 100084)

2)(南京理工大學能源與動力工程學院,南京 210094)

1 引 言

熱光伏(TPV)發電系統是一種新型能源轉換系統,核心部分為熱發射器和熱光伏電池兩部分,工作時,半導體電池通過光伏效應將高溫發射器的輻射電磁波轉化為電子空穴對,完成熱能到電能的轉換過程[1,2].熱光伏系統能量來源多樣,無運動部件,運行安靜穩定;與直接光伏轉換相比,熱光伏的工作時間不受晝間陽光照射時段所限,且由于熱源的輻射光譜可進行人為調控,因此可以設計窄帶熱輻射源代替寬譜帶太陽能輻射,使能量轉換效率打破單能帶太陽能光伏的肖克利-奎伊瑟極限(S-Q極限)[3].熱光伏系統的輸出功率受熱發射器輻射功率制約,而發射器溫度受材料耐熱性限制不能無限升高,因此為提升輸出功率,提出了近場熱光伏系統的思路,熱發射器與電池僅由亞微米量級真空間隙隔開,此時輻射電磁波的倏逝分量可通過隧穿效應通過間隙,承載由發射器到電池的能量傳遞,即近場輻射效應[4].通過利用倏逝波,近場輻射熱流能夠遠遠超過遠場輻射的黑體上限,實驗也證明了將發射器與吸收器的真空間隙縮減到2μm以下時,熱光伏系統的電流出現躍升[5],驗證了近場熱光伏系統的可行性.

為實現高效率的熱光伏能量轉化,需要有效操控系統中的輻射熱流光譜分布,這是由于被電池吸收的輻射電磁波中能量低于電池能帶的長波輻射無法轉化為光生電流,而能量過高的短波部分在激發電子-空穴對后還有剩余能量,兩種情況都直接導致能量轉換效率下降,并會產生熱量使電池性能惡化,因此熱光伏系統的設計與常規光伏系統一樣,都希望盡可能增強波長略短于電池能帶的熱流,而抑制其他波段的輻射傳熱.遠場光伏系統常在半導體電池表面添加減反射薄膜或陷光微結構以增加吸收率[6],但現有設計難以直接應用于近場熱光伏系統,原因如下:1)添加的表面結構或涂層可能與近場間距尺度接近,顯著改變發射器與吸收器間的等效間距和輻射傳遞機理;2)近場輻射傳熱多由倏逝波主導,而遠場減反射設計中未考慮倏逝分量,直接應用于近場會產生新的問題.實驗和理論均已揭示遠場增益輻射傳熱結構在近場間距下可能失去增強效果甚至削弱近場輻射傳熱,反之亦然[7,8].此外,熱光伏系統常用的III-V族半導體電池在能帶附近折射率高,吸收輻射電磁波能力差,增加了熱光伏系統提升功率和效率的難度[9].

尋找和設計適用于近場熱光伏應用的材料和結構,需要高效、嚴謹的近場輻射傳熱計算方法,不能簡單套用遠場的輻射和光學設計流程,一是因為基于黑體輻射定義的熱發射率在近場輻射傳熱中失去了意義,二是因為近場熱輻射可以激發平面內任意波數切向分量的行波或倏逝波模式,而基于遠場光學響應的計算方法只能考慮行波部分[10].近場輻射傳熱的直接計算需要以漲落耗散理論為基礎給出熱輻射的源項分布規律,再由系統的格林函數給出輻射電磁波的分布[11].對于一維多層結構,格林函數的形式可以被簡化,得到近場輻射熱流的解析解[12].對于更復雜結構參與的情況,目前研究多使用等效介質理論(EMT)[13?15]將非均勻結構通過等效法則近似為具有等效介電常數的均勻層,再使用一維結構解析解.Vongsoasup等[16]設計用鎢的矩形光柵作為熱發射器,利用EMT計算發現發射器具有雙曲超穎材料特性,相比無結構鎢發射器獲得了更高的熱流,但該設計使不能被利用的長波熱流被同時增強,能源轉換效率未能有效提高.Chang等[17]采用EMT計算了鎢納米線到InGaSb電池的近場傳熱,并使用薄膜電池,由電池的全內反射控制了無用的長波熱流,同時顯著提升了系統的熱流峰值與能源轉化效率.但事實上EMT對近場輻射的應用范圍并沒有定論[18],原因在于近場輻射中倏逝波的參與,且傳熱路徑尺度與結構尺度可能處于同一量級.此外也有研究證明對于金屬性組元,即使結構特征尺度遠小于特征波長,使用EMT仍可能與精確解偏差較大[19].

目前能夠嚴格表征復雜結構參與的近場輻射熱流的計算方法,一般是將漲落耗散理論與Maxwell方程的解相結合,改寫已有的電磁學計算方法以適用于近場輻射計算,使用的方法包括時域有限差分(FDTD)法[10,20]、傅里葉模態法(FMM,又稱嚴格耦合波分析)[21]、有限元法(FEM)等[22].其中,FDTD方法通過模擬系統對電磁波的時間域響應,通過傅里葉變換得到頻域響應規律,一次模擬即可獲得全光譜特性,編程簡單,適用于任何形狀的系統,且時域模擬易于表征近場輻射中常涉及的表面波傳遞等物理現象,是一種效率高、普適性強的近場輻射傳熱計算方法.

本文借助基于漲落耗散理論的FDTD近場輻射模擬方法,設計并驗證了一種用于近場熱光伏系統的輻射傳熱光譜控制表面結構,即III-V族電池表面的二維光柵結構.計算結果顯示帶表面光柵的GaSb電池與支持近紅外表面波傳播的發射器結合,可獲得電池能帶波長以下窄波段的選擇性高近場輻射熱流,峰值可達使用同溫度遠場黑體發射器時的2—3倍.使用FDTD方法模擬了系統的頻域和時域響應特性,揭示了一種新的近場輻射光譜控制機理,即表面柵格結構與發射器表面波的近場耦合,由此滿足高效率、高功率熱光伏系統的輻射傳熱設計要求.

2 模 型

2.1 物質與結構

熱光伏系統的輻射換熱組件包括電池與熱發射器兩部分.受熱發射器耐高溫能力限制,一般設計輻射熱流峰值和電池能帶在近紅外波段匹配以實現高功率與高轉換效率.電池通常選擇低能帶III-V族半導體,以最常用的GaSb為例,其能帶為0.726 eV,對應電磁波波長λ=1.7μm、頻率ω=1.1×1015rad/s.熱發射器方面,為激發窄波段內的高近場輻射熱流,可以使用支持表面等離子體激元(SPP)、表面聲子激元(SPhP)等表面波傳播的材料,但自然界中不存在支持近紅外表面波的材料,而能調制表面波頻率的納米結構難以滿足熱發射器的高溫熱穩定性需求.因此本文設計采用人工合成半導體——摻雜氧化鋅作為熱發射器,其SPP頻率位于近紅外區間,不需要附加表面結構即可激發SPP,且可通過摻雜組分類型與質量分數調制SPP頻率.考慮采用GaSb電池,故選擇摻雜比例為6%的ZnO:Ga(GZO),其SPP頻率(ωSPP=1.39×1015rad/s)恰好略高于電池能帶,以有選擇性地增加電池能帶之上的近場輻射熱流.圖1對比了幾種材料的近場輻射光譜熱流,使用GZO作為發射器獲得的輻射熱流比傳統耐高溫金屬(W,Mo,Ta等)更高,但與GZO同種材料間的輻射換熱相比,GZO到GaSb的輻射熱流峰值明顯較低,原因是GaSb在紅外波段的折射系數高(n≈4),導致其菲涅耳反射系數高,不能吸收GZO發射器的SPP窄帶輻射,光譜選擇性潛力沒有得到充分利用.

圖1 平行平面材料間的輻射光譜熱流qλ解析解(近場間距d=200 nm,兩側材料溫度為1800和300 K)Fig.1.Analytical solution of spectral radiative fl ux between bulk materials(1800 K to 300 K)separated by a gap of 200 nm.

圖2 計算模型示意圖 (a)物質、結構及尺寸參數定義;(b)數值計算域示意圖Fig.2. Illustration ofthecomputation model:(a)Materials,surface structures,and de fi nition of parameters;(b)numerical calculation domain.

為此,進一步設計對GaSb端進行改性以提升其近場輻射吸收能力.比擬遠場光伏系統中的減反射表面,在GaSb表面添加二維正方形柵格結構(圖2(a)),形狀由周期長度Λ、柱面邊長a與柱面高度h三個尺寸參數完全確定,相應地有光柵填充比例f=(a/Λ)2,槽寬w=Λ?a.已證明具有高尺寸比例的光柵、納米線等表面結構在遠場輻射問題中能有效增加III-V族半導體的吸收率[23],但對近場輻射傳熱的作用尚不清楚.此外,輻射傳熱計算中一般不需要針對異質結獨立建模,因為其光學特性與GaSb基體基本一致,對輻射傳熱計算的影響可以忽略[24].本文計算中除非另加說明,否則均假定近場間距d=200 nm,發射器溫度為1800 K,GaSb電池溫度為300 K.

2.2 模擬方法

使用FDTD方法計算GZO與帶表面結構GaSb間的輻射傳熱,模擬域如圖2(b)所示.普通FDTD方法計算發射率和吸收率時,一般使用遠場入射平面波作為源項研究系統的響應特性,該方法顯然無法直接模擬近場輻射傳熱.根據漲落耗散定律給出的熱輻射能量源項形式,在物體內部填充電流源項J,其分布規律服從[20,25,26]:

式中〈〉代表取統計平均值,?為取共軛,Im{}為取虛部;r為空間坐標矢量,i,j=x,y,z代表各坐標分量;δ為Kronecker符號,反映任意不同位置以及任意不同坐標分量電流源項的不相干性;ω為電磁波角頻率,ε(ω)為物質介電常數,ε0為真空介電常數;Θ= ?ω/{exp[?ω/(kBT)]?1}為普朗克諧振子能量,?為約化普朗克常數,kB為玻爾茲曼常數,T為熱力學溫度.Θ與Im{ε}的乘積揭示了電流源項的光譜強度受溫度和材料物性的影響規律.為避免電流源項的光譜相關性在FDTD算法中需要作額外傅里葉變換,本文使用Luo等[20]提出的Langevin法,以隨機白噪聲ˉJ作為計算中采用的實際源項形式,其特點是電流源光譜強度歸一化,彼此非相干,

在計算結束的后處理中,根據(1)式和(2)式的系數之比對熱流結果進行反歸一化,即得到正確結果.實際在離散網格時間域模擬中,在發射體內部填充的電流源項是單位振幅、相位隨機的高斯脈沖形式,使其在數值意義上逼近(2)式所要求的非相干性;計算時每個算例都運行多次模擬取平均,消除單次模擬中隨機噪聲的影響,如圖3(a)所示,本文算例在30—40次模擬后收斂.除源項處理部分外,本文FDTD方法架構與常規方法一致,計算隨機電流源激發的真空間隙中的電磁場強度,叉乘得到坡印廷矢量,即得到通過近場空隙的輻射熱流.物質的色散特性用Drude-Lorentz模型描述,并使用分段線性遞歸卷積(PLRC)方法完成頻域到時間域的變換處理.物質介電常數的Drude-Lorentz模型表達式為

圖3 近場FDTD算法驗證 (a)Langevin法過程與收斂速率展示;(b)FDTD與解析解對比計算GZO平板間的近場輻射熱流Fig.3.Validation of the near- fi eld FDTD calculation:(a)Illustration of Langevin approach and numerical convergence;(b)comparison of FDTD and analytical results for radiative transfer between GZO plates.

式中ωp為等離子振蕩頻率,ωl為Lorentz諧振頻率,Γp和Γl為阻尼系數.計算中使用的擬合參數如表1所示.

表1 計算中使用的物質Drude-Lorentz模型參數Table 1.Drude-Lorentz parameters for all materials used in calculation.

此外,雖然是周期性結構模擬,但由于電流源項具有不同波數分量,因此不能直接使用x,y方向的Bloch周期邊界條件;模擬域中包含數個結構周期,邊界均使用完美匹配層(PML)條件[30]吸收向外傳遞的電磁波.模擬域尺寸為x=y=2.5μm,z=1.2μm,經檢驗三個維度上尺寸分別加倍不會導致結果顯著變化.網格邊長為10 nm,計算結果取48次模擬平均值,單次模擬在Intel E5-2630 v3 CPU上單核運行時間48—60 h.圖3(b)中以GZO平板間近場輻射傳熱為算例,與解析解進行對比,驗證了近場FDTD結果的準確性.

3 計算結果與物理機理分析

3.1 半無限大GaSb電池

借助近場FDTD方法闡明添加表面光柵結構對GaSb吸收能力的影響.為使電池表面的光柵結構與發射體表面波激發的交替正負電荷聚集區域充分耦合,光柵周期Λ應設計為SPP的半波長.由SPP的色散關系得到其波長為[11]

式中SPP頻率ωSPP=1.39×1015rad/s,c為真空光速,ε0=1為真空相對介電常數. 由此可得λSPP=829 nm,為方便網格離散近似取Λ=400 nm.考慮結構加工的可行性,光柵高度不宜過大,計算中取h=200 nm,檢驗發現當h進一步增加時對結果影響不顯著.

在固定結構參數Λ,h的取值后,重點研究頂端邊長a的取值影響.對于有結構GaSb,不能用平行平板間的解析解驗證FDTD數值解的準確性,使用文獻[31]中的FMM半解析解作為對比.圖4所示為兩種典型參數條件下光柵結構對近場輻射熱流的影響,其中a=320 nm代表高填充比例(f=0.64)、光柵溝槽高寬比大(h/w=2.5)的情況,a=200 nm代表低填充比例(f=0.25)、光柵溝槽高寬比小(h/w=1.0)的情況.從圖4(a)可以看出,GaSb表面添加光柵結構有效提高了熱流的峰值和波長的選擇性,且FDTD結果與FMM符合良好,長波區FDTD熱流結果略偏低是因為模擬區域尺寸受內存限制,長波區計算精度降低.

圖4 表面光柵結構對GZO與GaSb間的輻射熱流影響(Λ=400 nm,h=200 nm) (a)光譜熱流結果;(b)相對使用遠場黑體發射器的熱流比值Fig.4. Radiative heat transfer between GZO and GaSb with and without surface structure(Λ =400 nm,h=200 nm):(a)Spectral radiative heat fl ux;(b)ratio over the condition using far- fi eld blackbody emitter.

為突出參數影響,圖4(b)展示了近場熱流與同溫度遠場黑體發射器到無結構GaSb熱流的比值.可見GZO到無結構GaSb的近場熱流已經在電池能帶內超過了黑體輻射情況,而添加GaSb表面結構又顯著增加了熱流峰值.其中,高填充比例的光柵參數對應較高的峰值熱流,a=320 nm時峰值熱流達到遠場黑體的2倍以上,是算例中峰值熱流最高的情況,但熱流增益光譜較寬,也增加了能帶以外的長波區的輻射傳熱,影響了能量轉化效率的提升.低填充比例(a=200 nm)時熱流峰值略低,但光譜選擇性更好,相比無結構情況的熱流增益全部位于電池能帶內,能夠實現更高的能量轉化效率.

3.2 表面減反射物理機理

相比其他計算方法,FDTD作為時域模擬方法能夠清晰展示電磁場的時間演化規律,適于解釋與表面波傳播相關的物理現象,這也是使用FDTD法的重要原因.為揭示GaSb表面光柵操控近場輻射熱流的物理機理,在FDTD的發射器表面中央添加單一電偶極子源項,設定其頻率為GZO的SPP頻率,觀察一個時間周期內的電場強度變化(圖5(a)).單點源項在GZO表面激發了正負電荷交替的聚集區,隨時間推移沿表面方向傳播,體現了GZO在近紅外波段的SPP特性,正負電荷的周期長度即為表面波波長.而在與GZO間隔200 nm的GaSb側,表面周期柵格結構的溝槽中同樣形成了電荷聚集,這是介電常數周期變化結構的陷光效應.除了正對偶極子的溝槽直接受點源影響外,兩側溝槽中電荷的符號與正對的GZO恰好相反,這是與GZO表面波跨越近場間隙耦合作用的結果,且間隙兩側電荷分布的空間周期一致,因此可隨時間演化保持穩定的近場間隙兩側電荷異號關系,輔助電磁波從發射器到吸收器的傳輸,由此增加SPP頻率附近的輻射熱流.物理機理簡化過程如圖5(b)所示.

圖5 單偶極子源項情形的時間域結果 (a)一個時間周期內的電場強度演化;(b)電荷分布簡化示意圖Fig.5. FDTD simulation with a single dipole source:(a)Evolution of electric fi eld intensity within one time period;(b)simpli fi ed illustration of charge distribution.

時域模擬結果說明了GaSb表面光柵調制近場輻射熱流的物理機理,即周期結構與發射器SPP耦合的陷光效應.為驗證光譜調制是發射器、吸收器兩端在近場間距耦合的結果,圖6給出了不同近場間距下的光譜熱流計算結果.GZO與GaSb距離很近(d=100—200 nm)時,輻射熱流在GZO的SPP位置附近有很高的峰值,而隨著間距增大,輻射熱流絕對值和光譜選擇性均迅速減弱,在d=1μm時熱流峰值已經完全消失.這證明了本文的光譜調控策略是一種特別適用于近場輻射傳熱的有效方法,只有發射器和吸收器間隔在近場尺度耦合才會誘發這種物理現象.

圖6 不同近場間距下GZO與GaSb間的輻射熱流(Λ=400 nm,a=320 nm,h=200 nm)Fig.6.Radiative heat fl ux between GZO and GaSb with various d(Λ=400 nm,a=320 nm,h=200 nm).

3.3 有限厚度GaSb電池

在添加GaSb表面結構的基礎上,可進一步將電池薄膜化,在高反射率的金屬基底上制備超薄膜III-V族半導體電池(圖7(a)),通過內全反射實現光譜選擇性吸收,提升光譜操控幅度[32].圖7(b)對比了無限厚與有限厚GaSb電池的近場輻射熱流,共同參數為Λ=400 nm,a=320 nm;對無限厚GaSb電池取h=200 nm,而對有限厚GaSb電池經優化計算取h1=100 nm,h2=250 nm.結果顯示使用超薄膜GaSb電池配合高反射率金屬底板進一步提升了輻射熱流的峰值和光譜選擇性,與遠場黑體熱流比值達到2.84,且增益區域全部位于電池能帶內,同時實現了高輻射功率與高轉化效率.

為說明基于嚴格結構的計算方法對揭示物理機理的意義,圖7(b)對比了兩種近場輻射傳熱的常用近似計算方法——EMT與Derjaguin近似[29].其中,EMT將非均勻介質層近似為具有等效介電常數的均勻介質,再應用一維結構解析解計算輻射熱流;對比采用各向異性的Maxwell等效法則,這是一種精度較高、適用范圍廣的EMT法則.Derjaguin理論則分別計算近場間距為d和d+h1的一維解,再根據光柵填充比例取加權平均.對比可見兩種近似方法均效果不佳,Derjaguin近似低估了全波段熱流,而EMT熱流結果過高且峰值位置有偏差,原因在于光柵層被等效為均勻的理想低折射率減反射層.這一結論與文獻[29]中金屬表面結構對近場輻射熱流的影響相似,也再次證明了GaSb表面光柵對近場輻射熱流的操控是由x,y方向介電常數的非均勻性引起的,使用近似方法不能描述電磁波與非均勻結構的作用機理,會對計算結果精度造成嚴重影響.

圖7 使用有限厚度GaSb作為吸收器的輻射熱流計算 (a)尺寸參數示意圖;(b)輻射熱流計算結果與近似方法的比較Fig.7.Radiative heat fl ux computation using GaSb absorber with fi nite thickness:(a)Illustration of size parameters;(b)radiative fl ux results compared with estimation methods.

4 結 論

設計了表面光柵結構用于增加GaSb對近場輻射傳熱的吸收能力,與具有近紅外SPP特性的GZO發射器搭配,作為提升近場熱光伏系統輸出功率和能源轉化效率的方案.使用結合漲落耗散定律的隨機FDTD方法,直接對有復雜結構參與的近場輻射傳熱進行嚴格計算,時域和頻域模擬結果揭示了一種新的近場輻射光譜調控機理,即與表面波耦合的陷光效應.將表面光柵結構與超薄膜電池結合,增強了輻射熱流光譜的調控幅度,使近場熱流峰值達到使用同溫度遠場黑體輻射源情況下的2.84倍,且熱流增益區域全部位于波長略短于電池能帶的窄波段,由此可有效提升近場熱光伏系統的輸出效率和轉換功率.

[1]Liu D,Yu H T,Yang Z,Duan Y Y 2015J.Eng.Thermophys.36 698(in Chinese)[劉東,于海童,楊震,段遠源2015工程熱物理學報36 698]

[2]Coutts T J 1999Renewable and Sustainable Energy Reviews3 77

[3]Lenert A,Bierman D M,Nam Y,Chan W R,Celanovi C I,Solja?i? M,Wang E N 2014Nat.Nanotechnol.9 126

[4]Basu S,Chen Y,Zhang Z M 2007Int.J.Energ.Res.31 689

[5]Hanamura K,Fukai H,Srinivasan E,Asano M,Masuhara T 2011ASME/JSME 8th Thermal Engineering Joint ConferenceHawaii,USA,March 2011

[6]Geng C,Zheng Y,Zhang Y Z,Yan H 2016Acta Phys.Sin.65 070201(in Chinese)[耿超,鄭義,張永哲,嚴輝2016物理學報65 070201]

[7]Ijiro T,Yamada N 2015Appl.Phys.Lett.106 23103

[8]Chalabi H,Hasman E,Brongersma M L 2015Phys.Rev.B91 14302

[9]Molesky S,Jacob Z 2015Phys.Rev.B91 205435

[10]Lu D,Das A,Park W 2017Opt.Express25 12999

[11]Zhang Z M 2007Nano/Microscale Heat Transfer(New York:McGraw-Hill)p377

[12]Francoeur M,Mengü? M P,Vaillon R 2009J.Quant.Spectr.Radiat.Transfer110 2002

[13]Li J Y,Xuan Y M,Li Q,Han Y G 2013J.Eng.Thermophys.34 1548(in Chinese)[李佳玉,宣益民,李強,韓玉閣2013工程熱物理學報34 1548]

[14]Wu H H,Huang Y,Zhu K Y 2016J.Eng.Thermophys.37 597(in Chinese)[吳會海,黃勇,朱克勇2016工程熱物理學報37 597]

[15]Zhu K Y,Huang Y,Wu H H 2016J.Eng.Thermophys.37 2393(in Chinese)[朱克勇,黃勇,吳會海 2016工程熱物理學報37 2393]

[16]Vongsoasup N,Francoeur M,Hanamura K 2017Int.J.Heat Mass Transfer115 326

[17]Chang J Y,Yang Y,Wang L 2015Int.J.Heat Mass Transfer87 237

[18]Zhang R Z,Zhang Z M 2017J.Quant.Spectr.Radiat.Transfer197 132

[19]Yu H T,Liu D,Duan Y Y,Zhen Y 2015Int.J.Heat Mass Transfer87 303

[20]Luo C,Narayanaswamy A,Chen G,Joannopoulos J D 2004Phys.Rev.Lett.93 213905

[21]Lussange J,Guérout R,Rosa F S S,Greffet J J,Lambrecht A,Reynaud S 2012Phys.Rev.B86 85432

[22]Bai Y,Jiang Y,Liu L 2015J.Quant.Spectr.Radiat.Transfer158 36

[23]Kanamori Y,Kobayashi K,Yugami H,Hane K 2003Jpn.J.Appl.Phys.42 4020

[24]Bernardi M P,Dupré O,Blandre E,Chapuis P O,Vaillon R,Francoeur M 2015Sci.Rep.5 11626

[25]Didari A,Mengü? M P 2017J.Quant.Spectr.Radiat.Transfer197 95

[26]Datas A,Hirashima D,Hanamura K 2013J.Therm.Sci.Tech.8 91

[27]Kim J,Naik G V,Emani N K,Guler U,Boltasseva A 2013IEEE J.Sel.Top.Quant.19 4601907

[28]Djuri?i? A B,Li E H,Raki C D,Majewski M L 2000Appl.Phys.A70 29

[29]Yang Y,Wang L 2016Phys.Rev.Lett.117 44301

[30]Wei B,Li X Y,Wang F,Ge D B 2009Acta Phys.Sin.58 6174(in Chinese)[魏兵,李小勇,王飛,葛德彪2009物理學報58 6174]

[31]Yu H,Liu D,Yang Z,Duan Y 2017Sci.Rep.7 1026

[32]Kim J,Hwang J,Song K,Kim N,Shin J C,Lee J 2016Appl.Phys.Lett.108 253101