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ì)構成的光學(xué)層狀結構的傳播。我們還開(kāi)發(fā)了一種基于快速傅里葉變換技術(shù)的算法，具有數值高效的采樣規則。通過(guò)將此算法與其它系統建模技術(shù)相結合，我們展示了一些仿真案例，如經(jīng)過(guò)一個(gè)各向同性法珀標準具的光場(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.引言
由平行平面構成的光學(xué)層在光學(xué)中廣泛應用。層狀結構可以用作許多不同情況的模型，像平板和標準具�；谶@個(gè)事實(shí)，光與層狀結構相互作用的主題一直引起大家的注意并且對此已經(jīng)進(jìn)行了大量的研究。
在這類(lèi)研究中，大多數觀(guān)點(diǎn)都側重于平面波，然而僅僅少數的研究使用了平面波譜方法（SPW）來(lái)考慮一般的電磁場(chǎng)。例如，參考文獻[1-6]中研究了各向同性-各向同性的界面上，高斯光束的反射率和透射率；在參考文獻[7-11]中研究了各向同性層或者平板的情況；參考文獻[12-22]討論了各向同性-各向異性界面的情況，在參考文獻[23-26]中則討論了各向異性層或者平板的情況。
上面所提到的許多研究都用于特定的研究主題，像[1,3,5]中研究了高斯光束全內反射的橫向偏移，并且他們常常關(guān)注于具體的配置。因此，將這些方法推廣到更一般的情況的可能性受到了限制。
在這篇文章中，我們從一個(gè)更一般的觀(guān)點(diǎn)來(lái)考慮此問(wèn)題。光學(xué)層幾乎不會(huì )單獨使用；相反，他們常常是一個(gè)光學(xué)系統的一部分并且和其他的元件一起使用，如圖1中所示�；诖耸聦�(shí)，我們遵循場(chǎng)追跡的概念[27]，并使用不同的場(chǎng)追跡算子組合[28-32]，如圖1中所示，以對一個(gè)包含了層介質(zhì)元件的系統進(jìn)行物理光學(xué)模擬�？紤]到模擬是對整個(gè)系統而不是單個(gè)元件，仿真層結構必須與系統的前后部分相連接。這要求我們傳播步驟（圖1中的P）進(jìn)行適當的考慮，將前一個(gè)元件的輸出連接到當前元件的輸入，并將當前元件的輸出傳遞到下一個(gè)元件。一般情況下，這樣的傳輸步驟會(huì )出現在平行或者非平行平面之間。在參考文獻[28,29]中已經(jīng)提到了平行平面間幾種有效的傳輸方法，在參考文獻[33]中則可以找到對非平行平面間傳輸的一個(gè)詳細的討論。在這篇文章中，我們不會(huì )研究傳輸步驟，但會(huì )關(guān)注層狀結構的元件算子C。
此外，從數值計算的觀(guān)點(diǎn)出發(fā)，為了執行一個(gè)連續且有效的系統模擬，要求元件算子C
正確地處理采樣場(chǎng)數據并和其他的算子以一種統一的格式傳遞場(chǎng)數據；
優(yōu)化數值計算的效率。
考慮到上述兩個(gè)標準，我們開(kāi)發(fā)了一種具有自動(dòng)數值采樣規則的SPW方法。與之前一些利用積分方法對空間和角譜相關(guān)的傅里葉變換進(jìn)評估的研究相比（如參考文獻[23]中的二維中點(diǎn)規則和參考文獻[12-14，20，25]中的Stamnes–Spjelkavik–Pedersen方法[34]），我們使用了快速傅里葉變換（FFT）技術(shù)，此技術(shù)在大部分數值軟件包中容易訪(fǎng)問(wèn)并且效率高。再加上在角譜域中經(jīng)過(guò)深入考慮的數值采樣規則，我們的方法具有一般適用性，對層元件和入射場(chǎng)沒(méi)有任何限制。因此，此算法可以直接包含在一個(gè)物理光學(xué)系統模擬之中。
圖1.結合使用不同的場(chǎng)追蹤算子來(lái)模擬光學(xué)系統： C是元件算子，P是相鄰元件之間的傳輸算子。
2.理論
如圖2所示，層狀結構分別由兩個(gè)位于  和  的平行平面構成。  和  的區域充滿(mǎn)了復折射率為  和  的均勻各向同性介質(zhì)。參考文獻[27]中表明使用橫向分量Ex和Ey已足夠表征均勻各向同性介質(zhì)中電磁場(chǎng)了。因此，我們可以使用以下表達式來(lái)描述此問(wèn)題：
其中，分別在平面 和 處定義輸入和輸出橫向電場(chǎng)矢量，（兩者位于界面的數學(xué)位置，但總是認為在均勻介質(zhì)的一側），由下式給出
其中  。方程（1）中的元件算子是一個(gè)2x2的矩陣形式，
圖2.層狀結構分別由兩個(gè)位于 和 的平行平面構成。 和 的區域由均勻各向同性介質(zhì)填充，其折射率分別是 和 。輸出場(chǎng)和輸出場(chǎng)在層表面進(jìn)行定義，但總是在相應的各向同性介質(zhì)的一側。
在這章節，我們的目標是找到C的精確的形式，以連接層介質(zhì)元件的輸入和輸出場(chǎng)。為了研究與層結構的相互作用，我們對輸入橫向場(chǎng)分量進(jìn)行了一個(gè)傅里葉變換，并獲得了
其中  , F表示二維傅里葉變換，
 。逆傅里葉變換定義如下
方程（6）中的積分可以解釋為將  分解為具有不同橫向波矢分量κ的平面波。因此，在我們的情況下，每個(gè)輸入平面波都可以單獨處理——我們首先計算每個(gè)輸入平面波的輸出，然后進(jìn)行求和從而獲得輸出場(chǎng)。
此外，根據邊界條件對電磁場(chǎng)施加的連續性要求，可以顯示出一個(gè)給定的輸入平面波在與層結構相互作用的過(guò)程中其橫向波矢分量κ必定保持不變。同樣可以顯示出，通過(guò)疊加原理的有效性，不同的κ之間沒(méi)有耦合。因此，對于輸出角譜，我們可以寫(xiě)下
其中
公式（8）中分別是透射和反射系數矩陣。為了計算T(κ) 或者R(κ)，我們選擇使用數值穩定S矩陣方法。為了計算S矩陣，首先必須確定每個(gè)各向異性層的平面波�；谖墨I[35]中Berreman的4x4矩陣公式，Landry和Maldonado開(kāi)發(fā)并展示了一種數值友好的形式，見(jiàn)參考文獻[23]。我們采用了他們的方法，對于每一層，求解了參考文獻[23]中由方程（28）所描述的特征系統的特征值和特征向量。
不同于[23，25]中直接使用本征解來(lái)構建一個(gè)轉換矩陣，另外，我們還需要根據他們的傳輸方向整理出平面波，這是為計算S矩陣所做的一個(gè)必需的準備。為此，我們遵循[36]中4.3部分由Li所提出的標準。
然后可以應用遞歸S矩陣公式。我們在這篇文中不再重復給出已發(fā)展成熟的S矩陣方法，讀者可以參考文獻[37]中的方程（5）-（8）以獲得更多的信息。在我們的情況中，由于沒(méi)有反向傳輸輸入場(chǎng)，我們僅對正向透射或者反向反射感興趣，因此這篇文章中的矩陣系數T(κ) 和R(κ)對應于[37]中方程（5c）或者（5d）中的子矩陣T_uu 或者R_du。
一旦獲得了矩陣系數，通過(guò)方程(7)即可獲得輸出角譜。對輸出角譜進(jìn)行一個(gè)逆傅里葉變換，我們獲得了輸出橫向場(chǎng)矢量
通過(guò)聯(lián)合方程（4），（7），（9），我們可以寫(xiě)出從輸入場(chǎng)到輸出場(chǎng)整個(gè)計算流程，如果
我們以透射的情況作為例子，則
因此，方程（1）中元件算子C的精確形式如下
通過(guò)使用系數矩陣R(κ)代替 T(κ)，可以獲得反射情況的表達式。
3.算法
按照方程（10）的順序，我們可以應用一種數值算法以計算場(chǎng)經(jīng)過(guò)分層介質(zhì)元件的傳輸。讓我們從橫向輸入場(chǎng)矢量  和  ，以圖3（a）中的均勻網(wǎng)格進(jìn)行采樣。這種網(wǎng)格定義為  ,其中  和  作為索引數，δx和δy是x方向和y方向的采樣距離。初始采樣參數應該受到合適的控制以使他們符合先前算子的奈奎斯特-香農采樣定理。方程（10）中的F為聯(lián)系空間域和頻率域的算子，可以使用不同的數值方法實(shí)現，像廣泛使用的快速傅里葉變換（FFT）技術(shù)，以及包含了能夠進(jìn)一步提高數值效率的半解析傅里葉變換[38]和啁啾z變換[39-41]的更高級的方法。在此篇文章中，我們使用了FFT技術(shù)，并以此獲得輸入角譜。但我們的算法不受限于該技術(shù)。
圖3.在角譜域定義均勻采樣網(wǎng)格。（a），初始網(wǎng)格定義，5x5采樣點(diǎn)，  和  為采樣距離；（b），使用5x5采樣點(diǎn)，沿垂直方向生成一個(gè)測試網(wǎng)格， 和  =0.5 為采樣距離；（c）中使用9x5個(gè)采樣點(diǎn)定義了沿水平方向的一個(gè)測試網(wǎng)格，  =0.5 和 為采樣距離。（（b）和（c）中的實(shí)心點(diǎn)是在（a）的初始網(wǎng)格中出現的共同的采樣點(diǎn)位置，而空心圓環(huán)與初始網(wǎng)格中不一樣的采樣點(diǎn)。
其中角譜域中的均勻網(wǎng)格定義為  ，其中 和 為Kx和Ky的采樣距離。利用FFT計算的結果，在兩個(gè)域中的采樣點(diǎn)是一樣的，因此我們有  。
接下來(lái)，將輸入角譜乘以系數矩陣，然后我們可以獲得輸出角譜場(chǎng)，例如，對于透射
角譜  的采樣通過(guò)傅里葉變換關(guān)系自動(dòng)確保。然而，  的采樣沒(méi)有必要進(jìn)行確保，因為他們包含了系數和輸入角譜之間的點(diǎn)乘。這個(gè)問(wèn)題將在第四部分的案例中進(jìn)行清楚地說(shuō)明。
一般來(lái)說(shuō)，采樣距離 ， 和帶寬  必須進(jìn)行適當調整以確保正確的采樣。在我們的情況中，方程（13）中的運算沒(méi)有改變譜范圍。因此，我們僅僅只需要找到合適的采樣距離即可。此問(wèn)題并沒(méi)有解析解，因為系數矩陣T或者R是通過(guò)遞歸S矩陣方法數值上獲得的。確保 采樣的唯一方法是實(shí)行數值測試。
讓我們在角譜域定義一個(gè)測試采樣網(wǎng)格，如如3（b）或者3（c）所示。測試網(wǎng)格定義為  ,其中  為指數；為了不改變角譜范圍，要求  。與初始網(wǎng)格相比，它也需要進(jìn)行細化，這意味著(zhù)  。
一方面，我們以一種嚴謹的方式來(lái)計算測試網(wǎng)格上的輸出橫向角譜分量。由于輸入橫向角譜分量的合適采樣，根據需要可以進(jìn)行插值，例如，在測試網(wǎng)格上，我們獲得
其中“INTRPL”代表的是數值插值運算。為了在一個(gè)更精細的網(wǎng)格上嚴格地獲得輸入場(chǎng)，如方程（14）中，我們總是采用基于FFT的Sinc插值方法，把它們代入到方程（13）中，可以嚴格地獲得輸出橫向角譜量
另一方面，通過(guò)對  的插值，我們也可以獲得測試網(wǎng)格上的輸出橫向角譜量，并將插值結果表示為  。我們將此過(guò)程描述為
值得強調的是，在方程（16）中無(wú)需再使用嚴格的Sinc插值。因為這些結果僅僅用于評估一下情況。在這篇文章中，我們使用三次插值。上面的插值沒(méi)必要給出正確的輸出值，因此，在方程（16）的左邊，我們對那些量使用上標“intrpl”。下一步，我們定義
作為插值結果和嚴格結果之間的相對偏差。只要  在初始網(wǎng)格上的采樣滿(mǎn)足奈奎斯特-香農采樣定理，基于他們的插值結果就不會(huì )對嚴格的插值方法表現出一個(gè)大的偏差。在這種情況下，根據所選的插值方法，兩種結果的偏差應該在數值誤差級次，由σ_0表示。使用上述標準，通過(guò)一步步的減小采樣間距來(lái)測試場(chǎng)采樣，直至σ<σ_0。
以上的測試過(guò)程對應著(zhù)一個(gè)循環(huán)，在每個(gè)測試循環(huán)中都會(huì )執行方程（14）到（17）的計算。為了在每個(gè)測試循環(huán)中充分地使用測試數據，如方程（15）中的結果，我們總是使用圖3（b）和3（c）中的測試網(wǎng)格。這樣規則的測試網(wǎng)格定義為
對每個(gè)循環(huán)內，需要對先前的采樣距離減半。使用這個(gè)測試網(wǎng)格具有如下好處：當計算  時(shí)，僅僅需要計算那些在空心圓圈位置處的值，而實(shí)心點(diǎn)位置處的值在之前的測試循環(huán)中已經(jīng)計算過(guò)了，并且讓這些值簡(jiǎn)單地進(jìn)行傳輸即可。在這種方式下，這些基于測試目的所進(jìn)行嚴格計算的值不會(huì )被丟掉，而會(huì )作為下一次循環(huán)中的起始點(diǎn)以重新使用。因此，在一個(gè)完整的計算中所執行的是嚴格的計算，從而驗證合適的采樣，并且這些計算的值能夠有效地用于構建最終輸出場(chǎng)。我們將上面的方法總結在算法1中。
算法1：通過(guò)合適的采樣控制來(lái)進(jìn)行輸出角譜計算的數值流程
1)在初始網(wǎng)格嚴格地計算
2)初始化相對誤差值σ=+∞
3)判斷σ>σ_0， ⊳如果為真，則采樣不合適
4)沿Kx（或者Ky）二等分采樣距離
5)應用新的采樣距離定義測試網(wǎng)格 
6)在測試網(wǎng)格中對輸入角譜插值，根據方程（14）獲得 
7)根據方程（15），在測試網(wǎng)格上嚴格地計算輸出角譜 
8)根據方程（16），通過(guò)插值獲得 
9)根據方程（17），計算誤差σ
10)如果σ>σ_(0 )  ，則⊳將目前的輸出場(chǎng)設置為下一次循環(huán)的起始點(diǎn)
11)返回
為了有效地處理非對稱(chēng)情況，例如，光束在x方向和y方向有不同的發(fā)散角，算法1中的測試需要按順序沿兩個(gè)方向進(jìn)行。開(kāi)始方向的選擇是任意的，在我們的情況中，我們是沿y方向開(kāi)始測試的。
4.示例
在VirtualLab Fusion[42]軟件中，我們將第三節中提出的算法實(shí)現在“可編程元件”的編程界面中。這個(gè)元件可以與VirtualLab Fusion中其它的物理光學(xué)仿真技術(shù)進(jìn)行聯(lián)合仿真。接下來(lái)，我們展示了四個(gè)案例：前兩個(gè)主要關(guān)注元件本身并以一種嚴格的數值方式檢查此算法；后兩個(gè)案例中，元件將用于光學(xué)系統中，例如，此算法與其它仿真技術(shù)一起使用。
在進(jìn)入實(shí)際的案例之前，設置方程（17）中的迭代終止標準σ_0很重要。對于在此文中所使用的三次插值，我們預先檢查了它在一般情況下的表現并在我們的數值環(huán)境中找到了一個(gè)0.01的基準值。
A.各向同性標準具
第一個(gè)案例模擬了一個(gè)線(xiàn)性偏振高斯光束經(jīng)過(guò)一個(gè)標準具的傳輸。通過(guò)這個(gè)案例，我們將清楚的指出第三節中所說(shuō)的采樣問(wèn)題并描述了算法1的工作原理。標準具由熔融石英制成，兩側有多層薄膜，如圖4所示。關(guān)于標準具的光學(xué)參數和結構的更多信息，請見(jiàn)表1。輸入場(chǎng)為波長(cháng)633nm，x方向線(xiàn)偏光的高斯光束。其在元件的輸入平面定義為E_⊥^in,且束腰半徑為（2um，2um）。在經(jīng)過(guò)傅里葉變換后，我們獲得其角譜，同樣具有高斯輪廓，如圖5所示。按照方程（10）中的操作算子序列，輸入角譜將乘以透射或者反射率系數。我們仍以透射作為例子，并且對于線(xiàn)性偏振輸入場(chǎng)，我們使用t_xx和t_yx乘以E ̃_x^in，以獲得輸出角譜分量。
圖4.由熔融石英制成，兩側有多層膜的標準具。其結構和光學(xué)參數如表1中所示。
表1 標準具的結構和光學(xué)參數
圖5 .（a）輸入高斯場(chǎng)分量的振幅；（b）對應的角譜分量。由于輸入光場(chǎng)為沿x方向的線(xiàn)性偏振光，因此僅顯示Ex分量。
如第三節中所指出，乘積  的采樣不能自動(dòng)得到保證，此案例中將顯示該現象。標準具由于其頻率選擇功能（頻譜或角頻率）而得到廣泛的使用。在我們的案例中，角頻率選擇可以解釋為系數 txx和 tyx以一種方式調制輸入角譜，以使特定的角頻率加強而其它的減弱。這種調制可以出現在一種非常精細的頻率水平上。因此，需要使用更精細的采樣以在輸出角譜中解析這樣一個(gè)精細的調制。為了獲得需要的采樣間距，我們遵循算法1，圖6中顯示了部分結果。
通過(guò)圖6和表2，我們根據算法1中的步驟描述了工作流程，如下：
第一步：從圖6（a）中所示的輸入角譜開(kāi)始，計算各個(gè)系數并乘上元件矩陣以生成初始化的輸出角譜；
第二步：初始化相對偏差σ=+∞；
第三步：開(kāi)始測試循環(huán)；
第四&五步：將采樣距離沿κ_x或者κ_y方向減半，以定義測試網(wǎng)格 ，對應的采樣點(diǎn)如表2中所示；
第六步：對輸入角譜在測試網(wǎng)格上插值；
第七步：在測試網(wǎng)格上嚴格的計算輸出角譜 ，在此例中，對應圖6（e）；
第八步：執行插值以獲得 ，在此例中對應圖6（b）-6（d）；
第九步：比較嚴格仿真和插值結果，并計算相對誤差；
第十步：對于較差的插值結果，如6（b）和6（c），表2中的0-3行，其結果是σ>σ_0，嚴格的結果 將會(huì )傳遞給下一個(gè)循環(huán)并用于輸入；否則，返回當前的結果。
圖6.算法1中不同步驟時(shí)沿κ_x方向的一維提取結果：（a）輸入角譜振幅，（b）-（d），在測試循環(huán)中的插值角譜振幅以及（e）在最后的循環(huán)中嚴格地計算輸出角譜。所有的子圖中的值都縮放到相同的范圍內。
從表2中我們也可以看到測試首先是沿y方向，之后再沿x方向執行的，如第三節的最后所提到的。從第0輪到第4輪測試，采樣距離δκx并沒(méi)有改變，因此采樣點(diǎn)Mx^，保持不變；在第四輪測試后σ<σ0，沿y方向的測試終止，意味著(zhù)場(chǎng)數據已經(jīng)可以從先前一輪的結果中恢復。因此，沿y方向的采樣點(diǎn)數是705，在第四輪中額外的704個(gè)數據僅僅是用于測試目的，對最終的輸出并沒(méi)有貢獻。然后以一種類(lèi)似的方式沿x方向開(kāi)始測試并在第11輪終止，同樣，場(chǎng)數據可以從先前一輪的結果中恢復。因此，最終輸出的采樣點(diǎn)數固定在2817x705。在表2的測試輪中，包含的總的采樣點(diǎn)數是5633x705+45x704，數據45x704來(lái)自于沿y方向的最后測試輪。再次強調一下，為了在角譜域進(jìn)行合適的采樣，有必要采用如此大的采樣點(diǎn)數。并且，除了沿x方向或者y方向最后一輪測試，用于測試目的而所有嚴格計算地值都用于構建最終輸出場(chǎng)。
表2 Etalon模擬中每個(gè)測試輪次中的采樣參數和誤差
文章中的模擬，使用的是一臺Intel Core i7-4910MQ處理器,2.9GHz，32Gb物理內存的電腦。表2中所顯示的示例所需的總的計算時(shí)間是118s。注意計算時(shí)間最多的部分是花在S矩陣計算上，而最后步驟的傅里葉逆變換僅花了約0.5s。
從表2中也可以看出沿x方向和沿y方向的測試輪次數也不一樣，即沿兩個(gè)方向所需的采樣間距不同。因此，我們在算法中更傾向沿每個(gè)方向分別執行采樣測試。
通過(guò)這個(gè)算法，我們獲得了標準具的透過(guò)場(chǎng)，如圖7所示。除了圖7中透射場(chǎng)的尺寸遠大于圖5中的輸入場(chǎng)，非零輸出場(chǎng)分量  也值得注意，盡管其振幅要遠小于  。Ex和Ey之間出現的偏振串擾是由于系數矩陣T或者R的非對角形式造成的。從 的振幅分布中，還可以看出精細的同心圓，這是由于標準具多層結構中發(fā)生的多重反射造成的。
圖7.以任意單位給出的標準具輸出場(chǎng)分量的振幅，（a）：|E_x^out |，（b）：|E_y^out |
B.單軸晶體板
我們的方法對各向同性介質(zhì)層和各向異性介質(zhì)層都同樣適用。在這個(gè)案例中，我們演示了聚光干涉儀的原理，其可以用于精確的測量單軸晶體的光軸傾斜角度[43]。當使用會(huì )聚單色光照射正交偏振器間的晶體時(shí)，即可產(chǎn)生聚光干涉。圖8中展示了聚光干涉儀的簡(jiǎn)要原理。在我們的模擬中，光源是一個(gè)Ex線(xiàn)性偏振，波長(cháng)為633nm，NA=0.25的球面波，在晶體板前20mm。光源將輸入場(chǎng)  傳遞到晶體板的前表面，后表面的輸出場(chǎng)將被分析。在這個(gè)子節中，我們沒(méi)有討論晶體板元件外部區域的傳播步驟。
圖8.使用聚光干涉儀測量光軸的角度。偏振片P1沿x方向，偏振片P2沿y方向，光軸表示為“o.a”。
我們使用聚光干涉儀，測試了一個(gè)6mm厚的向列液晶板，其光軸與z軸成θ角。其中n_e=1.7，n_0=1.5[43]，設θ=0°，2°以及5°，我們觀(guān)察了晶體板后y方向的輸出場(chǎng)，如圖9所示。當θ=0°，晶體的光軸垂直于晶體板表面，同心環(huán)的原點(diǎn)位于中心，如圖9（a）所示；當θ≠0°時(shí)，如此例中的圖9（b）和
圖9（c），同心環(huán)的原點(diǎn)相應的發(fā)生了橫向偏移。我們在模擬中獲得的值與參考文獻[43]中分析的結果一致。
圖9.針對不同的θ值，在晶體板后表面的輸出場(chǎng)分布（振幅平方）。使用沿y方向的偏振片P2分析輸出場(chǎng)，因此僅顯示了E_y分量。紅色的十字標出了同心圓環(huán)的原點(diǎn)，為（0，0），（210um，0）以及（525um，0），分別對應子圖（a），（b）和（c）。
為了獲得圖9中的透射場(chǎng)，算法1再次被用于確定采樣距離。為了完整性，表3中給出了θ=5°時(shí)此例中采樣距離的數值測試的細節。
在表3所顯示的測試輪次中，包含了3969x1985+63x1984個(gè)采樣點(diǎn)，對于每個(gè)采樣點(diǎn)，計算了各向異性S矩陣，耗時(shí)416s。
表3.光軸為5°的單軸晶體板在每個(gè)測試輪次中的采樣距離和偏差
C.傾斜的各向異性晶體板系統
如第一節所述，我們開(kāi)發(fā)的算法是系統建模的一部分。在這個(gè)例子中，我們模擬了一個(gè)包含方解石晶體板的光學(xué)系統以研究偏振轉換以及渦旋的生成[44]，如圖10所示。因為是針對整個(gè)系統進(jìn)行建模，因此在處理過(guò)程使用了幾種不同的場(chǎng)追跡算子[27]：偏振片使用瓊斯矩陣算子，透鏡被當作是理想透鏡，到傾斜方解石板表面的光場(chǎng)傳輸使用參考文獻[33]中的技術(shù)，方解石板中的光場(chǎng)傳輸則使用了此文中的方法。
根據參考文獻[44]，給出了方解石晶體折射率，其中
其中λ單位為微米。
圖10.包含傾斜晶體板的系統。晶體板材料為方解石，厚度為6mm，光軸在圖中以“o.a.”表示，垂直于晶體板表面。偏振片P1用于生成沿x方向的偏振光，而偏振片P2用于分析。透鏡L1和L2焦距都為30mm。
此例中我們從α=0°開(kāi)始，在系統中幾個(gè)位置的場(chǎng)分布如圖11所示。在偏振片P1前的輸入場(chǎng)為波長(cháng)633nm，x方向和y方向束腰半徑1.5mm的高斯光束。經(jīng)過(guò)偏振片P1之后，獲得了一個(gè)沿x方向的線(xiàn)性偏振場(chǎng)。透鏡L1后的聚焦場(chǎng)是方解石晶體板的輸入場(chǎng)。根據算法1，我們獲得了晶體板后表面的輸出場(chǎng)，如圖11（a）和11（b）。最后，在偏振片P2之后的平面上，我們獲得了系統的輸出場(chǎng)，如圖11（c）和11（d）所示。
圖11.圖10中方解石晶體板及P2偏振片后的場(chǎng)分布（振幅平方），此例中α=0°。
讓我們來(lái)仔細地觀(guān)察一下圖11中所展現的結果。數值測試細節在表4中給出。與4A部分的情況相似，由于T或者R的非對角線(xiàn)形式，沿x方向的線(xiàn)性偏振輸入場(chǎng)為輸出場(chǎng)增加了一個(gè)非零y偏振方向。但此例中單軸晶體板的偏振串擾相對更強，E(x )≈E(y ),而在標準具的例子中，E(x )和E(y )的比例大概是1000：1，如圖7（a）和7（b）。這是由于方解石晶體的雙折射性質(zhì)，導致E(x )和E(y )分量之間具有一個(gè)強的耦合。
表4 1.5mm高斯光束入射到單軸晶體板，在每個(gè)測試輪次中的采樣距離和偏差
P2后的場(chǎng)，如圖11（e）和11（f）所示，與參考文獻[44]中的測量結果吻合的很好。除了由雙折射引起的四極結構，我們同樣也看到了疊加在其上的同心圓環(huán)結構。與圖11的分析類(lèi)似，這些精細結構是由于板內的多重反射造成的�？梢栽趨⒖嘉墨I[44]中的圖4中看到同樣的影響，盡管實(shí)驗測量的對比度較低。
下一步，如參考文獻[44]中所展示，我們將晶體傾斜一定的角度α并將輸入高斯束腰設置為0.5mm，以生成單極光渦旋。模擬中，我們使用[33]中的方法處理非平行板之間的傳輸。通過(guò)這種方式，可以獲得傾斜晶體板表面的輸入場(chǎng)，且傾斜板表面的輸出場(chǎng)可以繼續傳輸至下一個(gè)元件。模擬結果如圖12和13所示。當α=1.2°時(shí)，可以在表5中獲得采樣距離的數值測試細節�？梢钥吹�，沿x方向和y方向的偏差σ分別在輪次2和5中滿(mǎn)足停止判據。我們沿兩個(gè)方向繼續往前迭代一次（標記有*）以獲得一個(gè)過(guò)采樣，以為了更好的與之后的非平行平面間的傳播相結合[33]，這之間的傳播通常需要過(guò)采樣因子為2.
圖12.當α=0.8°時(shí)，方解石晶體板和偏振片P2之后的輸出場(chǎng)（振幅平方）
圖12和13再次與參考文獻[44]中的結果一致。與之前分析α=0°時(shí)的情況類(lèi)似，光渦旋結構是由方解石晶體的雙折射產(chǎn)生的。如圖13（c）和13（d）所示，同樣也可以看到相位錯位引起的光渦旋。除了渦旋，圖13（b）中還可以看到環(huán)形的精細結構，這是由于晶體板內的多重反射造成的。
將表5中的數據測試與之前的進(jìn)行比較，我們可以看到一個(gè)更快的收斂，特別是沿x方向的測試。包括表5中的所有測試循環(huán)，經(jīng)過(guò)傾斜的方解石晶體板的傳輸僅使用了大約4s的計算時(shí)間就完成了，而整個(gè)系統的模擬，包含透鏡，大約是30s的時(shí)間。由于方解石板繞y軸傾斜1.2°，若將晶體板表面作為參考，則對應著(zhù)一個(gè)輕微的傾斜入射場(chǎng)。根據傅里葉變換偏移理論，沿y方向的傾斜將導致角譜沿此方向有一個(gè)偏移。因此，透射場(chǎng)的角譜同樣存在沿x方向和y方向的不同調制。
圖13.當α=1.2°時(shí)，方解石晶體板和偏振片P2后的輸出場(chǎng)，其中（a）和（b）為振幅平方的分布，（c）和（d）為相位分布
表5 0.5mm高斯光束入射到傾斜晶體板時(shí)，每輪測試中的采樣參數和偏差
D.在空氣-單軸表面的反射
對于反射情況，同樣能以類(lèi)似于算法中的透射情況的方法進(jìn)行處理，我們將在下文展示。我們選用文獻[45]的工作，即空氣-單軸晶體表面發(fā)生的自旋霍爾效應，來(lái)作為一個(gè)例子。系統原理圖如圖14所示，且自旋相關(guān)納米量級偏移可以通過(guò)檢查反射場(chǎng)的中心來(lái)進(jìn)行測量。
圖14.在空氣-單軸晶體表面的反射。晶體元件使用  坐標系統來(lái)表示，其光軸（o.a.）沿  方向，輸入和反射場(chǎng)分在  中給出
在我們的模擬中，使用了波長(cháng)為632.8nm的線(xiàn)性偏振高斯場(chǎng)作為輸入場(chǎng)，在此波長(cháng)下，對于o光和e光，LiNbO3晶體樣品的折射率分別為n0=2.232和ne=2.156。在文獻[45]中并沒(méi)有給出輸出光場(chǎng)的光束尺寸，我們將其設置為20um。由于反射場(chǎng)的橫向偏移是由空間頻率域中的局部線(xiàn)性相位決定的，因此不改變光束尺寸的大小。
為了研究此效應，反射矢量場(chǎng)投影在圓偏單位矢量坐標系，而后我們測量了右旋圓偏振分量的中心。此外，對Ex和Ey兩個(gè)線(xiàn)性偏振入射場(chǎng)，我們旋轉晶體元件并監控位移隨著(zhù)入射角度的變化。模擬結果如圖15中所示。
我們使用參數掃描獲得圖15中的結果。對于每個(gè)入射角度，我們使用我們的算法計算了反射場(chǎng)，且在一個(gè)大的角度范圍內顯示出了較好的適應性。模擬結果與[45]中的結果匹配的很好。
圖15.入射角度與右圓偏振量反射場(chǎng)的橫向偏移之間的關(guān)系，（a）Ex和（b）Ey分別為線(xiàn)性偏振輸入。
總結
我們從系統仿真的角度研究了電磁場(chǎng)經(jīng)過(guò)各向同性或各向異性層介質(zhì)元件的傳輸。利用SPW分析并將以方程（10）的形式概括的結果，以為后續的數值實(shí)現做準備。我們清楚地討論了在很多實(shí)際問(wèn)題中會(huì )遇到的采樣問(wèn)題，此外，我們提出了一種數值測試算法，以更有效地在角譜域確定采樣參數。通過(guò)標準具的案例，我們詳細地描述了我們的算法的工作流程，顯示了我們方法的一般有效性。此方法已經(jīng)成功地應用到了對激光晶體封裝技術(shù)中的應力誘導雙折射的分析[46]。
致謝
我們感謝Olga Baladron-Zorita女士對此文章的校正以及其日常的幫助。
參考文獻
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).
|
