光柵衍射問題數(shù)值計(jì)算的最小二乘方法

欒 天，王玉潔，鄭恩希，3

（1.北華大學(xué) 數(shù)學(xué)與統(tǒng)計(jì)學(xué)院，吉林 吉林132033；2.吉林大學(xué) 數(shù)學(xué)研究所，長(zhǎng)春130012；3.大連海事大學(xué) 理學(xué)院，遼寧 大連116026）

1 光柵衍射問題

在電磁場(chǎng)散射理論中，散射體是周期結(jié)構(gòu)表面的散射問題稱為光柵衍射［1］.光柵衍射在微光學(xué)領(lǐng)域（如光譜分析、納米尺度光學(xué)原件設(shè)計(jì)和光纖通信）應(yīng)用廣泛，關(guān)于這類問題的數(shù)值研究目前已取得了一些結(jié)果［2－4］，通常采用的數(shù)值方法為有限元方法［5－7］，適用于一般形狀的光柵結(jié)構(gòu).但當(dāng)問題的波數(shù)較大時(shí)，有限元方法將帶來(lái)較大的計(jì)算量［8］，不再適合實(shí)際應(yīng)用.為此，本文針對(duì)光柵衍射問題，提出一種更有效的算法——最小二乘方法，該方法不僅適用于一般形狀的衍射光柵，而且能克服大波數(shù)帶來(lái)的困難，應(yīng)用過程簡(jiǎn)單，所需剖分單元少，收斂速度快。

考慮一維光柵的時(shí)諧衍射問題.設(shè)光柵的表面為Γ，周期為d.帶狀區(qū)域S定義為

設(shè)Γ0?S為Γ在帶狀區(qū)域S內(nèi)的一個(gè)周期，Γ0將S分割為兩個(gè)部分，上半部分記為E.入射平面波為uI＝exp｛iαx－iβy｝，其中α＝ksinθ，β＝kcosθ，k為波數(shù)，－π／2＜θ＜π／2為入射角.在區(qū)域E 中，全場(chǎng)u滿足Helmholtz方程：Δu＋k2u＝0.考慮滿足擬周期條件的解，即u滿足

這里u（x，y）exp｛－iαx｝在x方向是以d為周期的周期函數(shù).為了保證數(shù)學(xué)問題的解存在唯一，并且符合物理要求，還要求衍射場(chǎng)ud＝u－uI在無(wú)窮遠(yuǎn)處滿足有界外行平面波條件，即衍射場(chǎng)在無(wú)窮遠(yuǎn)處僅由有限多個(gè)向外傳播的平面波構(gòu)成.

一維光柵時(shí)諧衍射問題的數(shù)學(xué)模型為：給定平面入射波uI，求擬周期解u，滿足如下邊值問題：

且衍射場(chǎng)ud＝u－uI滿足有界外行平面波條件，這里a，b∈?不同時(shí)為零.當(dāng)a＝0，b≠0時(shí)，對(duì)應(yīng)Dirichlet邊界條件；當(dāng)a≠0，b＝0時(shí)，對(duì)應(yīng)Neumann邊界條件；當(dāng)a≠0，b≠0時(shí)，對(duì)應(yīng)Robin邊界條件.

2 最小二乘方法

將一個(gè)周期內(nèi)的計(jì)算區(qū)域進(jìn)行簡(jiǎn)單剖分，在每個(gè)有界子域中選用平面波函數(shù)近似解u的局部性態(tài).關(guān)于解在無(wú)窮遠(yuǎn)處的性態(tài)，使用Rayleigh展開的前有限項(xiàng)截?cái)嗳ソ疲@種近似可以自然滿足有界外行平面波條件和擬周期條件.

考慮 子 區(qū) 域 Ej（j＝1，2，…，s）.由 于Helmholtz方程的解u∈H1（Ej）可以利用平面波函數(shù)去逼近［9］，因此可采用平面波函數(shù)近似解u的局部性態(tài).

在每個(gè)Ej內(nèi)部選取點(diǎn)xj＝（xj，yj）（j＝1，2，…，s），并選取Nj個(gè)方向d（j）l（l＝1，2，…，Nj），定義局部近似空間

圖1 符號(hào)示意圖Fig.1 Diagram of notations

在子區(qū)域Es＋1中，u的 Rayleigh展式［2］為

其中：

因此，可以采用Rayleigh展開的有限項(xiàng)截?cái)嘧鳛榻鈛的近似.在Es＋1中，定義近似空間

結(jié)合上述所有近似空間Vj（j＝1，2，…，s＋1），定義試探函數(shù)空間V為

于是，可定義誤差匹配泛函為

其中［v］表示函數(shù)v在子區(qū)域相交邊界處的躍度，定義為

利用對(duì)偶技巧可以得到最小二乘方法的一個(gè)基本估計(jì)：

命題1 若k2≠（αn＋α）2，n∈?，則存在一個(gè)與u和Nj無(wú)關(guān)的常數(shù)C，使得

命題1的證明過程與文獻(xiàn)［10］中定理3.1和文獻(xiàn)［11］中定理2.1的證明完全類似，故略.命題1表明，對(duì)任意非共振波數(shù)k，J（uN）1／2控制著解的內(nèi)部誤差，因此可以用于判斷算法的收斂性.

3 數(shù)值模擬

下面通過數(shù)值模擬驗(yàn)證算法在計(jì)算光柵衍射問題時(shí)的有效性.數(shù)值實(shí)驗(yàn)均使用Matlab軟件實(shí)現(xiàn).

為簡(jiǎn)單，取Nj＝p（j＝1，2，…，s＋1），p∈?，即每個(gè)單元上選取相同數(shù)目的平面波函數(shù)，且平面波的方向按如下方式選?。篸（j）l＝（cosθl，sinθl），θl＝2π（l－1）／p，l＝1，2，…，p.

例1 直線光柵Γ＝｛（x，y）∈?2，y＝0｝，入射波為平面波uI＝exp｛iαx－iβy｝.取波數(shù)k＝50，光柵周期d＝2，入射角θ＝π／4，a＝0，即Dirichlet邊界條件.將計(jì)算區(qū)域｛（x，y）；－1＜x＜1，y＞0｝剖分，選取點(diǎn)x1＝（－0.5，0.25），x2＝（0.5，0.25），如圖2所示.

由于例1中問題存在真解u＝uI－exp｛iαx＋iβy｝，因此本文計(jì)算了真解和數(shù)值解在區(qū)間Ω＝（－1，1）×（0，1）上的L2誤差，并分析了L2誤差和泛函J（uN）1／2的收斂性.數(shù)值結(jié)果表明，誤差隨基底數(shù)目的增加而快速減少，而且L2誤差可以被J（uN）1／2所控制，如圖3所示.

圖2 直線光柵與區(qū)域剖分Fig.2 Straight line grating and domain decomposition

圖3 L2 誤差與J（uN）1／2收斂結(jié)果Fig.3 L2 error and convergence result of J（uN ）1／2

圖5為泛函J（uN）1／2的收斂性結(jié)果.由圖5可見，隨基底數(shù)目的增加，泛函值隨之快速衰減.當(dāng)基底數(shù)目達(dá)到收斂性要求時(shí)，收斂速度也很快.此外，由于L2誤差可以被泛函J（uN）1／2所控制，所以即使在沒有真解的情況下，仍可以判定算法是關(guān)于基底數(shù)目p收斂的.

圖4 曲線光柵與區(qū)域剖分Fig.4 Curve grating and domain decomposition

圖5 J（uN）1／2收斂結(jié)果Fig.5 Convergence result of J（uN ）1／2

由例1和例2可見，本文算法在處理光柵衍射問題時(shí)是高效的.一方面，算法僅需較少的剖分單元，從而減小了計(jì)算量；另一方面，當(dāng)波數(shù)較大時(shí)，算法同樣可以達(dá)到較好的精度.

［1］MENG Pin－chao，JIANG Zhi－xia， LI Yan－zhong.Electromagnetic Scattering of Diffractive Grating in Homogeneous Medium ［J］.Journal of Jilin University：Science Edition，2012，50（6）：1151－1155.（孟品超，姜志俠，李延忠.均勻介質(zhì)中衍射光柵的電磁散射 ［J］.吉林大學(xué)學(xué)報(bào)：理學(xué)版，2012，50（6）：1151－1155.）

［2］BAO Gang，Cowsar L，Masters W.Mathematical Modeling in Optical Science［M］.Philadelphia：Frontiers in Applied Mathematics，2001.

［3］YIN Wei－shi，ZHANG De－yue，MA Fu－ming.Numerical Calculation of the Scattering Problem for Grating by Integral Equation Method［J］.Journal of Jilin University：Science Edition，2009，47（6）：1112－1121.（尹偉石，張德悅，馬富明.光柵散射問題數(shù)值計(jì)算的積分方程方法 ［J］.吉林大學(xué)學(xué)報(bào)：理學(xué)版，2009，47（6）：1112－1121.）

［4］LUAN Tian，MA Fu－ming.Well－Posedness of Anisotropic Layers Scattering above Rough Surfaces［J］.Journal of Jilin University：Science Edition，2012，50（2）：213－218.（欒天，馬富明.粗糙曲面上各向異性介質(zhì)層散射問題的適定性 ［J］.吉林大學(xué)學(xué)報(bào)：理學(xué)版，2012，50（2）：213－218.）

［5］BAO Gang.Numerical Analysis of Diffraction by Periodic Structures：TM Polarization［J］.Numer Math，1996，75（1）：1－16.

［6］BAO Gang，CAO Yan－zhao，YANG Hong－tao.Numerical Solution of Diffraction Problems by a Least－Squares Finite Element Method［J］.Math Methods Appl Sci，2000，23（12）：1073－1092.

［7］BAO Gang，CHEN Zhi－ming，WU Hai－jun.Adaptive Finite－Element Method for Diffraction Gratings［J］.J Opt Soc Am A，2005，22（6）：1106－1114.

［8］Deraemaeker A，Babu?ka I，Bouillard P.Dispersion and Pollution of the FEM Solution for the Helmholtz Equation in One，Two and Three Dimensions［J］.Int J Numer Meth Eng，1999，46（4）：471－499.

［9］Moiola A，Hiptmair R，Perugia I.Plane Wave Approximation of Homogeneous Helmholtz Solutions ［J］.Z Angew Math Phys，2011，62（5）：809－837.

［10］ZHENG En－xi.Application of Least－Squares Non－polynomial Finite Element Methods in Several Scattering Problems［D］.Changchun：Jilin University，2012.（鄭恩希.非多項(xiàng)式最小二乘有限元法在幾種散射問題中的應(yīng)用 ［D］.長(zhǎng)春：吉林大學(xué)，2012.）

［11］ZHENG En－xi，MA Fu－ming，ZHANG De－yue.A Least－Squares Non－polynomial Finite Element Method for Solving the Polygonal－Line Grating Problem ［J］.Journal of Mathematical Analysis and Applications，2013，397（2）：550－560.