SITE ZHANG,1,2,*CHRISTIAN HELLMANN,2 AND FRANK WYROWSKI1
1Applied Computational Optics Group, Institute of Applied Physics, Friedrich Schiller University Jena, Max-Wien-Platz 1, 07743 Jena, Germany
2Wyrowski Photonics UG, Kahlaische Straße 4, 07745 Jena, Germany
*Corresponding author: site.zhang@uni‑jena.de
Received 8 March 2017; revised 29 April 2017; accepted 3 May 2017; posted 3 May 2017 (Doc. ID 290298); published 19 May 2017
Received 8 March 2017; revised 29 April 2017; accepted 3 May 2017; posted 3 May 2017 (Doc. ID 290298); published 19 May 2017
通過(guò)使用平面波譜分析以及S矩陣方法，我們研究了一般電磁場(chǎng)經(jīng)過(guò)由各向同性介質(zhì)或者各向異性介質(zhì)構(gòu)成的光學(xué)層狀結(jié)構(gòu)的傳播。我們還開(kāi)發(fā)了一種基于快速傅里葉變換技術(shù)的算法，具有數(shù)值高效的采樣規(guī)則。通過(guò)將此算法與其它系統(tǒng)建模技術(shù)相結(jié)合，我們展示了一些仿真案例，如經(jīng)過(guò)一個(gè)各向同性法珀標(biāo)準(zhǔn)具的光場(chǎng)傳輸以及具有任意方位和光軸方向的單軸晶體平板的光場(chǎng)傳輸。
OCIS codes: (260.0260) Physical optics; (260.2110) Electromagnetic optics; (260.1440) Birefringence; (230.4170) Multilayers.
https://doi.org/10.1364/AO.56.004566
1.引言
由平行平面構(gòu)成的光學(xué)層在光學(xué)中廣泛應(yīng)用。層狀結(jié)構(gòu)可以用作許多不同情況的模型，像平板和標(biāo)準(zhǔn)具�；谶@個(gè)事實(shí)，光與層狀結(jié)構(gòu)相互作用的主題一直引起大家的注意并且對(duì)此已經(jīng)進(jìn)行了大量的研究。
在這類研究中，大多數(shù)觀點(diǎn)都側(cè)重于平面波，然而僅僅少數(shù)的研究使用了平面波譜方法（SPW）來(lái)考慮一般的電磁場(chǎng)。例如，參考文獻(xiàn)[1-6]中研究了各向同性-各向同性的界面上，高斯光束的反射率和透射率；在參考文獻(xiàn)[7-11]中研究了各向同性層或者平板的情況；參考文獻(xiàn)[12-22]討論了各向同性-各向異性界面的情況，在參考文獻(xiàn)[23-26]中則討論了各向異性層或者平板的情況。
上面所提到的許多研究都用于特定的研究主題，像[1,3,5]中研究了高斯光束全內(nèi)反射的橫向偏移，并且他們常常關(guān)注于具體的配置。因此，將這些方法推廣到更一般的情況的可能性受到了限制。
在這篇文章中，我們從一個(gè)更一般的觀點(diǎn)來(lái)考慮此問(wèn)題。光學(xué)層幾乎不會(huì)單獨(dú)使用；相反，他們常常是一個(gè)光學(xué)系統(tǒng)的一部分并且和其他的元件一起使用，如圖1中所示。基于此事實(shí)，我們遵循場(chǎng)追跡的概念[27]，并使用不同的場(chǎng)追跡算子組合[28-32]，如圖1中所示，以對(duì)一個(gè)包含了層介質(zhì)元件的系統(tǒng)進(jìn)行物理光學(xué)模擬。考慮到模擬是對(duì)整個(gè)系統(tǒng)而不是單個(gè)元件，仿真層結(jié)構(gòu)必須與系統(tǒng)的前后部分相連接。這要求我們傳播步驟（圖1中的P）進(jìn)行適當(dāng)?shù)目紤]，將前一個(gè)元件的輸出連接到當(dāng)前元件的輸入，并將當(dāng)前元件的輸出傳遞到下一個(gè)元件。一般情況下，這樣的傳輸步驟會(huì)出現(xiàn)在平行或者非平行平面之間。在參考文獻(xiàn)[28,29]中已經(jīng)提到了平行平面間幾種有效的傳輸方法，在參考文獻(xiàn)[33]中則可以找到對(duì)非平行平面間傳輸?shù)囊粋�(gè)詳細(xì)的討論。在這篇文章中，我們不會(huì)研究傳輸步驟，但會(huì)關(guān)注層狀結(jié)構(gòu)的元件算子C。
此外，從數(shù)值計(jì)算的觀點(diǎn)出發(fā)，為了執(zhí)行一個(gè)連續(xù)且有效的系統(tǒng)模擬，要求元件算子C
正確地處理采樣場(chǎng)數(shù)據(jù)并和其他的算子以一種統(tǒng)一的格式傳遞場(chǎng)數(shù)據(jù)；
優(yōu)化數(shù)值計(jì)算的效率。
考慮到上述兩個(gè)標(biāo)準(zhǔn)，我們開(kāi)發(fā)了一種具有自動(dòng)數(shù)值采樣規(guī)則的SPW方法。與之前一些利用積分方法對(duì)空間和角譜相關(guān)的傅里葉變換進(jìn)評(píng)估的研究相比（如參考文獻(xiàn)[23]中的二維中點(diǎn)規(guī)則和參考文獻(xiàn)[12-14，20，25]中的Stamnes–Spjelkavik–Pedersen方法[34]），我們使用了快速傅里葉變換（FFT）技術(shù)，此技術(shù)在大部分?jǐn)?shù)值軟件包中容易訪問(wèn)并且效率高。再加上在角譜域中經(jīng)過(guò)深入考慮的數(shù)值采樣規(guī)則，我們的方法具有一般適用性，對(duì)層元件和入射場(chǎng)沒(méi)有任何限制。因此，此算法可以直接包含在一個(gè)物理光學(xué)系統(tǒng)模擬之中。
圖1.結(jié)合使用不同的場(chǎng)追蹤算子來(lái)模擬光學(xué)系統(tǒng)： C是元件算子，P是相鄰元件之間的傳輸算子。
2.理論
如圖2所示，層狀結(jié)構(gòu)分別由兩個(gè)位于  和  的平行平面構(gòu)成。  和  的區(qū)域充滿了復(fù)折射率為  和  的均勻各向同性介質(zhì)。參考文獻(xiàn)[27]中表明使用橫向分量Ex和Ey已足夠表征均勻各向同性介質(zhì)中電磁場(chǎng)了。因此，我們可以使用以下表達(dá)式來(lái)描述此問(wèn)題：
其中，分別在平面 和 處定義輸入和輸出橫向電場(chǎng)矢量，（兩者位于界面的數(shù)學(xué)位置，但總是認(rèn)為在均勻介質(zhì)的一側(cè)），由下式給出
其中  。方程（1）中的元件算子是一個(gè)2x2的矩陣形式，
圖2.層狀結(jié)構(gòu)分別由兩個(gè)位于 和 的平行平面構(gòu)成。 和 的區(qū)域由均勻各向同性介質(zhì)填充，其折射率分別是 和 。輸出場(chǎng)和輸出場(chǎng)在層表面進(jìn)行定義，但總是在相應(yīng)的各向同性介質(zhì)的一側(cè)。
在這章節(jié)，我們的目標(biāo)是找到C的精確的形式，以連接層介質(zhì)元件的輸入和輸出場(chǎng)。為了研究與層結(jié)構(gòu)的相互作用，我們對(duì)輸入橫向場(chǎng)分量進(jìn)行了一個(gè)傅里葉變換，并獲得了
其中  , F表示二維傅里葉變換，
 。逆傅里葉變換定義如下
方程（6）中的積分可以解釋為將  分解為具有不同橫向波矢分量κ的平面波。因此，在我們的情況下，每個(gè)輸入平面波都可以單獨(dú)處理——我們首先計(jì)算每個(gè)輸入平面波的輸出，然后進(jìn)行求和從而獲得輸出場(chǎng)。
此外，根據(jù)邊界條件對(duì)電磁場(chǎng)施加的連續(xù)性要求，可以顯示出一個(gè)給定的輸入平面波在與層結(jié)構(gòu)相互作用的過(guò)程中其橫向波矢分量κ必定保持不變。同樣可以顯示出，通過(guò)疊加原理的有效性，不同的κ之間沒(méi)有耦合。因此，對(duì)于輸出角譜，我們可以寫下
其中
公式（8）中分別是透射和反射系數(shù)矩陣。為了計(jì)算T(κ) 或者R(κ)，我們選擇使用數(shù)值穩(wěn)定S矩陣方法。為了計(jì)算S矩陣，首先必須確定每個(gè)各向異性層的平面波。基于文獻(xiàn)[35]中Berreman的4x4矩陣公式，Landry和Maldonado開(kāi)發(fā)并展示了一種數(shù)值友好的形式，見(jiàn)參考文獻(xiàn)[23]。我們采用了他們的方法，對(duì)于每一層，求解了參考文獻(xiàn)[23]中由方程（28）所描述的特征系統(tǒng)的特征值和特征向量。
不同于[23，25]中直接使用本征解來(lái)構(gòu)建一個(gè)轉(zhuǎn)換矩陣，另外，我們還需要根據(jù)他們的傳輸方向整理出平面波，這是為計(jì)算S矩陣所做的一個(gè)必需的準(zhǔn)備。為此，我們遵循[36]中4.3部分由Li所提出的標(biāo)準(zhǔn)。
然后可以應(yīng)用遞歸S矩陣公式。我們?cè)谶@篇文中不再重復(fù)給出已發(fā)展成熟的S矩陣方法，讀者可以參考文獻(xiàn)[37]中的方程（5）-（8）以獲得更多的信息。在我們的情況中，由于沒(méi)有反向傳輸輸入場(chǎng)，我們僅對(duì)正向透射或者反向反射感興趣，因此這篇文章中的矩陣系數(shù)T(κ) 和R(κ)對(duì)應(yīng)于[37]中方程（5c）或者（5d）中的子矩陣T_uu 或者R_du。
一旦獲得了矩陣系數(shù)，通過(guò)方程(7)即可獲得輸出角譜。對(duì)輸出角譜進(jìn)行一個(gè)逆傅里葉變換，我們獲得了輸出橫向場(chǎng)矢量
通過(guò)聯(lián)合方程（4），（7），（9），我們可以寫出從輸入場(chǎng)到輸出場(chǎng)整個(gè)計(jì)算流程，如果
我們以透射的情況作為例子，則
因此，方程（1）中元件算子C的精確形式如下
通過(guò)使用系數(shù)矩陣R(κ)代替 T(κ)，可以獲得反射情況的表達(dá)式。
3.算法
按照方程（10）的順序，我們可以應(yīng)用一種數(shù)值算法以計(jì)算場(chǎng)經(jīng)過(guò)分層介質(zhì)元件的傳輸。讓我們從橫向輸入場(chǎng)矢量  和  ，以圖3（a）中的均勻網(wǎng)格進(jìn)行采樣。這種網(wǎng)格定義為  ,其中  和  作為索引數(shù)，δx和δy是x方向和y方向的采樣距離。初始采樣參數(shù)應(yīng)該受到合適的控制以使他們符合先前算子的奈奎斯特-香農(nóng)采樣定理。方程（10）中的F為聯(lián)系空間域和頻率域的算子，可以使用不同的數(shù)值方法實(shí)現(xiàn)，像廣泛使用的快速傅里葉變換（FFT）技術(shù)，以及包含了能夠進(jìn)一步提高數(shù)值效率的半解析傅里葉變換[38]和啁啾z變換[39-41]的更高級(jí)的方法。在此篇文章中，我們使用了FFT技術(shù)，并以此獲得輸入角譜。但我們的算法不受限于該技術(shù)。
圖3.在角譜域定義均勻采樣網(wǎng)格。（a），初始網(wǎng)格定義，5x5采樣點(diǎn)，  和  為采樣距離；（b），使用5x5采樣點(diǎn)，沿垂直方向生成一個(gè)測(cè)試網(wǎng)格， 和  =0.5 為采樣距離；（c）中使用9x5個(gè)采樣點(diǎn)定義了沿水平方向的一個(gè)測(cè)試網(wǎng)格，  =0.5 和 為采樣距離。（（b）和（c）中的實(shí)心點(diǎn)是在（a）的初始網(wǎng)格中出現(xiàn)的共同的采樣點(diǎn)位置，而空心圓環(huán)與初始網(wǎng)格中不一樣的采樣點(diǎn)。
其中角譜域中的均勻網(wǎng)格定義為  ，其中 和 為Kx和Ky的采樣距離。利用FFT計(jì)算的結(jié)果，在兩個(gè)域中的采樣點(diǎn)是一樣的，因此我們有  。
接下來(lái)，將輸入角譜乘以系數(shù)矩陣，然后我們可以獲得輸出角譜場(chǎng)，例如，對(duì)于透射
角譜  的采樣通過(guò)傅里葉變換關(guān)系自動(dòng)確保。然而，  的采樣沒(méi)有必要進(jìn)行確保，因?yàn)樗麄儼讼禂?shù)和輸入角譜之間的點(diǎn)乘。這個(gè)問(wèn)題將在第四部分的案例中進(jìn)行清楚地說(shuō)明。
一般來(lái)說(shuō)，采樣距離 ， 和帶寬  必須進(jìn)行適當(dāng)調(diào)整以確保正確的采樣。在我們的情況中，方程（13）中的運(yùn)算沒(méi)有改變譜范圍。因此，我們僅僅只需要找到合適的采樣距離即可。此問(wèn)題并沒(méi)有解析解，因?yàn)橄禂?shù)矩陣T或者R是通過(guò)遞歸S矩陣方法數(shù)值上獲得的。確保 采樣的唯一方法是實(shí)行數(shù)值測(cè)試。
讓我們?cè)诮亲V域定義一個(gè)測(cè)試采樣網(wǎng)格，如如3（b）或者3（c）所示。測(cè)試網(wǎng)格定義為  ,其中  為指數(shù)；為了不改變角譜范圍，要求  。與初始網(wǎng)格相比，它也需要進(jìn)行細(xì)化，這意味著  。
一方面，我們以一種嚴(yán)謹(jǐn)?shù)姆绞絹?lái)計(jì)算測(cè)試網(wǎng)格上的輸出橫向角譜分量。由于輸入橫向角譜分量的合適采樣，根據(jù)需要可以進(jìn)行插值，例如，在測(cè)試網(wǎng)格上，我們獲得
其中“INTRPL”代表的是數(shù)值插值運(yùn)算。為了在一個(gè)更精細(xì)的網(wǎng)格上嚴(yán)格地獲得輸入場(chǎng)，如方程（14）中，我們總是采用基于FFT的Sinc插值方法，把它們代入到方程（13）中，可以嚴(yán)格地獲得輸出橫向角譜量
另一方面，通過(guò)對(duì)  的插值，我們也可以獲得測(cè)試網(wǎng)格上的輸出橫向角譜量，并將插值結(jié)果表示為  。我們將此過(guò)程描述為
值得強(qiáng)調(diào)的是，在方程（16）中無(wú)需再使用嚴(yán)格的Sinc插值。因?yàn)檫@些結(jié)果僅僅用于評(píng)估一下情況。在這篇文章中，我們使用三次插值。上面的插值沒(méi)必要給出正確的輸出值，因此，在方程（16）的左邊，我們對(duì)那些量使用上標(biāo)“intrpl”。下一步，我們定義
作為插值結(jié)果和嚴(yán)格結(jié)果之間的相對(duì)偏差。只要  在初始網(wǎng)格上的采樣滿足奈奎斯特-香農(nóng)采樣定理，基于他們的插值結(jié)果就不會(huì)對(duì)嚴(yán)格的插值方法表現(xiàn)出一個(gè)大的偏差。在這種情況下，根據(jù)所選的插值方法，兩種結(jié)果的偏差應(yīng)該在數(shù)值誤差級(jí)次，由σ_0表示。使用上述標(biāo)準(zhǔn)，通過(guò)一步步的減小采樣間距來(lái)測(cè)試場(chǎng)采樣，直至σ<σ_0。
以上的測(cè)試過(guò)程對(duì)應(yīng)著一個(gè)循環(huán)，在每個(gè)測(cè)試循環(huán)中都會(huì)執(zhí)行方程（14）到（17）的計(jì)算。為了在每個(gè)測(cè)試循環(huán)中充分地使用測(cè)試數(shù)據(jù)，如方程（15）中的結(jié)果，我們總是使用圖3（b）和3（c）中的測(cè)試網(wǎng)格。這樣規(guī)則的測(cè)試網(wǎng)格定義為
對(duì)每個(gè)循環(huán)內(nèi)，需要對(duì)先前的采樣距離減半。使用這個(gè)測(cè)試網(wǎng)格具有如下好處：當(dāng)計(jì)算  時(shí)，僅僅需要計(jì)算那些在空心圓圈位置處的值，而實(shí)心點(diǎn)位置處的值在之前的測(cè)試循環(huán)中已經(jīng)計(jì)算過(guò)了，并且讓這些值簡(jiǎn)單地進(jìn)行傳輸即可。在這種方式下，這些基于測(cè)試目的所進(jìn)行嚴(yán)格計(jì)算的值不會(huì)被丟掉，而會(huì)作為下一次循環(huán)中的起始點(diǎn)以重新使用。因此，在一個(gè)完整的計(jì)算中所執(zhí)行的是嚴(yán)格的計(jì)算，從而驗(yàn)證合適的采樣，并且這些計(jì)算的值能夠有效地用于構(gòu)建最終輸出場(chǎng)。我們將上面的方法總結(jié)在算法1中。
算法1：通過(guò)合適的采樣控制來(lái)進(jìn)行輸出角譜計(jì)算的數(shù)值流程
1)在初始網(wǎng)格嚴(yán)格地計(jì)算
2)初始化相對(duì)誤差值σ=+∞
3)判斷σ>σ_0， ⊳如果為真，則采樣不合適
4)沿Kx（或者Ky）二等分采樣距離
5)應(yīng)用新的采樣距離定義測(cè)試網(wǎng)格 
6)在測(cè)試網(wǎng)格中對(duì)輸入角譜插值，根據(jù)方程（14）獲得 
7)根據(jù)方程（15），在測(cè)試網(wǎng)格上嚴(yán)格地計(jì)算輸出角譜 
8)根據(jù)方程（16），通過(guò)插值獲得 
9)根據(jù)方程（17），計(jì)算誤差σ
10)如果σ>σ_(0 )  ，則⊳將目前的輸出場(chǎng)設(shè)置為下一次循環(huán)的起始點(diǎn)
11)返回
為了有效地處理非對(duì)稱情況，例如，光束在x方向和y方向有不同的發(fā)散角，算法1中的測(cè)試需要按順序沿兩個(gè)方向進(jìn)行。開(kāi)始方向的選擇是任意的，在我們的情況中，我們是沿y方向開(kāi)始測(cè)試的。
4.示例
在VirtualLab Fusion[42]軟件中，我們將第三節(jié)中提出的算法實(shí)現(xiàn)在“可編程元件”的編程界面中。這個(gè)元件可以與VirtualLab Fusion中其它的物理光學(xué)仿真技術(shù)進(jìn)行聯(lián)合仿真。接下來(lái)，我們展示了四個(gè)案例：前兩個(gè)主要關(guān)注元件本身并以一種嚴(yán)格的數(shù)值方式檢查此算法；后兩個(gè)案例中，元件將用于光學(xué)系統(tǒng)中，例如，此算法與其它仿真技術(shù)一起使用。
在進(jìn)入實(shí)際的案例之前，設(shè)置方程（17）中的迭代終止標(biāo)準(zhǔn)σ_0很重要。對(duì)于在此文中所使用的三次插值，我們預(yù)先檢查了它在一般情況下的表現(xiàn)并在我們的數(shù)值環(huán)境中找到了一個(gè)0.01的基準(zhǔn)值。
A.各向同性標(biāo)準(zhǔn)具
第一個(gè)案例模擬了一個(gè)線性偏振高斯光束經(jīng)過(guò)一個(gè)標(biāo)準(zhǔn)具的傳輸。通過(guò)這個(gè)案例，我們將清楚的指出第三節(jié)中所說(shuō)的采樣問(wèn)題并描述了算法1的工作原理。標(biāo)準(zhǔn)具由熔融石英制成，兩側(cè)有多層薄膜，如圖4所示。關(guān)于標(biāo)準(zhǔn)具的光學(xué)參數(shù)和結(jié)構(gòu)的更多信息，請(qǐng)見(jiàn)表1。輸入場(chǎng)為波長(zhǎng)633nm，x方向線偏光的高斯光束。其在元件的輸入平面定義為E_⊥^in,且束腰半徑為（2um，2um）。在經(jīng)過(guò)傅里葉變換后，我們獲得其角譜，同樣具有高斯輪廓，如圖5所示。按照方程（10）中的操作算子序列，輸入角譜將乘以透射或者反射率系數(shù)。我們?nèi)砸酝干渥鳛槔�，并且�?duì)于線性偏振輸入場(chǎng)，我們使用t_xx和t_yx乘以E ̃_x^in，以獲得輸出角譜分量。
圖4.由熔融石英制成，兩側(cè)有多層膜的標(biāo)準(zhǔn)具。其結(jié)構(gòu)和光學(xué)參數(shù)如表1中所示。
表1 標(biāo)準(zhǔn)具的結(jié)構(gòu)和光學(xué)參數(shù)
圖5 .（a）輸入高斯場(chǎng)分量的振幅；（b）對(duì)應(yīng)的角譜分量。由于輸入光場(chǎng)為沿x方向的線性偏振光，因此僅顯示Ex分量。
如第三節(jié)中所指出，乘積  的采樣不能自動(dòng)得到保證，此案例中將顯示該現(xiàn)象。標(biāo)準(zhǔn)具由于其頻率選擇功能（頻譜或角頻率）而得到廣泛的使用。在我們的案例中，角頻率選擇可以解釋為系數(shù) txx和 tyx以一種方式調(diào)制輸入角譜，以使特定的角頻率加強(qiáng)而其它的減弱。這種調(diào)制可以出現(xiàn)在一種非常精細(xì)的頻率水平上。因此，需要使用更精細(xì)的采樣以在輸出角譜中解析這樣一個(gè)精細(xì)的調(diào)制。為了獲得需要的采樣間距，我們遵循算法1，圖6中顯示了部分結(jié)果。
通過(guò)圖6和表2，我們根據(jù)算法1中的步驟描述了工作流程，如下：
第一步：從圖6（a）中所示的輸入角譜開(kāi)始，計(jì)算各個(gè)系數(shù)并乘上元件矩陣以生成初始化的輸出角譜；
第二步：初始化相對(duì)偏差σ=+∞；
第三步：開(kāi)始測(cè)試循環(huán)；
第四&五步：將采樣距離沿κ_x或者κ_y方向減半，以定義測(cè)試網(wǎng)格 ，對(duì)應(yīng)的采樣點(diǎn)如表2中所示；
第六步：對(duì)輸入角譜在測(cè)試網(wǎng)格上插值；
第七步：在測(cè)試網(wǎng)格上嚴(yán)格的計(jì)算輸出角譜 ，在此例中，對(duì)應(yīng)圖6（e）；
第八步：執(zhí)行插值以獲得 ，在此例中對(duì)應(yīng)圖6（b）-6（d）；
第九步：比較嚴(yán)格仿真和插值結(jié)果，并計(jì)算相對(duì)誤差；
第十步：對(duì)于較差的插值結(jié)果，如6（b）和6（c），表2中的0-3行，其結(jié)果是σ>σ_0，嚴(yán)格的結(jié)果 將會(huì)傳遞給下一個(gè)循環(huán)并用于輸入；否則，返回當(dāng)前的結(jié)果。
圖6.算法1中不同步驟時(shí)沿κ_x方向的一維提取結(jié)果：（a）輸入角譜振幅，（b）-（d），在測(cè)試循環(huán)中的插值角譜振幅以及（e）在最后的循環(huán)中嚴(yán)格地計(jì)算輸出角譜。所有的子圖中的值都縮放到相同的范圍內(nèi)。
從表2中我們也可以看到測(cè)試首先是沿y方向，之后再沿x方向執(zhí)行的，如第三節(jié)的最后所提到的。從第0輪到第4輪測(cè)試，采樣距離δκx并沒(méi)有改變，因此采樣點(diǎn)Mx^，保持不變；在第四輪測(cè)試后σ<σ0，沿y方向的測(cè)試終止，意味著場(chǎng)數(shù)據(jù)已經(jīng)可以從先前一輪的結(jié)果中恢復(fù)。因此，沿y方向的采樣點(diǎn)數(shù)是705，在第四輪中額外的704個(gè)數(shù)據(jù)僅僅是用于測(cè)試目的，對(duì)最終的輸出并沒(méi)有貢獻(xiàn)。然后以一種類似的方式沿x方向開(kāi)始測(cè)試并在第11輪終止，同樣，場(chǎng)數(shù)據(jù)可以從先前一輪的結(jié)果中恢復(fù)。因此，最終輸出的采樣點(diǎn)數(shù)固定在2817x705。在表2的測(cè)試輪中，包含的總的采樣點(diǎn)數(shù)是5633x705+45x704，數(shù)據(jù)45x704來(lái)自于沿y方向的最后測(cè)試輪。再次強(qiáng)調(diào)一下，為了在角譜域進(jìn)行合適的采樣，有必要采用如此大的采樣點(diǎn)數(shù)。并且，除了沿x方向或者y方向最后一輪測(cè)試，用于測(cè)試目的而所有嚴(yán)格計(jì)算地值都用于構(gòu)建最終輸出場(chǎng)。
表2 Etalon模擬中每個(gè)測(cè)試輪次中的采樣參數(shù)和誤差
文章中的模擬，使用的是一臺(tái)Intel Core i7-4910MQ處理器,2.9GHz，32Gb物理內(nèi)存的電腦。表2中所顯示的示例所需的總的計(jì)算時(shí)間是118s。注意計(jì)算時(shí)間最多的部分是花在S矩陣計(jì)算上，而最后步驟的傅里葉逆變換僅花了約0.5s。
從表2中也可以看出沿x方向和沿y方向的測(cè)試輪次數(shù)也不一樣，即沿兩個(gè)方向所需的采樣間距不同。因此，我們?cè)谒惴ㄖ懈鼉A向沿每個(gè)方向分別執(zhí)行采樣測(cè)試。
通過(guò)這個(gè)算法，我們獲得了標(biāo)準(zhǔn)具的透過(guò)場(chǎng)，如圖7所示。除了圖7中透射場(chǎng)的尺寸遠(yuǎn)大于圖5中的輸入場(chǎng)，非零輸出場(chǎng)分量  也值得注意，盡管其振幅要遠(yuǎn)小于  。Ex和Ey之間出現(xiàn)的偏振串?dāng)_是由于系數(shù)矩陣T或者R的非對(duì)角形式造成的。從 的振幅分布中，還可以看出精細(xì)的同心圓，這是由于標(biāo)準(zhǔn)具多層結(jié)構(gòu)中發(fā)生的多重反射造成的。
圖7.以任意單位給出的標(biāo)準(zhǔn)具輸出場(chǎng)分量的振幅，（a）：|E_x^out |，（b）：|E_y^out |
B.單軸晶體板
我們的方法對(duì)各向同性介質(zhì)層和各向異性介質(zhì)層都同樣適用。在這個(gè)案例中，我們演示了聚光干涉儀的原理，其可以用于精確的測(cè)量單軸晶體的光軸傾斜角度[43]。當(dāng)使用會(huì)聚單色光照射正交偏振器間的晶體時(shí)，即可產(chǎn)生聚光干涉。圖8中展示了聚光干涉儀的簡(jiǎn)要原理。在我們的模擬中，光源是一個(gè)Ex線性偏振，波長(zhǎng)為633nm，NA=0.25的球面波，在晶體板前20mm。光源將輸入場(chǎng)  傳遞到晶體板的前表面，后表面的輸出場(chǎng)將被分析。在這個(gè)子節(jié)中，我們沒(méi)有討論晶體板元件外部區(qū)域的傳播步驟。
圖8.使用聚光干涉儀測(cè)量光軸的角度。偏振片P1沿x方向，偏振片P2沿y方向，光軸表示為“o.a”。
我們使用聚光干涉儀，測(cè)試了一個(gè)6mm厚的向列液晶板，其光軸與z軸成θ角。其中n_e=1.7，n_0=1.5[43]，設(shè)θ=0°，2°以及5°，我們觀察了晶體板后y方向的輸出場(chǎng)，如圖9所示。當(dāng)θ=0°，晶體的光軸垂直于晶體板表面，同心環(huán)的原點(diǎn)位于中心，如圖9（a）所示；當(dāng)θ≠0°時(shí)，如此例中的圖9（b）和
圖9（c），同心環(huán)的原點(diǎn)相應(yīng)的發(fā)生了橫向偏移。我們?cè)谀M中獲得的值與參考文獻(xiàn)[43]中分析的結(jié)果一致。
圖9.針對(duì)不同的θ值，在晶體板后表面的輸出場(chǎng)分布（振幅平方）。使用沿y方向的偏振片P2分析輸出場(chǎng)，因此僅顯示了E_y分量。紅色的十字標(biāo)出了同心圓環(huán)的原點(diǎn)，為（0，0），（210um，0）以及（525um，0），分別對(duì)應(yīng)子圖（a），（b）和（c）。
為了獲得圖9中的透射場(chǎng)，算法1再次被用于確定采樣距離。為了完整性，表3中給出了θ=5°時(shí)此例中采樣距離的數(shù)值測(cè)試的細(xì)節(jié)。
在表3所顯示的測(cè)試輪次中，包含了3969x1985+63x1984個(gè)采樣點(diǎn)，對(duì)于每個(gè)采樣點(diǎn)，計(jì)算了各向異性S矩陣，耗時(shí)416s。
表3.光軸為5°的單軸晶體板在每個(gè)測(cè)試輪次中的采樣距離和偏差
C.傾斜的各向異性晶體板系統(tǒng)
如第一節(jié)所述，我們開(kāi)發(fā)的算法是系統(tǒng)建模的一部分。在這個(gè)例子中，我們模擬了一個(gè)包含方解石晶體板的光學(xué)系統(tǒng)以研究偏振轉(zhuǎn)換以及渦旋的生成[44]，如圖10所示。因?yàn)槭轻槍?duì)整個(gè)系統(tǒng)進(jìn)行建模，因此在處理過(guò)程使用了幾種不同的場(chǎng)追跡算子[27]：偏振片使用瓊斯矩陣算子，透鏡被當(dāng)作是理想透鏡，到傾斜方解石板表面的光場(chǎng)傳輸使用參考文獻(xiàn)[33]中的技術(shù)，方解石板中的光場(chǎng)傳輸則使用了此文中的方法。
根據(jù)參考文獻(xiàn)[44]，給出了方解石晶體折射率，其中
其中λ單位為微米。
圖10.包含傾斜晶體板的系統(tǒng)。晶體板材料為方解石，厚度為6mm，光軸在圖中以“o.a.”表示，垂直于晶體板表面。偏振片P1用于生成沿x方向的偏振光，而偏振片P2用于分析。透鏡L1和L2焦距都為30mm。
此例中我們從α=0°開(kāi)始，在系統(tǒng)中幾個(gè)位置的場(chǎng)分布如圖11所示。在偏振片P1前的輸入場(chǎng)為波長(zhǎng)633nm，x方向和y方向束腰半徑1.5mm的高斯光束。經(jīng)過(guò)偏振片P1之后，獲得了一個(gè)沿x方向的線性偏振場(chǎng)。透鏡L1后的聚焦場(chǎng)是方解石晶體板的輸入場(chǎng)。根據(jù)算法1，我們獲得了晶體板后表面的輸出場(chǎng)，如圖11（a）和11（b）。最后，在偏振片P2之后的平面上，我們獲得了系統(tǒng)的輸出場(chǎng)，如圖11（c）和11（d）所示。
圖11.圖10中方解石晶體板及P2偏振片后的場(chǎng)分布（振幅平方），此例中α=0°。
讓我們來(lái)仔細(xì)地觀察一下圖11中所展現(xiàn)的結(jié)果。數(shù)值測(cè)試細(xì)節(jié)在表4中給出。與4A部分的情況相似，由于T或者R的非對(duì)角線形式，沿x方向的線性偏振輸入場(chǎng)為輸出場(chǎng)增加了一個(gè)非零y偏振方向。但此例中單軸晶體板的偏振串?dāng)_相對(duì)更強(qiáng)，E(x )≈E(y ),而在標(biāo)準(zhǔn)具的例子中，E(x )和E(y )的比例大概是1000：1，如圖7（a）和7（b）。這是由于方解石晶體的雙折射性質(zhì)，導(dǎo)致E(x )和E(y )分量之間具有一個(gè)強(qiáng)的耦合。
表4 1.5mm高斯光束入射到單軸晶體板，在每個(gè)測(cè)試輪次中的采樣距離和偏差
P2后的場(chǎng)，如圖11（e）和11（f）所示，與參考文獻(xiàn)[44]中的測(cè)量結(jié)果吻合的很好。除了由雙折射引起的四極結(jié)構(gòu)，我們同樣也看到了疊加在其上的同心圓環(huán)結(jié)構(gòu)。與圖11的分析類似，這些精細(xì)結(jié)構(gòu)是由于板內(nèi)的多重反射造成的�？梢栽趨⒖嘉墨I(xiàn)[44]中的圖4中看到同樣的影響，盡管實(shí)驗(yàn)測(cè)量的對(duì)比度較低。
下一步，如參考文獻(xiàn)[44]中所展示，我們將晶體傾斜一定的角度α并將輸入高斯束腰設(shè)置為0.5mm，以生成單極光渦旋。模擬中，我們使用[33]中的方法處理非平行板之間的傳輸。通過(guò)這種方式，可以獲得傾斜晶體板表面的輸入場(chǎng)，且傾斜板表面的輸出場(chǎng)可以繼續(xù)傳輸至下一個(gè)元件。模擬結(jié)果如圖12和13所示。當(dāng)α=1.2°時(shí)，可以在表5中獲得采樣距離的數(shù)值測(cè)試細(xì)節(jié)�？梢钥吹�，沿x方向和y方向的偏差σ分別在輪次2和5中滿足停止判據(jù)。我們沿兩個(gè)方向繼續(xù)往前迭代一次（標(biāo)記有*）以獲得一個(gè)過(guò)采樣，以為了更好的與之后的非平行平面間的傳播相結(jié)合[33]，這之間的傳播通常需要過(guò)采樣因子為2.
圖12.當(dāng)α=0.8°時(shí)，方解石晶體板和偏振片P2之后的輸出場(chǎng)（振幅平方）
圖12和13再次與參考文獻(xiàn)[44]中的結(jié)果一致。與之前分析α=0°時(shí)的情況類似，光渦旋結(jié)構(gòu)是由方解石晶體的雙折射產(chǎn)生的。如圖13（c）和13（d）所示，同樣也可以看到相位錯(cuò)位引起的光渦旋。除了渦旋，圖13（b）中還可以看到環(huán)形的精細(xì)結(jié)構(gòu)，這是由于晶體板內(nèi)的多重反射造成的。
將表5中的數(shù)據(jù)測(cè)試與之前的進(jìn)行比較，我們可以看到一個(gè)更快的收斂，特別是沿x方向的測(cè)試。包括表5中的所有測(cè)試循環(huán)，經(jīng)過(guò)傾斜的方解石晶體板的傳輸僅使用了大約4s的計(jì)算時(shí)間就完成了，而整個(gè)系統(tǒng)的模擬，包含透鏡，大約是30s的時(shí)間。由于方解石板繞y軸傾斜1.2°，若將晶體板表面作為參考，則對(duì)應(yīng)著一個(gè)輕微的傾斜入射場(chǎng)。根據(jù)傅里葉變換偏移理論，沿y方向的傾斜將導(dǎo)致角譜沿此方向有一個(gè)偏移。因此，透射場(chǎng)的角譜同樣存在沿x方向和y方向的不同調(diào)制。
圖13.當(dāng)α=1.2°時(shí)，方解石晶體板和偏振片P2后的輸出場(chǎng)，其中（a）和（b）為振幅平方的分布，（c）和（d）為相位分布
表5 0.5mm高斯光束入射到傾斜晶體板時(shí)，每輪測(cè)試中的采樣參數(shù)和偏差
D.在空氣-單軸表面的反射
對(duì)于反射情況，同樣能以類似于算法中的透射情況的方法進(jìn)行處理，我們將在下文展示。我們選用文獻(xiàn)[45]的工作，即空氣-單軸晶體表面發(fā)生的自旋霍爾效應(yīng)，來(lái)作為一個(gè)例子。系統(tǒng)原理圖如圖14所示，且自旋相關(guān)納米量級(jí)偏移可以通過(guò)檢查反射場(chǎng)的中心來(lái)進(jìn)行測(cè)量。
圖14.在空氣-單軸晶體表面的反射。晶體元件使用  坐標(biāo)系統(tǒng)來(lái)表示，其光軸（o.a.）沿  方向，輸入和反射場(chǎng)分在  中給出
在我們的模擬中，使用了波長(zhǎng)為632.8nm的線性偏振高斯場(chǎng)作為輸入場(chǎng)，在此波長(zhǎng)下，對(duì)于o光和e光，LiNbO3晶體樣品的折射率分別為n0=2.232和ne=2.156。在文獻(xiàn)[45]中并沒(méi)有給出輸出光場(chǎng)的光束尺寸，我們將其設(shè)置為20um。由于反射場(chǎng)的橫向偏移是由空間頻率域中的局部線性相位決定的，因此不改變光束尺寸的大小。
為了研究此效應(yīng)，反射矢量場(chǎng)投影在圓偏單位矢量坐標(biāo)系，而后我們測(cè)量了右旋圓偏振分量的中心。此外，對(duì)Ex和Ey兩個(gè)線性偏振入射場(chǎng)，我們旋轉(zhuǎn)晶體元件并監(jiān)控位移隨著入射角度的變化。模擬結(jié)果如圖15中所示。
我們使用參數(shù)掃描獲得圖15中的結(jié)果。對(duì)于每個(gè)入射角度，我們使用我們的算法計(jì)算了反射場(chǎng)，且在一個(gè)大的角度范圍內(nèi)顯示出了較好的適應(yīng)性。模擬結(jié)果與[45]中的結(jié)果匹配的很好。
圖15.入射角度與右圓偏振量反射場(chǎng)的橫向偏移之間的關(guān)系，（a）Ex和（b）Ey分別為線性偏振輸入。
總結(jié)
我們從系統(tǒng)仿真的角度研究了電磁場(chǎng)經(jīng)過(guò)各向同性或各向異性層介質(zhì)元件的傳輸。利用SPW分析并將以方程（10）的形式概括的結(jié)果，以為后續(xù)的數(shù)值實(shí)現(xiàn)做準(zhǔn)備。我們清楚地討論了在很多實(shí)際問(wèn)題中會(huì)遇到的采樣問(wèn)題，此外，我們提出了一種數(shù)值測(cè)試算法，以更有效地在角譜域確定采樣參數(shù)。通過(guò)標(biāo)準(zhǔn)具的案例，我們?cè)敿?xì)地描述了我們的算法的工作流程，顯示了我們方法的一般有效性。此方法已經(jīng)成功地應(yīng)用到了對(duì)激光晶體封裝技術(shù)中的應(yīng)力誘導(dǎo)雙折射的分析[46]。
致謝
我們感謝Olga Baladron-Zorita女士對(duì)此文章的校正以及其日常的幫助。
參考文獻(xiàn)
B. R. Horowitz and T. Tamir, “Lateral displacement of a light beam at a dielectric interface*,” J. Opt. Soc. Am. 61, 586–594(1971).
M. McGuirk and C. K. Carniglia, “An angular spectrum representation approach to the Goos–Hänchen shift,” J. Opt. Soc. Am. 67, 103–107(1977).
S. Kozaki and H. Sakurai, “Characteristics of a Gaussian beam at a dielectric interface,” J. Opt. Soc. Am. 68, 508–514 (1978).
C. C. Chan and T. Tamir, “Angular shift of a Gaussian beam reflected near the Brewster angle,” Opt. Lett. 10, 378–380 (1985).
F. I. Baida, D. V. Labeke, and J.-M. Vigoureux, “Numerical study of the displacement of a three-dimensional Gaussian beam transmitted at total internal reflection near-field applications,” J. Opt. Soc. Am. A17, 858–866 (2000).
S. Zhang, D. Asoubar, F. Wyrowski, and M. Kuhn, “Efficient and rigorous propagation of harmonic fields through plane interfaces,”Proc. SPIE 8429, 84290J (2012).
M. Tanaka, K. Tanaka, and O. Fukumitsu, “Transmission and reflection of a Gaussian beam at oblique incidence on a dielectric slab,” J. Opt. Soc. Am. 67, 819–825 (1977).
S. Kozaki and H. Harada, “Beam displacement of a reflected beam at an interface between an inhomogeneous medium and free space,”J. Opt. Soc. Am. 68, 1592–1596 (1978).
C. W. Hsue and T. Tamir, “Lateral displacement and distortion of beams incident upon a transmitting-layer configuration,” J. Opt.Soc. Am. A 2, 978–988 (1985).
R. P. Riesz and R. Simon, “Reflection of a Gaussian beam from a dielectric slab,” J. Opt. Soc. Am. A 2, 1809–1817 (1985).
T. Tamir, “Nonspecular phenomena in beam fields reflected by multilayered media,” J. Opt. Soc. Am. A 3, 558–565 (1986).
J. J. Stamnes and D. Jiang, “Focusing of electromagnetic waves into a uniaxial crystal,” Opt. Commun. 150, 251–262 (1998).
J. J. Stamnes and V. Dhayalan, “Transmission of a two-dimensional Gaussian beam into a uniaxial crystal,” J. Opt. Soc. Am. A 18,1662–1669 (2001).
J. J. Stamnes and V. Dhayalan, “Double refraction of a Gaussian beam into a uniaxial crystal,” J. Opt. Soc. Am. A 29, 486–497(2012).
L. I. Perez, “Reflection and non-specular effects of 2d Gaussian beams in interfaces between isotropic and uniaxial anisotropic media,” J. Mod. Opt. 47, 1645–1658 (2000).
L. I. Perez, “Nonspecular transverse effects of polarized and unpolarized symmetric beams in isotropic-uniaxial interfaces,” J. Opt. Soc.Am. A 20, 741–752 (2003).
S. Stallinga, “Axial birefringence in high-numerical-aperture optical systems and the light distribution close to focus,” J. Opt. Soc. Am.A 18, 2846–2859 (2001).
S. Stallinga, “Light distribution close to focus in biaxially birefringent media,” J. Opt. Soc. Am. A 21, 1785–1798 (2004).
M. Jain, J. Lotsberg, J. Stamnes, and Ø. Frette, “Effects of aperture size on focusing of electromagnetic waves into a biaxial crystal,” Opt.Commun. 266, 438–447 (2006).
M. Jain, J. K. Lotsberg, J. J. Stamnes, Ø. Frette, D. Velauthapillai, D.Jiang, and X. Zhao, “Numerical and experimental results for focusing of three-dimensional electromagnetic waves into uniaxial crystals,”J. Opt. Soc. Am. A 26, 691–698 (2009).
S. Zhang and F. Wyrowski, “Simulations of general electromagnetic fields propagation through optically anisotropic media,” Proc. SPIE 9630, 96300A (2015).
D. Asoubar, S. Zhang, and F. Wyrowski, “Simulation of birefringenceeffects on the dominant transversal laser resonator mode caused by anisotropic crystals,” Opt. Express 23, 13848–13865(2015).
G. D. Landry and T. A. Maldonado, “Gaussian beam transmission and reflection from a general anisotropic multilayer structure,” Appl. Opt.35, 5870–5879 (1996).
V. Dhayalan and J. J. Stamnes, “Focusing of electromagnetic waves into a dielectric slab, I: exact and asymptotic results,” Pure Appl. Opt.7, 33–52 (1998).
V. Dhayalan and J. Stamnes, “Focusing of electromagnetic waves into a dielectric slab, II: numerical results,” J. Eur. Opt. Soc. 6,11036 (2011).
A. Turpin, Y. V. Loiko, T. K. Kalkandjiev, and J. Mompart, “Light propagation in biaxial crystals,” J. Opt. 17, 065603 (2015).
F. Wyrowski and M. Kuhn, “Introduction to field tracing,” J. Mod. Opt.58, 449–466 (2011).
D. Asoubar, S. Zhang, F. Wyrowski, and M. Kuhn, “Parabasal field decomposition and its application to non-paraxial propagation,”Opt. Express 20, 23502–23517 (2012).
D. Asoubar, S. Zhang, F. Wyrowski, and M. Kuhn, “Efficient semianalytical propagation techniques for electromagnetic fields,”J. Opt. Soc. Am. A 31, 591–602 (2014).
H. Zhong, S. Zhang, and F. Wyrowski, “Parabasal thin-element approximation approach for the analysis of microstructured interfaces and freeform surfaces,” J. Opt. Soc. Am. A 32, 124–129 (2015).
F. Wyrowski, H. Zhong, S. Zhang, and C. Hellmann, “Approximate solution of Maxwell’s equations by geometrical optics,” Proc. SPIE9630, 963009 (2015).
F. Wyrowski and C. Hellmann, “Geometrical optics reloaded,”Optik Photonik 10, 43–47 (2015).
S. Zhang, D. Asoubar, C. Hellmann, and F. Wyrowski, “Propagation of electromagnetic fields between non-parallel planes: a fully vectorial formulation and an efficient implementation,” Appl. Opt. 55, 529–538 (2016).
34. J. J. Stamnes, B. Spjelkavik, and H. M. Pedersen, “Evaluation of diffraction integrals using local phase and amplitude approximations,”Opt. Acta 30, 207–222 (1983).
D. W. Berreman, “Optics in stratified and anisotropic media:4 × 4-matrix formulation,” J. Opt. Soc. Am. 62, 502–510 (1972).
L. Li, “Reformulation of the Fourier modal method for surface-relief gratings made with anisotropic materials,” J. Mod. Opt. 45,1313–1334 (1998).
L. Li, “Note on the s-matrix propagation algorithm,” J. Opt. Soc. Am. A20, 655–660 (2003).
Z. Wang, “Analytical handling of optical wavefront,” Master’s thesis(Friedrich Schiller University, 2015).
L. Rabiner, R. Schafer, and C. Rader, “The chirp z-transform algorithm,” IEEE Trans. Audio Electroacoust. 17, 86–92 (1969).
J. L. Bakx, “Efficient computation of optical disk readout by use of the chirp z transform,” Appl. Opt. 41, 4897–4903 (2002).
M. Leutenegger, R. Rao, R. A. Leitgeb, and T. Lasser, “Fast focus field calculations,” Opt. Express 14, 11277–11291 (2006).
Physical optics design software, “Wyrowski virtuallab fusion,” developed by Wyrowski Photonics UG, distributed by LightTrans GmbH. http://www.lighttrans.com.
B. L. V. Horn and H. H. Winter, “Analysis of the conoscopic measurement for uniaxial liquid-crystal tilt angles,” Appl. Opt. 40, 2089–2094(2001).
Y. Izdebskaya, E. Brasselet, V. Shvedov, A. Desyatnikov, W.Krolikowski, and Y. Kivshar, “Dynamics of linear polarizationconversion in uniaxial crystals,” Opt. Express 17, 18196–18208(2009).
Y. Qin, Y. Li, X. Feng, Z. Liu, H. He, Y. Xiao, and Q. Gong, “Spin Hall effect of reflected light at the air-uniaxial crystal interface,” Opt.Express 18, 16832–16839 (2010).
P. Ribes-Pleguezuelo, S. Zhang, E. Beckert, R. Eberhardt, F.Wyrowski, and A. Tünnermann, “Method to simulate and analyse induced stresses for laser crystal packaging technologies,” Opt.Express 25, 5927–5940 (2017).
|
