【正文】 ingly thought as a valuable numerical method，and will be widely applied the design of optoelectronics ponent.附錄B 外文翻譯 FDTD法在電磁場中的應(yīng)用1 簡介FDTD法是1966由Yee氏首次引入的一種求解麥克斯韋方程組的數(shù)值方法，并由Taflove將這一技術(shù)在微波領(lǐng)域中做了開創(chuàng)性的應(yīng)用。由于其簡單，健全性和多功能性，F(xiàn)DTD已經(jīng)在光頻率范圍內(nèi)得到廣泛的應(yīng)用，而且適合用這種方法分析的問題也已經(jīng)出現(xiàn)。對于平面光波導(dǎo)和梯度折射率光波導(dǎo)來說，用FDTD分析是一種尤其好的分析方法。其結(jié)構(gòu)非常簡單，是一種可以有效加以利用的模式，而且由于是平面的，在許多情況下只需要兩維分析即可。用超級計(jì)算機(jī)，可以分析該結(jié)構(gòu)中的數(shù)以百計(jì)的波長。時域有限差分法，實(shí)質(zhì)上是一個初值問題，在電磁場中允許用指明的來源來說明，在離散時間域可以用框架結(jié)構(gòu)來分析。該場的發(fā)展是由在每一個小單元中的復(fù)雜的介電常數(shù)決定的。在不同介電常數(shù)的邊界上，可以觀察到反射，折射和衍射。時間步長通常是幾個各不相同的完整的正弦波源周期實(shí)現(xiàn)的，場分量的最大值出現(xiàn)在存儲的最后一個周期的半個周期之間。通過這種方式，一個穩(wěn)定的解決辦法得以實(shí)現(xiàn)（允許波在整個模型空間傳播） 。在分析中，必須十分小心地選擇相應(yīng)的邊界條件。由于時域有限差分法是一種初值時域方法，傳播波將在格的邊界上反射，除非在場的邊界格上強(qiáng)加特殊條件，使它們不反射（ “吸收” ） 。這些吸收邊界條件的吸收場發(fā)生在邊界上，就如反射沒有發(fā)生。雖然在這個時候有沒有完善的吸收邊界算法，許多近似邊界條件算法已經(jīng)讓網(wǎng)格邊界上的反射最小。一階和二階近似差分形式已經(jīng)發(fā)展為包括Taflove和Brodwin ，Enquist和Majda，和Bayliss和Turkel所提出的邊界條件 。這可以防止網(wǎng)格邊界上的大量源反射，雖然吸收邊界中的網(wǎng)格可以進(jìn)行擴(kuò)展，進(jìn)行分析的結(jié)構(gòu)可以減少空間，提出了FDTD方法節(jié)省內(nèi)存和執(zhí)行時間。而其他的方法需要存儲和計(jì)算時間與和成正比，其中N是網(wǎng)格總數(shù)，而FDTD的計(jì)算時間僅需要與N成正比。這是一個更直接有效的時域計(jì)算方法。2公式表達(dá)用Yee氏網(wǎng)格表達(dá)的Maxwell旋度方程的有限差分表達(dá)形式。用FDTD或Yee氏網(wǎng)格計(jì)算該旋度方程的規(guī)則如下（我們假定這個媒介是無磁性的i .e.） （1）這里是磁導(dǎo)率，是介電常數(shù)，是電導(dǎo)率，在直角坐標(biāo)系中，Maxwell旋度方程表達(dá)式如下： (2a)(2b)(2c)圖 1 三維場分量模型(2d)(2e)(2f) 各組成部分E和H的位置結(jié)構(gòu)的單位網(wǎng)格這一方法的目標(biāo)是建立一個電磁波在介質(zhì)空間和時間空間進(jìn)行傳播的模型，在相應(yīng)空間網(wǎng)格上的每個單元中多次用有限差分法模擬了旋度方程。在入射波的傳播路徑是它傳播到電介質(zhì)結(jié)構(gòu)中與它所經(jīng)過的介質(zhì)表面發(fā)生作用激發(fā)了色散和衍射。當(dāng)正弦穩(wěn)態(tài)行為在每個網(wǎng)格中完成后使波傳播完成。E或H在直角坐標(biāo)系中任一分量，在時間和空間中的離散表示為（其中n為時間步長）：(3)FDTD方法利用Yee網(wǎng)格將六個E和H的分別放置在有限差分網(wǎng)格的每個節(jié)點(diǎn)上。式中，和分別是長方體網(wǎng)格沿x、y、z方向的空間步長。通過展開三個E和H部分來描述如圖1所示的有限差分近似的方法，使用中央差分近似法，在空間和時間上，利用===，就可以得到如式（4）和（5）所表示的Yee氏網(wǎng)格中的場分量，但是計(jì)算機(jī)所執(zhí)行計(jì)算是計(jì)算機(jī)網(wǎng)格中的節(jié)點(diǎn)或者是單元（i，j，k）。從式（5a）—（5c）表明了每個網(wǎng)格上的E分量在n+1上的時間步長，僅僅依賴于他的前一個（第n個）時間步長和他周圍的H分量。H場分量采用了一個類似的程序，除了磁場分量僅僅計(jì)算在1/2時間步長處。這就提供了一個通過Yee網(wǎng)格解決程序的常規(guī)的執(zhí)行方案。(4a)(4b)(4c) (5a)(5b) (5c) 為了避免可能造成數(shù)值色散和數(shù)值不穩(wěn)定性，以及誤差大的計(jì)算，對于時間和空間的離散化必須要有適當(dāng)?shù)倪x擇。在本文中，我們選擇1/20的空間步長來解決精確尺寸問題。時間不常是有網(wǎng)格大小來決定，而且必須滿足下列穩(wěn)定條件 (6)在這個空間模型中，哪里的波速度最大？這個立體的空間模型中===，通常會用到下面的關(guān)系(7)在相互作用空間中，電磁波最大。我們可以得到電磁波在邊界或其他空間位置中的相應(yīng)的。在本文中，我們以在i=1處的條件為例，能夠得到(8)這里所說的在網(wǎng)格截?cái)嗵幩玫降氖窍辔凰俣龋ǎ? 結(jié)論本文主要介紹了時域有限差分法的計(jì)算規(guī)則，并介紹了麥克斯韋方程組在電磁場領(lǐng)域中的數(shù)值計(jì)算結(jié)果。時域有限差分算法可以用來進(jìn)行電磁場中的六個場分量的全波分析。它不僅能夠適用于弱指導(dǎo)光波導(dǎo)的計(jì)算，而且更大的不同的波導(dǎo)結(jié)構(gòu)指數(shù)，并且由于它所占的存儲容量小，這種方法是一種非常有前景的寶貴的數(shù)值方法，它將廣泛應(yīng)用于光電元件的設(shè)計(jì)當(dāng)中。附錄C 部分程序(1) fortran計(jì)算程序c**********自動剖分邊界的邊界元通用程序 *********************c*********求解二維靜電場邊值問題（特性阻抗拋物線角點(diǎn)處理法）***************c**********自動剖分邊界的線性邊界通用程序 *********************c*********求解二維靜電場邊值問題（特性阻抗）(n2=3)***************c 求解任意形狀區(qū)域的電位分布及特性阻抗c******n2：閉合邊界數(shù)；l：內(nèi)部待求節(jié)點(diǎn)數(shù)；nn：1電場，2磁標(biāo)元；n0：1電場，2磁場；n1：1:Lapalace方程，2:Possion方程；c******x(2,ni):y某段單元的起點(diǎn)、終點(diǎn)坐標(biāo)；xm,ym：中點(diǎn)坐標(biāo)；n3(j):第j條閉合邊界的分段數(shù)；cx,cy：待求場點(diǎn)的場點(diǎn)坐標(biāo)；c******m(j,i):第J條第I段邊界的單元數(shù)目。x1(j,i),y1(j,i):第J條第I段邊界的起點(diǎn)坐標(biāo)。n4:計(jì)算次數(shù)。c******sx:第J條閉合邊界第I段X方向的單位長度；sy:第J條閉合邊界第I段Y方向的單位長度。x0(j,i),y0(j,i):第J條第I段邊界的圓心坐標(biāo)；c******r0(j,i):第J條閉合邊界第I段的曲率半徑；n5(j):第J條閉合邊界的最后一段線元編號；ss1=ds/dr。c******k0(j,i):1:邊界值待求，0：一類邊界條件 f0(j,i)=u的已知值，1：二類邊界條件 f0(j,i)=du/dn的已知值。c******m0(j,i):2:其他形狀，1：直線邊界，0：逆時針方向的圓弧，1：順時針方向的圓?。籲nc:0:電容，1：阻抗；nn(j):c******ll(j)：每一封閉曲線的周長；nn_4(j):每一封閉曲線的總段數(shù)；dat_l(j):每個封閉曲線的單元長；cc(j):電容。z0(j):阻抗；ccccccccc交聯(lián)聚乙烯 介電常數(shù)：ee0r= parameter(ni=310,nh=20,nnc=0,ee0r=) dimension k0(nh,nh),f0(nh,nh),m(nh,nh),x0(nh,nh),y0(nh,nh) dimension x(2,ni),y(2,ni),g(ni,ni),kode(ni),ex(ni),nn_4(ni), $h(ni,ni),fi(ni),dfi(ni),r(nh,nh),m0(nh,nh),cx(ni),cy(ni), $n3(nh),x1(nh,nh),y1(nh,nh),sol(ni),ey(ni),xm(ni),ym(ni),ll(nh), $dat_l(nh),cc(nh),z0(nh),n5(nh),nn(nh),ffi(nh) real ll,dat_l mon n,l,n2 mon uu pi= n=0 open(6,file=39。39。) read(6,*)l,n2 do 36 j=1,n2 read(6,*)n3(j) do 30 i=1,n3(j) read(6,*)m(j,i),k0(j,i),f0(j,i),x1(j,i),y1(j,i),m0(j,i),r(j,i), $x0(j,i),y0(j,i)30 n=n+m(j,i) x1(j,n3(j)+1)=x1(j,1)36 y1(j,n3(j)+1)=y1(j,1) read(6,*) read(6,*)uu do 40 k=1,l read(6,*)cx(k),cy(k)40 continue close(6) do 44 j=1,n2 nn4=0 do 41 i=1,n3(j) nn4=nn4+m(j,i)41 nn_4(j)=nn4c write(*,*)nn_4(j)44 continuec open(6,file=39。39。)c do 41 k=1,lc41 write(6,*)cx(k)c close(6)c open(6,file=39。39。)c do 42 k=1,lc42 write(6,*)cy(k)c close(6) n4=0 do 260 j=1,n2 ll(j)=0c write(*,*)39。j=39。,j,39。ll(j)=39。,ll(j) do 256 i=1,n3(j) if (m0(j,i))45,100,100c******************45語句為直線邊45 if(m0(j,i).) then sx=(x1(j,i+1)x1(j,i))/m(j,i) sy=(y1(j,i+1)y1(j,i))/m(j,i) do 50 k=1,m(j,i) n4=n4+1 x(1,n4)=x1(j,i)+sx*(k1) y(1,n4)=y1(j,i)+sy*(k1) x(2,n4)=x1(j,i)+sx*k y(2,n4)=y1(j,i)+sy*k kode(n4)=k0(j,i)50 fi(n4)=f0(j,i) ll(j)=ll(j)+sqrt((x1(j,i+1)x1(j,i))**2+(y1(j,i+1)y1(j,i))**2) goto 256 else if(m0(j,i).) then sy=(y1(j,i+1)y1(j,i))/m(j,i) do 80 k=1,m(j,i) n4=n4+1 y(1,n4)=y1(j,i)+sy*(k1) x(1,n4)=y(1,n4)**2 y(2,n4)=y1(j,i)+sy*k x(2,n4)=y(2,n4)**2 kode(n4)=k0(j,i)80 fi(n4)=f0(j,i) go to 256 end if100 if(m0(j,i).)then ll(j)=2*pi*r(j,i) if(y1(j,i+1)y1(j,i)) 110,130,120110 ss1=(x1(j,i)x0(j,i))/r(j,i) ss2=(x1(j,i+1)x0(j,i))/r(j,i) if()then ss1= else if (+)then ss1= end if if()then ss2= else if(+)then ss2= end if(2) 耦合電容計(jì)算程序程序1clc。clear all。n=::c0_1 = [,,]。 plot(n,c0_1,39。rp39。,39。linewidth39。,2,39。markerfacecolor39。,39。w39。)。 hold off。 h = legend(39。耦合電容39。,2)。 set(h,39。FontSize39。,9)。xlabel(39。列車左邊距隧道左邊距離與隧道寬度比39。,39。FontSize39。,9),ylabel(39。列車與隧道間耦合電容的值（單位：pF）39。,39。FontSize39。,9)。 程序2 clc。clear all。n = ::。c12 = [,,,]。 plot(n,c12,39。rp39。,39。linewidth39。,2,39。markerfacecolor39。,39。w39。)。 hold off。 h = legend(39。列車間耦合電容c1239。,2)。 %標(biāo)注不同的二維曲線并將標(biāo)注放在西北角處set(h,39。FontSize39。,9)。xlabel(39。兩列車高度與隧道高度比39。,39。FontSize39。,9),ylabel(39。兩列車間的耦合電容的值（單位：pF）39。,39。FontSize39。,9)。 %設(shè)置x、y坐標(biāo)軸標(biāo)