半導體器件模擬及數值分析第一章

引

言什么叫半導體器件計算機模擬為什么要學習本課程怎樣學好本課程本課程的內容參考書考試方式1.什么叫半導體器件計算機模擬1．1

定義半導體器件的數值模擬是指這樣的過程：對所研究的半導體器件建立或選用合適的物理模型，并對其抽象得到相應的數學表述，然后利用適當的數值方法開發計算機軟件，并賦以器件的工藝、幾何尺寸、電學方面的模型參數，進行計算機計算，得到器件的特性及其內部的物理圖象。1.2

發展半導體器件模擬的概念起源于肖克萊（Shockley）1949年發表的論文，這篇文章奠定了結型二級管和晶體管的基礎。但這是一種局部分析方法，不能分析大注入情況以及集電結的擴展。1964年，古默爾（H.K.Gummel）首先用數值方法代替解析方法模擬了一維雙極晶體管，從而使半導體器件模擬向計算機化邁進。1969年，D．P．Kennedy和R.R.O’Brien第一個用二維數值方法研究了JFET。1969年，J．W．Slotboom用二維數值方法研究了晶體管的DC特性。從此以后，大量文章報導了二維數值分析在不同情況和不同器件中的應用。相應地也有各種成熟的模擬軟件，如CADDET和MINIMOS等。1.什么叫半導體器件計算機模擬1.什么叫半導體器件計算機模擬1

.3 分類

表1

.

1 半導體器件計

算機模擬分類

類 別M

O

S

器件

雙極型

器件

Ga

As

M

E

S

F

E

T傳感器

件

及

其

它基本模型

方

程泊

松方程

、電流

連續性

方程或

玻爾茲

曼方程泊

松方程

、電流連續性

方程

泊

松

方程

、電流

連

續性

方程或

玻

爾茲

曼方程麥

克斯韋

方程組

、電流

連續性

方程以

及其它模型方

程維 數1-3

維1-3

維1-2

維1-2

維數值處理

方

法有

限元法

、有限

差分法

或蒙特卡羅法

有

限元法

、有限差分法

有

限

元法

、有限

差

分法

或蒙特卡羅法

有

限元法

、有限差分法

非

線性偏

微分方

程組求

解方法耦

合方法

、非耦

合方法

或粒子模擬方

法耦

合方法

、非耦合方法

耦

合

方法

、非耦

合

方法

或粒子模擬方

法耦

合方法

、非耦合方法

物理模型

經

典模型

或半經典模型

經典模型

經

典

模型

或半經典模型

經典模型

1.什么叫半導體器件計算機模擬器件模擬小尺寸器件大尺寸器件二維或三維模擬一維模擬二維模擬半經典模型(波松方程、玻爾茲曼方程)特征尺寸小于一定值特征尺寸大于一定值蒙特卡羅方法耦合方法有限元法有限差分法非耦合方法大電流流小電經典模型(波松方程、連續性方程)圖1.1各劃分方法之間的內在聯系2．1

器件研制的周期性和效益性（CAD）2.為什么要學習本課程用戶要求器件橫向設計器件縱向設計光刻版制備工藝設計工藝流程參數測試圖1.2

傳統的器件試制過程用戶要求器件橫向設計器件縱向設計工藝設計器件模擬參數顯示圖1.3

器件模擬得到電參數的方法工藝模擬2．1

器件研制的周期性和效益性（CAD）2.為什么要學習本課程圖1.4

器件模型化的概念2．2

器件結構的復雜性（CAA）隨著半導體工藝的飛速發展，器件尺寸越來越小，器件的二維和三維效應已非常明顯，同時各種新結構器件也不斷出現?，F在，用傳統的解析方法已不能準確分析器件的特性，而必須采用數值方法。通過數值模擬不但可以得到器件的外部特性，同時可以顯示器件內部的物理圖象，從而加深對器件的理解和分析。2.為什么要學習本課程理論學習半導體物理，半導體器件物理，量子阱，半導體激光器數值分析上機實踐FORTRAN語言，C語言，MATLAB語言理論計算與實際器件相結合3.怎樣學好本課程4．1

雙極器件的模擬一維PN結：穩態特性瞬態特性一維晶體管：穩態特性頻率特性二維晶體管：穩態特性4．2

單極型器件的模擬有限元模擬方法蒙特卡羅模擬方法4.本課程的內容4．3

光增益的計算普通量子阱應變量子阱4．4

光波導的模擬一維實數光波導一維復數光波導二維光波導4．4

半導體激光器的特性模擬穩態特性瞬態特性小信號特性何野

魏同立

半導體器件的計算機模擬方法

科學出版社趙鴻麟

半導體器件計算機模擬

天津大學出版社張義門

任建民

半導體器件計算機模擬

電子工業出版社倉田衛著

張光華譯

半導體器件的數值分析

電子工業出版社S．賽爾勃赫著

阮剛等譯

半導體器件的分析與模擬

上海科學技術文獻出版社吉利久著

計算微電子學

科學出版社葉良修著

小尺寸半導體器件的蒙特卡羅模擬

科學出版社5.參考書第二章

器件模擬的模型方程主要內容器件模擬的基本方程組方程中變量的表示式半導體器件的分級模擬在低頻或靜電場條件下，有如下本構關系：其中

和

分別是半導體材料的介電常數和磁導率。介質中的Maxwell微分方程為：其中E為電場強度，D為電位移矢量，H為磁場強度，B為磁感應強度，J為傳導電流密度，

為電荷密度。2.1

器件模擬的基本方程組半導體器件中的電磁場運動規律應遵循Maxwell方程，從Maxwell方程出發，可以導出半導體器件的基本方程。(2.1.3)(2.1.4)(2.1.2)(2.1.1)

H

J

D

J

E

-

B

D

B

0total

t

t(2.1.5)(2.1.6)D

εEB

H2.1

器件模擬的基本方程組2.1.1 Poisson方程

(2.1.5)

D

E

D

(2.1.3)

(

E

)

E

(

)

Poisson方程通式:ad

N

)在半導體中:

q(

p

n

N

,

則對同一材料的半導體器件,

Poisson方程可簡為:2

ad(

p

n

N

N

)q

2.1

器件模擬的基本方程組2.1.2

電流連續性方程

H

J

D (2.1.2)

t

H

(J

D

)

0

t(2.1.3)

D

0

t

J

上式的物理意義是單位體元中流出的電流等于單位時間內該體元中凈電荷密度的減少全電流連續性方程:2.1

器件模擬的基本方程組2.1.2電流連續性方程在半導體器件中，傳導電流密J由電子電流密度Jn和空穴電流密度Jp兩部分組成，則

t

tn p

(J

J)

q

p

q

n

0

t

tp n

J

q

p

(

J

q

n)

t

J

0上式的物理意義是：單位時間、單位體元內流出以及增加的正電荷數等于單位時間、單位體元內流入以及減少的的負電荷數。構成這一事實的原因在于在該體元內存在著增加或減少電荷的“源”或“漏”。這就是存在著產生或復合中心。=-qR2.1

器件模擬的基本方程組2.1.2電流連續性方程定義R為單位時間、單位體積內復合的電子和空穴對數，即復合率，則有，

p

t

n

tnp

J

R

J

R1q1q上面第一式的物理意義為：單位時間、單位體元內增加的空穴數等于單位時間流入該體元的空穴數減去該體元內復合掉的空穴數。2.1

器件模擬的基本方程組n p式中q是單位荷量，

、

分別是電子和空穴遷移率，n、p是電子和空穴濃度，

n、

p是電子空穴的準費米勢.

由載流子濃度的玻爾茲曼分布：載流子輸運的基本方程電流連續性方程中的載流子密度Jn和Jp應根據器件的尺寸有兩種形式常規半導體器件的載流子輸運方程J

n

q

nn

nJ

p

q

pp

p對于傳統的常規半導體器件，電子和空穴的電流密度方程為:

(

)

pi

kTqp

n exp

q

n

niexp

kT (

n)

ip

i

npqn

n

kT

qkTlnln

n

pn nn

q

ppE-

qDp

pJ

q

nE

qD

nJE

=

-

，Dn=q qkT kT

n

，Dp=

p其中,2.1

器件模擬的基本方程組2.1.3載流子輸運的基本方程2.1.3.2小尺寸半導體器件的載流子輸運方程(a) (b)圖2.1

半導體中的載流子過沖.(a)GaAs材料,(b)Si材料2.1

器件模擬的基本方程組2.1.3載流子輸運的基本方程2.1.3.2小尺寸半導體器件的載流子輸運方程當器件的特征長度達到一定值（Si器件0.1

m，GaAs器件為0.5

m）時，電流的輸運需用玻爾茲曼方程來描述，

f

t

{f(k')Pk(',k)

f(k)P(k,k')}k'

1F

f

υ

fk rh式中f是載流子分布函數，h=h/2

，h是Planck常數，P（k,k’）是單位時間載流子波矢從k變到k’的幾率，

是載流子速度（不是漂移速度）矢量，F是作用在載流子上的力。分布函數隨時間的變化外場引起的分布函數變化溫度場不同引起的分布函數變化散射引起的分布函數隨時間變化2.1

器件模擬的基本方程組2.1.3載流子輸運的基本方程2.1.3.2小尺寸半導體器件的載流子輸運方程在求得分布函數之后，可通過下式確定電流密度：其中“+”、“-”分別對應于空穴和電子，波矢為k的電子的速度由下式確定，2qJ

(x

,

t)

3

υ(k

)

f

(k

,

x

,

t)dk(2

)kh

(k

)

1

E(k

)式中E(k)代表電了能量。如果作準靜態局部電場近似，則可以將六維問題轉化為單一的三維實空間問題，半經典模型轉化為漂移擴散模型。準靜態局部電場近似的基本假設是：載流子隨電場的變化速率比電場的有效變化速率快得多，即與載流子的基本響應相比較，器件內部的介電弛豫過程和器件響應要慢得多。2.1

器件模擬的基本方程組2.1.4熱流方程用于熱計設的熱流方程為：

t

C

T

Q

[

(T)

T

]其中

是材料的密度，C是材料的比熱，T為絕對溫度，Q為單位時間單位體積內產生的熱，

(T)為熱導率。上式的物理意義是：單位時間單位體積內凈吸收的熱量等于該體積上熱流的散度。單位時間單位體積內所產生的熱來源于能量的消耗，而能量的消耗來源于兩個方面：一是電流在材料中傳輸消耗的能量，二是載流子的產生和復合所消耗的能量?？杀硎緸椋篞=（Jn+Jp）

E

+R

Eg2.1

器件模擬的基本方程組2.1.5

光子速率方程對于半導體激光器，光子速率方程可寫為，其中S是光子密度，c

是真空光速，n是材料的折射率，g是有源區的光增益，

有源區的光限制因子，N是有源區內的載流子濃度，

是自發發射光子壽命，

是自發發射因子，A是模式損耗，包括腔損耗Ai和端面損耗Am，A=

Ai+Am其中，

dt ndS

c(g

A)S

N

1 212L R

RmA

1 ln式中，L是腔長，R1、R2分別是兩個端面的能量反射率。2.1

器件模擬的基本方程組2.1.6光波導方程由Maxwell方程組同樣可以導出在半導體材料中傳輸的光波的電場分量E所滿足的方程：0

2E

n2k2E

0式中n為材料的折射率，k0

=2

/

，

是波長。對于沿z方向傳播的波，E

(

x,

y,

z,

t

)

(

Ex

,

E

y

,

Ez

)

exp j

(

t

z

)所以,

/

z

j

,

2

/

z

2

2式中

是波沿z方向的傳播常數，可得到Helmholtz方程為，02 22 2

T E

(n k

)E

02 2 2 2 2式中,

T

/

x

/

y2.1

器件模擬的基本方程組半導體器件模擬中的模型方程歸納如下:

(

)

p

n

1

J

R

t q n

p

1

J

R

t q

dt ndS

c(g

A)S

N

T E

(n k0

)E

02 2 2 2

t

C

T

Q

[

(T

)

T

]J

n

q

nnE

qDn

nJ

p

q

p

pE

-

qDp

pk'

f

1F

f

υ

f

{f(k')P(k',k)

f(k)P(k,k')}k r

t

h2qJ

(x

,

t)

3

υ(k

)

f

(k

,

x

,

t)dk(2

)常規小尺寸2.2

方程中變量的表示式2.2.1載流子濃度考慮到重摻雜情況,

引入有效本征載流子濃度,

則有

(

)(

)pienie

kT

q

p

n exp

kT

q

n

n expnie稱為有效本征載流子濃度，可用與摻雜濃度有關的經驗表達式來表示。

ie i2kTln

d a

C

n

n(T)expN

0

N

N

2

N

N

N

0

qV1

ln

d a

對Si器件:其中：V1=9

10-3V，N0=1017cm-3，C=0.5，T是絕對溫度，ni(T)是本征Si的載流子濃度。。2.2

方程中變量的表示式2.2.1載流子濃度圖2.2

室溫下,

Si材料的本征載流子濃度與摻雜濃度的關系1510501016 1017 1018 1019 1020 10212025n

/nie iD AN+N

(cm-3)當摻雜濃度大于1017cm-3時,

重摻雜效應開始出現2.2

方程中變量的表示式ppnq

p

pJ

q

nn

nJ

(

)

pie(

n)

kTqp

n exp

qn

nieexp

kT

ie

p

ie

np

qnnqkTkTlnln

n

定義有效電場強度,2.2.1載流子濃度

ieppppiennnnq

kTq

kT

ln

nn

q n

lnn

J

q

pE

-

qD

p

q

pJ

q

nE

qD

kTqkTq

lnn

ie

lnn

ieE

E

En

E

p則有,n nqD

nJ

p

q

p

pE

p

-

qDp

pJ

n

q

nnE2.2

方程中變量的表示式2.2.2遷移率定義:單位電場強度下載流子的漂移速度

d

(E)

/

E影響遷移率的主要因素有：1.

晶格振動散射

~T

-a，

2.電離雜質散射~

T

3/2/N，

3.載流子-載流子散射，

4.中性雜質散射

~極低溫下起主要作用，

5.速度飽和效應

，

6.表面散射。當考慮遷移率與雜質濃度、溫度及速率飽和效應的關系時，可有如下經驗公式:(1

(E

/

Ec)

)1/

(T

/

T )

0 0 0(1

0.8

exp(T

/

600))

0

(T

/

T

) 2.319

*Vm

Ec

1

(Nt/Nr)

0

min

max

min

其中,2.2

方程中變量的表示式2.2.2遷移率表2.1

遷移率經驗參數取值載流子類型參數電 子空 穴

min(cm2/V.s)65.047.7

max(cm2/V.s)Nr(cm-3)1330.08.5

1016495.06.3

1016

0.720.76Vm(cm/s)1.1

1079.5

106

2.422.2

2.02.0T0(K)3003001200Nt=1016cm-3Electron

Mobility(cm2/V.s)1000800Nt=1017cm-3600400Nt=1018cm-3200Nt=1019cm-30102 103 104 105 106Electric

Field(V/cm)1400Nt=1015cm-3T=300K圖2.3

溫度等于300K時，不同摻雜濃度下，電子的遷移率與電場強度的關系2.2

方程中變量的表示式2.2.2遷移率(a) (b)圖2.3電場強度等于1000V/cm時，不同摻雜濃度下，電子(a)和空穴(b)的遷移率與溫度的關系3003203804000200400600800100012001400Nt=1019cm-3Nt=1018cm-3Nt=1017cm-3Electron

Mobility(cm2/V.s)340 360Temperature(K)Nt=1015cm-3Nt=1016cm-3E=103V/cm3003203804000100200300400500Nt=1018cm-3Nt=1019cm-3Nt=1017cm-3Nt=1016cm-3Hole

Mobility(cm2/V.s)340 360Temperature(K)Nt=1015cm-3E=103V/cm2.2

方程中變量的表示式2.2.3載流子的產生和復合

在半導體中存在著五種產生和復合機構：SRH（Shockley-Read-

Hall）復合，輻射復合，俄歇復合，表面復合以及碰撞電離

1.

SRH復合間接復合機構，在復合過程中放出聲子，對間接帶隙半導體材料，這是一種主要復合機構。RSRHnp

n2

ie

p

(n

n1

)

n

(

p

p1

)n1，p1是費米能級和復合中心能級重合時導帶和價帶中的載流子濃度。對于深能級復合中心，可以取n1=p1=nie。在高摻雜下，會有和摻雜有關的復合中心能級，所以載流子的壽命就是摻雜的函數，

1

N

n

n0

1

N0n

1

0

pN

N

p

p0

1

對于Si材料的一些經驗參數取值為：

p0=

n0=0.5

10-6s，N0n=N0p=5

1016cm-32.2

方程中變量的表示式2.2.3載流子的產生和復合

在半導體中存在著五種產生和復合機構：SRH（Shockley-Read-

Hall）復合，輻射復合，俄歇復合，表面復合以及碰撞電離

2.

輻射復合這是一種直接復合，在復合過程中放出光子，在直接帶隙材料中，這是一種主要的復合機構r

為電子空穴的復合幾率3.

俄歇復合是碰撞離化的逆過程。在俄歇復合中，電子-空穴對復合，并將過剩的能量傳遞給第三個電子或空穴，所以，是一個三粒子過程。對于重摻雜半導體，俄歇復合作用比較明顯cn和cp是由實驗決定的常數，對于Si，cn=

2.8

10-31cm6.s-1，cp=9.9

10-32cm6.s-1Rsp=r

(np-nie2)RAUG=cnn(np-nie)+cpp(np-nie

)2 22.2

方程中變量的表示式半導體中載流子的產生和復合統一寫為，R=RSRH+RSP+RAUG+RSUR+

RA2.2.3載流子的產生和復合4.表面復合通過表面復合中心能級形成的復合。它對于載流子在靠近表面附近傳輸的器件中很重要的。采用和SRH復合類似的公式來計算。式中sn

和sp

分別為電子和空穴的表面復合速度，其值與表面狀況關系很大，通常在100

~

106cm/s，

(y)為

函數。5.碰撞電離這是一種在高場作用下的凈產生過程。設在強電場下的碰撞電離的凈產生率為GA，那么-RA=GA=

nn

n+

pp

p或者，1 1ssp nRSUR

1

(n

n)

1

(p

p)np

n

2

ie

(

y)qppnnAJqJ

G

l

n pE

b

m

(

)

Aexp

離化系數，2.3

半導體器件的分級模擬2.3.1問題目的提出判斷一個半導體器件模擬軟件優劣的指標是功能全、精度高、速度快和便于用戶使用。功能全主要指能處理問題面廣，便于用戶使用則主要指程序輸入參數形式簡單，并以交互或對話方式工作。實際開發半導體器件模擬軟件時要考慮這兩點，但這不是衡量半導體器件模擬方法本身優劣的指標。衡量半導體器件模擬方法優劣的指標是速度快、精度高。在半導體器件的計算機模擬中，除了從指標要求出發選取好的方法外，在給定精度的條件下，還經常使用分級模擬技術以減少計算時間和提高計算速度。2.3

半導體器件的分級模擬2.3.2

MOS器件二維穩態模擬的分級級別應用條件求解的模型方程I亞閾區，

電子濃度比溝道摻雜低

2

(二維)

n

和

p

恒定II亞閾區，

需計算電流或閾電壓

2

(二維)

J

n

0 (一維)Jp

0III正常區，

忽略碰撞離化

2

(二維)

Jn

qRn (二維)Jp

0

2

(二維

)IV正常區或雪崩倍增區，

考慮碰撞離化

Jn

q(Rn

Gn)(二維

)

Jp

q(Rp

Gp)(二維

)2.3

半導體器件的分級模擬2.3.3分級模擬的意義隨著工件條件的變化，模型方程的復雜性越來越高，相應地，模擬的復雜性也越來越高。對于復雜的模擬問題，往往需要采用分級模擬的方法，該方法包括兩點：（1）根據具體的工作條件，選用級別較低的模型方程，以在保證精度的條件下大大減少計算時間。（2）利用低一級的解作為初值。由于低一級的解是本級的很好近似，這樣做將有效減少計算時間。這種方法尤其適宜于半導體器件整個工作范圍的模擬。第三章 一維PN結的模擬2主要內容基本方程及邊界條件網格劃分及差分離散PN結穩態特性求解的非耦合法PN結穩態特性求解的耦合法變量初值的選取PN結瞬態特性的求解PN結特性的模擬結果33．1．1

器件模型圖3.1給出了PN結二極管的一維結構，在這種結構中，摻雜濃度是坐標

x

的函數，即，N(x)=Nd(x)

-

Na(x) (3.1)為了簡化分析和模擬計算，我們特作如下假設：半導體材料滿足非簡度條件，載流子濃度由玻爾茲曼條件決定。在0

~

L內半導體的溫度是一常數。雜質全部電離。3.1基本方程及邊界條件NP0L圖3.1

PN結二極管的一維結構x43．1．2

基本方程3.1基本方程及邊界條件(3.2)p

(

N

p

n)

(3.4)

J

G

R

(3.3)

t

q

xG

R

t

q

x

p 1

n

1

Jn

2

x

2

q(

3

.8

)(

3

.7

)(

3

.6

)(

3

.5

)

x

J

0

E

x

x

x

nJ

qD l

tpnpppnnn

xJ

J

J

p

J

qD

p

q

n

q

(3.9)np

n

2ip

(n

ni

)

n(

p

pi

)R

53．1．3

邊界條件在兩個歐姆接觸邊界處，載流子濃度滿足電中性和熱平衡關系N(0)+p(0)–n(0)=

0N(L)+p(L)–n(L)=

0n(0)p(0)=n

2in(L)p(L)=

ni2(3.10a)(3.10b)(3.10c)(3.10d)上述邊界條件與時間無關，通過聯立求解，即可確定在邊界處的載流子濃度。定義邊界處的準費米勢為，

n(0)=

p(0)=0

n(L)=

p(L)=

VA(3.11a)(3.11b)AV

為外加在PN結兩端的電壓3.1基本方程及邊界條件(3.12b)kT(3.12

a)kTln

q n

i

p

i

n

p

n

q ln

n

ln(3.13b)nq(3.13a)A

i

kT

p(L)

i

kT

n(0)

(L)

V

(0)

ln

q n得邊界處的電位為，63.1基本方程及邊界條件Dq2

n

kTiL

其中

是材料的介電常數。載流子濃度和摻雜濃度：以本征載流子濃度

ni

為單位。電位、費米勢：以熱電壓Vt(kT/q)為單位。其它各量：根據上述單位依次類推，見表3.1。方程的歸一化定義：在方程

(3.2)-(3.8)

中，各物理量的數值變化范圍通常很大，為了能可靠地進行半導休器件的數值模擬，保證計算過程中不至于因數據過高或過低發生“溢出”現象，需要將模型方程中的各物理量除以適當的比例常數。這個過程被稱為歸一化，相應的比例常數被稱為歸一化因子，最歸一化因子的選取有一定任意性，只要能保證模型方程成立就行。歸一化因子的選取：(1)空間坐標以本征德拜長度LD為單位：73.1基本方程及邊界條件3．1．4方程的歸一化表

3.1

歸一化因子表參數名稱歸一化參數因子Si

材料歸一化因子值位置坐標xLD3.405

10-3cm時間坐標tLD2/D04.483

10-4s靜電勢

Vt0.025875V施加電壓VAVt0.025875V電場強度ElVt/LD7.295V/cm載流子濃度n、pni1.45

1010cm-3雜質濃度N,Nd,

Nani1.45

1010cm-3電流密度J,Jn,JpqD0ni/LD1.76

10-8A/cm2產生、復合率G,

RD0ni/LD23.234

1013cm-3.s-1載流子擴散系數Dn,

DpD00.025875cm2/s載流子遷移率

n,

pD0/Vt1cm2/V.s83．1．4方程的歸一化3.歸一化方程3.1基本方程及邊界條件2p

(N

p

n)

(3.17)

2

x

J

G

R

(3.16)

t

x

G

R (3.15)

n

t

x

p

n

J

(3.22)(3.21)(3.20)(3.19)(3.18)np

1

x

J

0

tJ

J

J

x

x

x

x

x

x

ElpnppppnnnnR

(n

1)

(

p

1)p n

)

p

pJ

(

p

)

n

nJ

(n盡管式(3.15)-(3.22)與式(3.2)-(3.9)的符號相同，但它們已經歸一化了。9(3.23a)(3.23b)(3.23c)(3.23d)3.1基本方程及邊界條件(3.23e)(3.23f

)3．1．4方程的歸一化4.

歸一化邊界條件：歸一化以后的邊界條件可寫為，N(0)+p(0)–n(0)=

0N(L)+p(L)–n(L)=

0n(0)p(0)=

1n(L)p(L)=

1

(0)

ln

n(0)

(L)

VA

ln

p(L)10有限差分法基本概念離散數值分析法：對函數所在區間分離成小區間后求值，故稱離散分析過程。其中包括有限差分法和有限元法。有限差分法具有如下優點：

差分格式簡潔。典型的是三點、五點和七點差分格式。

適于處理邊界不太復雜的問題。大多數半導體器件（如MOS器件、雙極器件）的邊界是簡單的。

編制程序方便。差分格式是規格化的，對半導體器件而言，只須列出器件內部格點和邊界格點的為數有限的差分格式。3.2網格劃分及差分離散113．2．1

有限差分法基本概念2．一階導數的有限差分近似為了求解函數f(x)在區間0

x

L上的值，將待求解區域以N個格點劃分為N-1段，其中每段的長l度為

x=L/(N-1)。在格點x

處的一階偏導數可近似表示為：前差分近似：3.2網格劃分及差分離散x1x2xlxl-1 xl+1xN-1xNx0L

x

x

x

x圖3.2

區間網格化中心差分近似：(3.24c)可以證明，前差分近似和后差分近擬的誤差都與

x成正比，而中心差分近似的誤差卻與

x2成正比，因此一階導數多用中心差分來近似。顯然，

x越小，近似程度越好。dxdf f

f

x

x

l l

1 l

xdx后差分近似：

dff

ll

x

x

xf

l

1dfdx

x

xlf

l

1

f

l

12

x123．2．1

有限差分法基本概念3．二階導數的有限差分近似二階導數d2f/dx2在點

xl

處的差分近似為：3.2網格劃分及差分離散2dx

2d f

x

22fl

x

x

fl

fl

df

/

dx

df/

dx

xfl

1

fl

1

fl

1

xfl

1x

xlxl

1

/

2xl

1

/

2可以證明，上式的誤差也與

x2成正比。x1 x2xlxl-1 xl+1xN-1xNxL

x0

x

x

x圖3.2

區間網格化(3.25)133.2網格劃分及差分離散x=0M+1M+2x=L1K-1 K K+1圖3.3

中線網格劃分方案K+2

N’M-1h(M-1)Mh(M)3．2．2

PN結的網格劃分方案基于上述分析，我們用“中線網格”方案對PN結進行網格劃分，如圖3.3所示。對求解區域0

x

L

用N’個格點進行網格劃分，在K=1至K=N’的N’個格點中，每兩個格點的中點設置一個中線網格點，中線網格點從1變化到N’-1。格點的選取遵循如下原則：在摻雜濃度變化陡峭區和PN結耗盡區選擇較細的網格，而在中性區，間距可以適當地取大。根據這種網格劃分，在中線網格點上計各物理量的導數，在原始格點上計算變量n、p和

的值。這種網格嵌套技術有利于提高差分近似的精度。h’(K)143．2．3

模型方程的差分離散對于穩態情況，載流子濃度不隨時間變化，電流連續性方程寫為，3.2網格劃分及差分離散pn

x

R

0 (3

.

26 a

)

R

0 (3

.

26 b

)

x

J

J根據圖3.3的網格劃分方案，寫出穩態下電流連續性方程和泊松方程的差分離散式為：(3.27a)Jn

(M)

Jn

(M

1)

R(K)

0(3.27b)h'(K)h'(K)Jp

(M)

Jp

(M

1)

R(K)

0(3.27c)

h'(K)

h(M) h(M

1)

1

(K

1)

(K)

(K)

(K

1)

[p(K)

n(K)

N(K)]

n

[n(

K

)

1]

p

[

p(

K

)

1]21h'(K

)

[h(M

)

h(M

1)]其中, R

(

K

)

p(

K

)n(

K

)

1 17 21153.2網格劃分及差分離散(3.30)h(M) h(M)nn nJ

(M)

(M)

n(M)

(K

1)

(K)

(M)

n(K

1)

n(K)式中的n(M)必須用相鄰的原始格點K和K+1上的n值來表示，這將影響求解的精度。Scharfetter-Gummel方法：認為在格點K和K+1之間，電場、遷移率、電流密度的變化可以忽略。則：n n

x

xJ

(n

n)

0nnlnE

dxdnJ

3．2．4

電流密度方程的差分離散根據式

(3.18) ，按照通常的差分方法，中點M處的電子電流密度的差分形式為，163.2網格劃分及差分離散Jn[1

exp(

Elx

)]3．2．4

電流密度方程的差分離散n

(

x

)

n

(0

)

exp(

El

x

)

{1

exp[

(K

1)

(K)]}Jn

(M)h(M)

n[

(K

1)

(K)]n(K

1)

n(K)exp[

(K

1)

(K)]

整理上式得到電子的電流密度離散為，

1)

exp[

1)]

n(K

)h(M

)J(M)

nnexp[

(K

)

(K

1)]

1

(K

)

(K

[

(K

)

(K

1)]

(M

) n(K

Jn

[n(K

),

n(K

1),

(K

),

(K

1)]同理得空穴的電流密度離散式為，h(M

)J (M

)

pp

J

p

[

p(K

),

p(K

1),

(K

),

(K

1)]exp[

(K

)

(K

1)]

1

p(K

)

exp[

(K

)

(K

1)]

p(K

1)[

(K

)

(K

1)]

(M

)（3.34a)（3.34b)n

El重新回到

圖3.3 坐標系，以格點K為原點，得173.23．3

PN結穩態特性求解的非耦合法

(K

1)]

0F

[(n(K),p(K),(3.35a)(3.35b)(3.35c)

(K

1),

(K),

Fp[(n(K),p(K

1),p(K),p(K

1),

(K

1),

(K),

(K

1)]

0

Fn[(n(K

1),n(K),n(K

1),p(K),

(K

1),

(K),

(K

1)]

0

其中2

K

N’-1。由于在K=1和N’處的n、p和

值可由邊界條件確定，因此共有3

(N’-2)個未知量需要確定，而我們現在正好有3

(N’-2)個方程，因此解是可以確定的。由于方程的數目較多，用傳統的解析方法是無法求解的，而必須借助計算機來求解。這就需要研究與計算機相結合的數值求解方法。在本節和下節我們分別給出用非耦合法和耦合法求解方程組(3.35)的方法。將式（3.34）代入式((3.277))，則可得到3

(N’-2)個方程，21 25183．3．1

非耦合法求解的基本思想3．3

PN結穩態特性求解的非耦合法設有一組含有三個未知數的線性方程組a11x1+a12x2+a13x3=b1a21x1+a22x2+a23x3=b2a31x1+a32x2+a33x3=

b3(3.36a)(3.36b)(3.36c)設系數矩陣非奇異，且對角線元素11 22 33a 、a 、a 皆不為零，上列各式可分別改寫為，x1=(b1-a12x2-a13x3

)/a11x2=(b2-a21x1-a23x3

)/a22x3=(b3-a31x1-a32x2

)/a33(3.37a)(3.37b)(3.37c)1 2選擇一組試探性的初始解x

(0)，x

(0)，3x

(0)，x(1)=(b

-

a x

(0)-

a x(0)

)/a1 1 12

2 13

3 11x(1)=(b

-

a x

(1)-

a x(0)

)/a2 2 21

1 23

3 22x(1)=(b

-

a x

(1)-

a x(1)

)/a3 3 31

1 32

2 331 2 3 1用新得到的x

(1)，x

(1)，x

(1)置換x

(0)，2 3x

(0)，x

(0)開始新的循環。對第k次，x(k)=(b-

a1 1 12

2x(k-1)-

ax(k-1)

)/a133 11x2 =(b2-a21x1 -a23x3 )/a22(k) (k) (k-1)x3 =(b3-a31x1 -a32x2 )/a33(k) (k) (k)(3.38a)(3.38b)(3.37c)通過多次循環，直到滿足下述條件，i

1,2,

3x(k

)ix(k

)

x(k

1)(k

)M

max

i i

193．3

PN結穩態特性求解的非耦合法3．3．1

非耦合法求解的基本思想例:設方程組為4x1-x2+x3

=4x1+6x2+2x3

=9-x1-2x2+5x3

=2不難證明，這組方程組的準確解是x1=

x2=x3

=1。用上述迭代方法，令初始解為x1(0)=

x2(0)=x3(0)

=0，則因x1=(4+x2-x3)/4x2=(9-x1-2x3)/6x3=(2+x1+2x3

)/5經過5次迭代的結果如表3.2所示，可見，如果取小數點后4位，則經過五次迭代可求得其解。表 3

.

2 非耦合法解線性方程的迭代解結

果迭代次

數x

1x

2x

3000010

.

1000

10

10

.

1333

10

10

.

1

1

33

10

120

.

1050

10

10

.

9473

10

00

.

9889

10

030

.

9896

10

00

.

1005

10

10

.

9999

10

040

.

1001

10

10

.

9999

10

00

.

1000

10

150

.

1000

10

10

.

1000

10

10

.

1000

10

1203．3．2

非耦合法求解程序框圖—Gummel方法)3．3

PN結穩態特性求解的非耦合法開

始網格劃分求h(M),h’(K)賦初值n0、p0、

0解泊松方程得

解電子連續性方程

、p固定得n解空穴連續性方程

、n固定得p

—

0<

?N輸出結果下一偏壓初值

Y 圖3.4

非耦合法求解PN結特程序框圖213.2.353．3

PN結穩態特性求解的非耦合法(3.40)

aK

(K

1)

bK

(K

)

cK

(K

1)

p(K

)

n(K

)

N

(K

)

0aKKK

1 h'(K

)h(M

)b

1 [

1

1 ], ch'(K

)h(M

1) h'(K

) h(M

) h(M

1)

1 ,其中，考慮載流子濃度與電位的關系p(K)=exp[

p(K)-

(K)]，n(K)=exp[

(K)-

n(K)]后，方程(3.40)是

的非線性方程，其求解可用牛頓迭代法。泊松方程的求解非線性泊松方程離散化泊松方程

(3.27c)

的求解可以有兩種方法。（1）線性化泊松方程：在求解方程時，認為n、p是固定的，這時方程是線性的。（2）非線性泊松方程：即計入n、p與

的關系，這時的泊松方程是非線性的。非線性泊松方程加快了Gummel方法的收斂速度，因此比前者應用更廣泛。下面對這種方法進行說明。式

((33.35cc))

的具體形式為，F

[(n(K

),

p(K

),

(K

1),

(K

),

(K

1)]22（3.42）3．3．3

泊松方程的求解2．牛頓迭代法求非線性方程的思想為了求非線性方程f(x)

= 0的解，將f(x)在x0處進行Talor展開，并忽略二階以上高次項，f(x0)+df0/dx

x=

0其中

x

=

x

–

x0，導數的上標0表示在x0處導數，上式是關于

x的線性3．3

PN結穩態特性求解的非耦合法因此式(3.42)的解可近似為，x=

x0+

x之后用上式替換x0，重新由式(3.43)求

x，并由上式修下x，如此循環，直到

x小于規定的精度。df /

dx0f(x0

)方程，可很容易地求得，

x

233.43．3

PN結穩態特性求解的非耦合法上式是關于

的線性方程組，其系數矩陣是三對角形式，因此可用三對角陣追趕法求解。

由式(3.44)解得

后，由下式對初值

0

進行修正，(3.46)

0（K）=

0（K）+

(2

K

N’-1)同時也對n0、p0進行如下修正，(3.47a)(3.47b)n0(K)=n0(K)

exp[

(K)]p0(K)=p0(K)

exp[-

(K)]通過反復求解式(3.44)，直到

小于規定的精度。(3.44)0

F

0

F

0

F

0F

F

(K

1)

(K

1)

(K

)

(K

)

(K

1)

(K

1)

03．3．3

泊松方程的求解3．非線性泊松方程的牛頓迭代解將方程

(3.40)

在

0(K-1)，

0(K)和

0(K+1)處進行Talor展開，并忽略二階以上高次項，得243．3．3

泊松方程的求解3．非線性泊松方程的牛頓迭代解式(3.44)中的各常數項為，3．3

PN結穩態特性求解的非耦合法000 0 0 0 00(3.45

f

)n

(K)(3.45e)p

(K)(3.45d)(3.45c)(3.45b)(3.45a)KK

n0(K)

p0(K)

ck

p0

(K)

n0

(K)

bk

ak

(K)

(K)

F0

(K

1)

(K)

(K)

F0

(K)

F0

(K

1)F

a

(K

1)

b

(K)

c

(K

1)

p(K)

n(K)

N(K)

K25353．3

PN結穩態特性求解的非耦合法3．3．4

電流連續性方程的求解將電位和空穴濃度視為固定值時，離散電子電流連續性方程((33..35aa))重寫為，Fn

[(n(K

1),

n(K

),

n(K

1)]

0 (3.48)由于產生復合項的影響，上式是關于n的非線性方程，同樣可以通過牛頓迭代法求解，n的初值由非線性泊松方程的解根據式(3.47)得到。類似地，由空穴電流連續性方程可求得空穴的分布。26耦合法求解的基本思想基本思想非耦合的Gummel方法是建立在方程之間的弱耦合假設基礎上的。如果問題本身的非線性嚴重，方程之間的耦合強列，例如器件工作于大電流、高電壓的工作條件下或者載流子的產生-

復合項很顯著，Gummel方法收斂很慢，甚至失效。在這類問題的模擬中，應該采用耦合法。耦合法的基本思想是認為方程(3.35)之間的耦合強烈，各方程必須同時求解。3．4

PN結穩態特性求解的耦合法273．4．1

耦合法求解的基本思想耦合法的基本過程是：將F

、Fn和Fp在

0、n0和p0處進行Talor展開，忽略二階以上高次項，得到各離散網格點的線性化矩陣方程，在所有離散點上同時求解上述方程和進行牛頓迭代，直至

、

n和

p滿足收斂條件。顯然，對于N’個離散點，耦合法需要處理3N’

3N’階的系數矩陣，因此它要比非耦合法占用更多的存貯單元，且計算量大，在半導體器件的多維模擬中，有時甚至無法實現。3．4

PN結穩態特性求解的耦合法

F

F

F

Fn

n

Fn

p

Fn

Fp

n

Fp

p

Fp

n

p(3.49)0pnF

0F

0F

0

p

n

283．4．1

耦合法求解的基本思想2．耦合法和非耦合法的比較耦合法與非耦合相比，各有自已的優缺點，主要體現在：存貯單元：非耦法特別有利于節省存貯單元，對于N’個離散點，只涉及N’

N’階的系數矩陣，而非耦合法的系數矩陣為3N’

3N’階。收斂性：耦合法的收斂速度快，而非耦合法的收斂速度慢。初值的選?。厚詈戏▽Τ踔档倪x取較嚴格，而非耦合法則可在較廣的初值范圍內保證收斂。對方程的要求：非耦合法適于處理方程之間是弱耦合的情況，而當方程之間耦合強列時，則必須用耦合方法。3．4

PN結穩態特性求解的耦合法293．4．2

方程的線性化處理由于Jn、Jp的非線性，方程(3.35a)和(3.35b)也是非線性的，這對求解高階方程組(3.35)帶來了困難，所以，必須將這些方程線性化，通過解線性化方程并進行牛頓迭代來求解。將差分離散方程(3.35)在n0、p0和

0處進行Talor展開，并忽略二階以上高次項，有3．4

PN結穩態特性求解的耦合法(3.50a)n

p(K)

p(K)

Fn

n(K

1)

n(K

1)

Fn

n(K)

n(K

1)

n(K

1)

n(K)

Fn

FnF0

(K

1)

0

Fn

(K

1)

(K)

(K

1)

Fn

(K)

Fn

(K

1)30以及,3．4

PN結穩態特性求解的耦合法0(3.50c)(3.50b)F0p

p(K)

p(K)

FF

p(K

1)

p(K

1)

Fp

p(K)

p(K)

p(K

1)

Fp

p(K

1)

Fp

Fp

n(K)

n(K)

Fp

(K

1)

0

(K

1)

(K)

(K

1)

0

Fp

(K

1)

(K)

(K)

(K

1)

Fp

(K

1)

F

(K

1)

(K

1)

F

(K)

n(K)

F

n(K)

F

其中，2

K

N’-1，Fn0、Fp0和F

0分別對應于變量n0、p0和

0的Fn、Fp和F

的函數值。31進一步將上述三個方程整理成統一的矩陣形式：3．4

PN結穩態特性求解的耦合法其中矩陣下標“0”表示所有矩陣元是在n0、p0和

0處的偏導數。各矩陣元可由式(3.27)-(3.29)和式(3.34)求出

Fn

Fn

Fn

n(K)

Fp

p(K

)

Fp

(K

)

Fp

n(K)

F

p(K

)

F

(K

)

F

n(K

)

p(K

)

(K

)(3.51)00000

pn

(K

1)

F

0

p(K

1)

F

0

n(K

1)

F

0

F

(K

1)

p(K

1)00

Fp

(K

1)

Fp

(K

1)

Fn

Fn

n(K

1)

0

(K

)

p(K

)

n(K

)

(K

1)

p(K

1)

n(K

1)

F

(K

1)

Fn

(K

1)

p(K

1)00

Fp

Fn

(K

1)

n(K

1)0

Fn32例如，對泊松方程，式(3.51)的各矩陣元為，3．4

PN結穩態特性求解的耦合法11(3.52f

)(3.52e)(3.52d

)h'(K

)h(M

)(3.52c)(3.52b)(3.52a)h'(K

)h(M

)

1

n(K

)

1

1

h'(K

)

h(M

1) h(M

)

h'(K)h(M

1)1

1 1

[N

(K

)

p0

(K

)

n0

(K

)]

1 1

1 1h'(K

)h(M

1) h'(K

)

h(M

1) h(M

)

F0

p(K

)

F

F

F

(K

1)

F

(K