MÔ PHỎNG BA CHIỀU LINH KIỆN BÁN DẪN<br />
SỬ DỤNG THUẬT TOÁN BICGSTAB TIỀN ĐIỀU KIỆN<br />
VỚI TIỀN ĐIỀU KIỆN JACOBI<br />
CHO LỜI GIẢI PHƯƠNG TRÌNH POISSON<br />
ĐINH NHƯ THẢO<br />
Trường Đại học Sư phạm - Đại học Huế<br />
NGUYỄN TIẾN NGỌC<br />
Trường THPT Tây Trà - Quảng Ngãi<br />
Tóm tắt: Bài báo này trình bày các kết quả mô phỏng ba chiều động lực<br />
học hạt tải trong đi-ốt p-i-n bán dẫn GaAs bằng phương pháp Monte<br />
Carlo tập hợp tự hợp sử dụng tiền điều kiện Jacobi vào thuật toán<br />
BICGSTAB tiền điều kiện để giải phương trình Poisson. Chúng tôi đã<br />
tính toán vận tốc trôi dạt của điện tử cũng như sự phân bố điện thế<br />
trong linh kiện ứng với điện trường ngoài là 100 kV/cm và 150 kV/cm.<br />
Kết quả tính toán chỉ ra rằng thuật toán này không những cho kết quả<br />
tốt mà còn rút ngắn thời gian tính toán so với thuật toán BICGSTAB.<br />
<br />
1. GIỚI THIỆU<br />
Việc nghiên cứu các linh kiện bán dẫn na-nô có vai trò rất lớn để cải thiện hiệu quả<br />
và tốc độ làm việc của các sản phẩm công nghệ cao. Ngày nay, với sự phát triển của<br />
phương pháp mô phỏng, việc nghiên cứu các linh kiện bán dẫn này trở nên dễ dàng<br />
và khả thi hơn so với phương pháp giải tích. Trong nhiều phương pháp mô phỏng<br />
đã được sử dụng, phương pháp Monte Carlo tập hợp tự hợp với ưu điểm là tính ổn<br />
định và độ chính xác cao nên sớm được quan tâm và sử dụng rộng rãi [1], [2], [4], [5].<br />
Mô phỏng Monte Carlo tập hợp tự hợp ba chiều đã được thực hiện thành công trên<br />
nhiều linh kiện bán dẫn [2], [3]. Gần đây nhất, năm 2009, Lê Hoài Linh và Trần<br />
Thiện Lân đã mô phỏng khá thành công đi-ốt p-i-n bán dẫn GaAs. Tuy nhiên, việc<br />
sử dụng thuật toán BICGSTAB để giải phương trình Poisson trong các công trình<br />
này có nhược điểm là thời gian tính toán dài, tốc độ hội tụ chậm và trong nhiều<br />
Tạp chí Khoa học và Giáo dục, Trường Đại học Sư phạm Huế<br />
ISSN 1859-1612, Số 04(16)/2010: tr. 34-41<br />
<br />
MÔ PHỎNG BA CHIỀU LINH KIỆN BÁN DẪN SỬ DỤNG THUẬT TOÁN...<br />
<br />
35<br />
<br />
trường hợp nghiệm không hội tụ. Để khắc phục những nhược điểm của thuật toán<br />
BICGSTAB, người ta đã đưa ra thuật toán BICGSTAB tiền điều kiện [3], [5]. Trong<br />
các tiền điều kiện đã được sử dụng thì tiền điều kiện Jacobi có ưu điểm là tốc độ<br />
hội tụ nhanh và sử dụng đơn giản hơn [5]. Thuật toán này được kì vọng sẽ cải thiện<br />
đáng kể tốc độ hội tụ và rút ngắn thời gian mô phỏng. Chúng tôi đã xây dựng một<br />
chương trình mô phỏng chạy trên máy tính cá nhân có cấu hình Intel(R) Pentium(R)<br />
Dual, CPU T2330 1.60GHz, 1GB DDR II với sai số 10−12 . Bộ công cụ mới này cho<br />
kết quả tốt và thời gian mô phỏng cũng được rút ngắn đáng kể so với thuật toán<br />
BICGSTAB.<br />
2. MÔ PHỎNG MONTE CARLO TẬP HỢP TỰ HỢP BA CHIỀU<br />
Phương pháp mô phỏng Monte Carlo là phương pháp lặp đi lặp lại việc tính toán<br />
kiểu không tất định (nghĩa là kết quả thu được là như nhau bất kể chúng ta lặp lại<br />
việc tính toán bao nhiêu lần) và dùng một tập số ngẫu nhiên như các dữ liệu đầu<br />
vào. Phương pháp Monte Carlo tập hợp tự hợp ba chiều là sự kết hợp đồng thời<br />
phương pháp Monte Carlo tập hợp với việc giải phương trình Poisson trong ba chiều<br />
không gian.<br />
Trong quá trình mô phỏng, việc giải phương trình Poisson ba chiều là công việc quan<br />
trọng nhất và cũng khó khăn nhất. Có nhiều phương pháp để giải phương trình<br />
Poisson, sau đây chúng tôi sẽ trình bày phương pháp sai phân hữu hạn. Phương<br />
trình Poisson cho vật liệu là đồng nhất có dạng:<br />
ρ<br />
∇2 ϕ = − ,<br />
s<br />
<br />
(1)<br />
<br />
với ϕ là điện thế, ρ là mật độ điện tích, s là hằng số điện môi của vật liệu.<br />
Trong trường hợp ba chiều phương trình trên được viết tường minh như sau:<br />
∂ 2ϕ ∂ 2ϕ ∂ 2ϕ<br />
ρ<br />
+<br />
+<br />
=<br />
−<br />
.<br />
∂x2<br />
∂y 2<br />
∂z 2<br />
s<br />
<br />
(2)<br />
<br />
Sau khi thực hiện phép sai phân hữu hạn, ta nhận được hệ phương trình sau (với i,<br />
j, k là các chỉ số chạy)<br />
ϕi−1,j,k − 2ϕi,j,k + ϕi+1,j,k ϕi,j−1,k − 2ϕi,j,k + ϕi,j+1,k<br />
+<br />
∆x2<br />
∆y 2<br />
ϕi,j,k−1 − 2ϕi,j,k + ϕi,j,k+1<br />
ρi,j,k<br />
=−<br />
.<br />
+<br />
2<br />
∆z<br />
s<br />
<br />
(3)<br />
<br />
36<br />
<br />
ĐINH NHƯ THẢO - NGUYỄN TIẾN NGỌC<br />
<br />
Nếu ta chọn kích thước không gian mắt lưới đồng nhất: ∆x = ∆y = ∆z thì hệ<br />
phương trình (3) trở thành:<br />
ϕi−1,j,k + ϕi,j−1,k + ϕi,j,k−1 − 6ϕi,j,k + ϕi,j,k+1 + ϕi+1,j,k + ϕi,j+1,k = −<br />
<br />
ρi,j,k<br />
∆x2 . (4)<br />
s<br />
<br />
Hệ phương trình (4) là một hệ phương trình đại số tuyến tính và có thể được viết<br />
lại dưới dạng một phương trình ma trận:<br />
Aϕ = b,<br />
<br />
(5)<br />
<br />
trong đó A là ma trận vuông cấp M N P × M N P , B và ϕ là các ma trận cột cấp<br />
M N P với M , N , P lần lượt là số điểm lưới theo các phương x, y, z.<br />
Như vậy, bằng phương pháp sai phân hữu hạn chúng ta đã biến đổi phương trình<br />
Poisson ba chiều thành một hệ phương trình đại số tuyến tính. Các phương trình này<br />
sẽ được giải bằng thuật toán BICGSTAB tiền điều kiện với tiền điều kiện Jacobi.<br />
Thuật toán BICGSTAB tiền điều kiện được phát triển từ thuật toán BICGSTAB<br />
nhằm tăng cường tốc độ hội tụ, theo đó giảm thời gian tính toán. Thuật toán này<br />
gồm các bước sau [1], [7]:<br />
e1 = e<br />
1. Chọn x<br />
ˆ1 , s1 bất kỳ; tính e<br />
r1 = Aˆ<br />
x1 − b và đặt u<br />
r1 , i = 1.<br />
2. Bắt đầu vòng lặp: trong khi chưa hội tụ nghiệm.<br />
ei<br />
- Tính u<br />
ˆ từ Pˆ<br />
u=u<br />
e i = Aˆ<br />
- Tính w<br />
u và αi =<br />
<br />
ri<br />
sT<br />
1e<br />
Te<br />
s1 wi<br />
<br />
ei<br />
- Tính e<br />
zi = e<br />
ri − αi w<br />
ei = Aˆ<br />
- Tính ˆ<br />
z từ Pˆ<br />
z=e<br />
zi và q<br />
z<br />
T<br />
e<br />
ei<br />
z q<br />
ei<br />
- Tính ωi = iT và e<br />
ri+1 = e<br />
zi − ωi q<br />
ei q<br />
ei<br />
q<br />
- Tính βi+1 =<br />
<br />
αi s T<br />
ri+1<br />
1e<br />
Te<br />
ωi s1 ri<br />
<br />
e i+1 = rei+1 + βi+1 (e<br />
e i)<br />
- Tính u<br />
ui − ωi w<br />
- Tính x<br />
ˆi+1 = x<br />
ˆi − (αi u<br />
ˆ + ωi ˆ<br />
z)<br />
- Đặt i = i + 1<br />
3. Kết thúc vòng lặp.<br />
4. Kết thúc thuật toán BICGSTAB tiền điều kiện.<br />
<br />
MÔ PHỎNG BA CHIỀU LINH KIỆN BÁN DẪN SỬ DỤNG THUẬT TOÁN...<br />
<br />
37<br />
<br />
So với thuật toán BICGSTAB thì thuật toán BICGSTAB tiền điều kiện có chứa<br />
thêm một số bước tính toán làm cho thời gian tính toán của mỗi vòng lặp tăng lên<br />
nhưng số vòng lặp cũng được giảm đáng kể. Kết quả là tổng thời gian tính toán<br />
được rút ngắn so với thuật toán BICGSTAB, đó là ưu điểm chính của thuật toán<br />
BICGSTAB tiền điều kiện.<br />
Ma trận P trong thuật toán BICGSTAB tiền điều kiện ở trên được gọi là ma trận<br />
tiền điều kiện. Ma trận tiền điều kiện P được chọn sao cho xấp xỉ với ma trận A<br />
hay P −1 A xấp xỉ ma trận đơn vị. Tuy nhiên, người ta chưa đưa ra một lý thuyết<br />
chung nào cho việc xây dựng một ma trận tiền điều kiện tối ưu [7]. Có rất nhiều kĩ<br />
thuật xây dựng ma trận tiền điều kiện, các tiền điều kiện được sử dụng phổ biến<br />
là: Jacobi, SSOR, ILU, MILU, AMG [3], [4], [9]. Hiệu quả của chúng phụ thuộc vào<br />
dạng cụ thể của A và chỉ có thể kiểm tra bằng kết quả áp dụng. Trong nhiều kĩ thuật<br />
kể trên, tiền điều kiện Jacobi (còn gọi là tiền điều kiện đường chéo) có ưu điểm là<br />
việc xây dựng và áp dụng đơn giản, thuận tiện. Cụ thể tiền điều kiện Jacobi được<br />
xây dựng như sau [8]:<br />
(<br />
Aij nếu i bằng j<br />
(6)<br />
Pij =<br />
0 nếu i khác j<br />
Trong đó Pij và Aij lần lượt là phần tử thứ ij của ma trận P và A. Nghĩa là ma<br />
trận tiền điều kiện P ở đây sẽ là ma trận đường chéo của A.<br />
3. KẾT QUẢ MÔ PHỎNG VÀ THẢO LUẬN<br />
Trong bài báo này chúng tôi sử dụng mô hình<br />
đi-ốt p-i-n bán dẫn GaAs được chỉ ra ở Hình<br />
1. Mật độ điện tử và lỗ trống thuần lần lượt<br />
100nm<br />
là 1011 cm−3 và 3 × 1014 cm−3 , mật độ pha tạp<br />
phần p và n tương ứng là 0.5 × 1017 cm−3 và<br />
2.5 × 1017 cm−3 . ở đây, chúng tôi sử dụng mô<br />
100nm p<br />
i<br />
n<br />
hình ba thung lũng là thung lũng Γ, thung<br />
lũng L và thung lũng X. Bằng phương pháp<br />
440nm<br />
Monte Carlo tập hợp tự hợp ba chiều, chúng<br />
Hình 1. Mô hình đi-ốt p-i-n GaAs.<br />
tôi tiến hành mô phỏng động lực học của 580<br />
hạt tải trong đi-ốt cho hai trường hợp cường độ điện trường bằng 100 kV /cm và<br />
150 kV /cm, mật độ hạt tải được kích thích quang khoảng 1017 cm−3 .<br />
Sau khi thực hiện mô phỏng động lực học ba chiều của hạt tải trong đi-ốt p-i-n bán<br />
dẫn GaAs với cường độ điện trường lần lượt là 100 kV /cm và 150 kV /cm, chúng tôi<br />
thu được vận tốc trôi dạt của điện tử và phân bố điện thế trong mô hình linh kiện.<br />
Các kết quả này được thể hiện trên các đồ thị tương ứng.<br />
<br />
38<br />
<br />
Hình 2. Vận tốc trôi dạt toàn phần và theo phương<br />
x, y, z của điện tử theo thời gian ứng với<br />
E = 100 kV /cm.<br />
<br />
ĐINH NHƯ THẢO - NGUYỄN TIẾN NGỌC<br />
<br />
Hình 3.Vận tốc trôi dạt toàn phần của điện tử theo<br />
thời gian ứng với E = 100 kV /cm, 150 kV /cm.<br />
<br />
Đồ thị vận tốc trong Hình 2 cho thấy vận tốc các hạt tăng nhanh vượt qua giá trị<br />
bão hòa và sau đó mới giảm mạnh về giá trị bão hòa. Hiện tượng này được gọi là<br />
hiện tượng vượt quá vận tốc. Vì trường ngoài đặt theo phương x nên trên Hình 2, ta<br />
cũng thấy rằng vận tốc trôi dạt của điện tử theo phương y và vận tốc trôi dạt của<br />
điện tử theo phương z biến thiên nhỏ. Vận tốc trôi dạt của điện tử theo hai phương<br />
này ảnh hưởng không nhiều đến vận tốc toàn phần nên vận tốc toàn phần gần như<br />
bằng vận tốc theo phương x.<br />
Hình 3 cho thấy vận tốc trôi dạt toàn<br />
phần của điện tử ứng với các giá<br />
trị điện trường ngoài là 100 kV /cm<br />
và 150 kV /cm. Chúng tôi thấy rằng<br />
ứng với điện trường càng cao thì hiện<br />
tượng vượt quá vận tốc của các hạt<br />
xảy ra càng nhanh và sau đó giảm<br />
nhanh về giá trị bão hòa. Nguyên<br />
nhân của hiện tượng này là do sự dịch<br />
chuyển liên thung lũng của điện tử từ<br />
thung lũng Γ đến thung lũng L. Khi<br />
điện trường càng cao thì các điện tử<br />
được gia tốc càng mạnh nên sẽ có rất<br />
nhiều điện tử tham gia vào quá trình<br />
dịch chuyển liên thung lũng này. Bởi<br />
vì khối lượng hiệu dụng của điện tử<br />
<br />
Hình 6. Giản đồ điện thế 3D ứng với cường độ điện trường<br />
bằng 100 kV /cm và 150 kV /cm.<br />
<br />