InAs/GaAs量子點生長的KMC模擬

陳 龍,徐凱宇

(1.上海大學上海市應用數學和力學研究所,上海200072;2.上海大學理學院,上海200444)

量子點結構是指在空間3個維度內尺寸均小于電子波長(載流子在空間各方向上均受到限制)、具有量子限域效應的3維量子受限系統[1].量子點是由半導體材料制成的納米粒子,具有許多獨特的納米性質.由于量子點具有明顯的量子約束、隧穿、干涉等光、電學效應,使其在微電子、光電子器件、超大規模集成電路以及量子計算機等方面具有較大的潛在應用優勢[2-3].目前,利用分子束外延技術[4],基于Stranski-Krastanov(S-K)生長模式[5-6]的自組織生長是制備量子點的一種有效途徑,制備高度有序的量子點以及研究量子點的定位生長成為研究熱點.

制備空間排列有序、量子點密度高的手段有很多,而在應變弛豫圖形襯底上生長量子點的方法具有可控性好、雜質較少等優點[7-8],從而受到越來越多的關注.有2種方法可以得到應變弛豫的圖形襯底:①通過物理化學腐蝕在襯底表面形成溝槽;②采用埋藏量子點方法,從而在襯底表面產生應變弛豫圖形.在GaAs襯底上生長InAs量子點系統中,一般會采用第2種方法,即在GaAs襯底中預先埋藏1層InAs量子點方法,而后在此襯底上進行自組織外延生長,此時島成核位置會更傾向于埋藏量子點的正上方,從而得到有序排列的島[9].

半導體量子點的生長通常采用的是自組織外延生長方式,其生長過程是一個復雜的熱動力學過程,因此影響量子點生長的因素比較多,主要包括生長過程中的溫度、沉積速率、應變能等[10].在實際制備中量子點的生長溫度不宜太低,因為過低的生長溫度有可能導致生成物為多晶和非晶結構[11];其次,在利用外延技術制備InAs/GaAs系統量子點的工藝過程中,原子沉積到襯底上的速率不同將對生長制備結果造成影響[12].由于沉積速率不一樣,單個原子沉積到襯底表面的時間也就不一樣,襯底上原子的擴散時間也不一樣,沉積速率越低,擴散的時間也越充足.因此,擴散是否充分對生長過程顯得猶為重要.另外,從熱力學角度來看,量子點生長過程是系統能量達到最小化的過程,也是表面能與應變能相互作用和相互競爭的過程,故應變能的分布是影響半導體量子點生長的因素之一.研究結果表明,在量子點的生長過程中原子趨向于向應變能較高的位置聚集,應變能不同的分布形式將決定量子點的成型位置與形狀[13],因而通過改變應變能分布可以有效控制量子點的定位生長.

考慮到量子點生長過程是一個隨時間變化的隨機過程,本工作運用動力學蒙特卡羅(kinetic Monte Carlo,KMC)方法,通過定義能量勢壘來模擬量子點早期生長階段單原子層成核(2D島形成)過程[14],分析量子點的生長溫度、沉積速率和埋層深度對InAs/GaAs量子點系統生長的影響.該階段將直接影響在S-K生長模式下量子點的最終位置和大小.

1 計算模型以及溫度和沉積速率對量子點生長的影響

1.1 KMC模型的建立

為了模擬量子點生長早期第1層單原子層階段,本工作將采用全擴散SOS(solid-onsolid)模型的KMC模型,此模型中將襯底劃分成離散的網格,每個沉積原子只可占據1個格點.該模型中沉積過程是任意的,也就是說原子沉積位置是隨機的,而表面擴散是該模型的主要機制,此過程由原子躍遷概率控制(躍遷概率可由Arrhenius公式得到).KMC模擬過程如下:首先,沉積原子按照沉積速率隨機沉積到襯底(襯底大小為LxLy,其中Lx為長度,Ly為寬度)上任意位置,如果該位置已有原子占據則重新選擇1個沉積位置;然后計算出每個位置原子相對應的躍遷概率p和總躍遷概率,由此確定原子是否躍遷以及躍遷方向.取總躍遷概率的倒數為每一步的時間步長,當時間滿足(式中,F為沉積速率,即單位時間內沉積到襯底表面的原子層數;t為沉積時間)時,新的原子開始沉積,重復以上步驟直至全部原子完成沉積(其模擬過程流程圖見圖1).

由于模擬過程為早期單原子層階段,故模擬對象都是小覆蓋率條件下的生長過程.表面覆蓋率

圖1 KMC模擬過程流程圖Fig.1 Flow chart of KMC simulation process

原子在襯底表面的運動為熱遷移過程,在本模型中原子每次只能躍遷到鄰近的空位,如果鄰近位置都已有原子則不躍遷.如果原子要從一個點躍遷到另一個點,由于這個過程包含了一個簡單的態的轉變,故其躍遷就需要克服相應的能量勢壘.由躍遷理論可知,在躍遷過程中所需克服的能量勢壘只取決于其初始位置,而與躍遷方向和終點位置無關.襯底表面的原子是否躍遷由其躍遷概率(只與初始位置相關)決定,躍遷概率滿足Arrhenius公式:

式中,v0為表面原子的嘗試躍遷頻率(約為1013s?1),kB為玻爾茲曼常數,T為襯底溫度(生長溫度),E為在躍遷過程中所需克服的能量勢壘,

式中,Es為原子與襯底表面的結合能,通常為一固定常數;En為與周圍原子之間的束縛能之和,

式中,n(n≤4)為周圍最鄰近位置的原子數,m(m≤4)為次近鄰位置上的原子數.在InAs/GaAs系統中,E0=0.3 eV,表示目標原子與最近鄰單個原子之間的束縛能[15];而Estr(x,y)為圖形化襯底表面的應變能,應變能的分布與位置坐標相關.

1.2 溫度對量子點生長的影響

對于InAs/GaAs系統而言,常見的生長制備得到的量子點通常為金字塔形,過高的溫度會出現塌陷現象,量子點形狀則變成圓屋頂形.另外,從熱力學角度看,當溫度過高時原子的熱運動比較活躍,擴散比較快速,島間距離比較分散,島的密度會減小,不利于控制;但溫度太低又不利于原子的結晶成核.同時,在實驗中也可以觀察到,小幅提高生長溫度可獲得更寬的尺寸分布和較少缺陷的量子點材料;而過高的生長溫度將會導致量子點面密度降低[16].因此,本工作將通過計算機模擬研究溫度對量子點生長的早期階段的影響.

為了研究GaAs應變弛豫圖形襯底上生長InAs量子點過程中溫度的影響,本工作模擬了溫度從350 K到700 K之間變化(每次溫度變化間隔為50 K)的量子點前期單原子層生長情況.模擬中選取的模擬參數[7,15,17]為選取的圖形襯底為正方形的漸變應變弛豫襯底,圖形襯底大小為100×100網格化格點.根據應變能的分布襯底又劃分成10×10正方形小區域,每個區域中應變能分布均為正四棱錐形式,在中間處取得最大值Estr2,在邊界處取最小值Estr1.生長時間t=1.5 s,沉積速率F=0.1 mL/s,表面覆蓋率為15%,即模擬中沉積原子數為1 500個,能量勢壘取值為Estr1=0.35 eV,Estr2=0.55 eV,Es=0.3 eV.

圖2顯示了在不同溫度條件下襯底上島的生長形態.從圖中可以很直觀地看到,在溫度很低時島的分布比較凌亂,且島的尺寸也較小;隨著溫度的升高,島的尺寸逐步增大,逐漸出現定位清晰整齊的較大島;當溫度達到450 K時,可以看到原子分布依舊比較混亂,但已經有原子開始聚集在一起;當溫度升高至600 K時,島的分布已經比較均勻,在每個網格之內基本上都有且僅有一個島,各個位置的島的大小也比較一致;當溫度達到650 K或更高時,雖然此時每個網格中僅有一島,但卻出現了部分尺寸較小的島,這些小島中的原子被一些較大島吸引而離開,這與奧斯特瓦爾德熟化現象相符.隨著進一步遷移,這些小島可能會逐漸消失,從而使得襯底上的2D島的分布不再有規律,島的尺寸也會相差較大,這對于工藝上定位生長量子點或者控制制備量子點尺寸是不利的.

圖2 在不同溫度下2D島的形態模擬結果Fig.2 2D island morphologies results under diあerent temperatures

通過模擬結果可以看出,溫度對于襯底表面2D島的生長影響顯著.在生長制備過程中所選取的溫度不宜過低,因為溫度過低會導致原子的躍遷概率減小,原子在襯底表面的擴散不活躍,就不利于原子聚集;而溫度過高會使原子擴散行為非常活躍,一些較小的島上的原子容易被更大的島所吸引,從而使得島的形態大小分布不再均勻,同樣也不利于量子點的生長,這一結論與實驗結果基本吻合[15].模擬結果顯示,溫度在600 K左右時,2D島的分布更趨于均勻,這表明600 K是最有利于量子點生長的溫度.

1.3 沉積速率對量子點生長的影響

從溫度條件的模擬結果中可以看出,擴散是否充分對生長非常重要,尤其是對于低溫條件,因為此時原子運動不夠活躍,很難擴散充分.那么針對低溫條件,如果給予更充足的擴散時間,即更低的沉積速率是否能得到更好的結果？為此本工作將進一步研究沉積速率對生長的影響.

為了研究GaAs應變弛豫圖形襯底上生長InAs量子點過程中沉積速率對其影響,本工作模擬了溫度在400 K和600 K的量子點前期單原子層生長情形.模擬中選取的模擬參數:能量勢Estr1=0.35 eV,Estr2=0.55 eV,Es=0.3 eV.表面覆蓋率為15%,即模擬中沉積原子數為1 500個,在保持同樣覆蓋率條件下根據沉積速率的不同考慮2種情況:①生長時間t=1.5 s,沉積速率F=0.1 mL/s;②生長時間t=1.5 s,沉積速率F=0.01 mL/s.

由于僅考慮沉積速率的影響,故本模擬中除了沉積速率外的其他參量保持不變,使系統在沉積結束時能達到相同的表面覆蓋率,即要求在不同的沉積速率下都沉積了相同數量的原子.這樣,低沉積速率就需要更長的時間來完成沉積過程,而在沉積過程中擴散過程同時進行.因此,在低沉積速率的情況下沉積原子的擴散將進行得更加充分.考慮到當溫度在400 K、沉積速率為0.1 mL/s時原子分布還比較混亂,只有少量原子聚集成核形成島.圖3為當溫度保持400 K不變時,只改變沉積速率的2維生長模擬結果.結果表明,在沉積剛剛結束時低沉積速率情況下成島的尺寸要比高沉積速率下成島的尺寸更大、更均勻,且島的分布也比較有規律,同一區域內幾乎只形成一個島.因此,當溫度保持在400 K時,如果不考慮由沉積速率太慢引起的沉積過程時間過長可能造成的外在因素的影響(沉積時間過長可能使生長過程中有更多雜質成分進入,從工藝生產制備角度出發,太長的生長時間會導致效率問題),較低的沉積速率有利于2D島聚集成核.相應地,實驗結果也表明較高的沉積速率會降低量子點的質量,較低的沉積速率可以得到分布更均勻、純度更高的量子點[12],本模擬所得到的結果與實驗結果基本一致.

圖3 不同沉積速率下2D島的形態模擬結果(T=400 K)Fig.3 2D island morphologies results under diあerent deposition rates(T=400 K)

相較于圖3結果,圖4為在相同條件下只將生長中襯底的溫度從400 K改成600 K時的情況.當溫度為600 K時,無論是在低沉積速率還是高沉積速率生長結束時都出現了分布均勻的島.從直觀上看二者差別并不大,島的分布情形基本一致,基本上每個島都占據一個區域,二者的差別只是島的大小略有出入.造成以上結果的原因在于,低沉積速率量子點生長比高沉積速率情形下能得到更多的沉積時間,使原子在表面的擴散能進行得更徹底,也就越容易聚集成核形成2D島.在400 K時,高沉積速率下沉積結束時原子在表面的擴散行為并未徹底完成,而低沉積速率給了原子足夠時間去擴散遷移,從而更有利于其生長;當溫度在600 K時,原子的表面擴散進行得比較充分,因此無論沉積速率高還是低均能得到分布均勻的生長情形.可見,在溫度較低或其他原因導致的不利于量子點生長的情形下,可以通過適當地降低沉積速率的方法得到尺寸更均勻的量子點結構.

圖4 不同沉積速率下2 D的島的形態模擬結果(T=600 K)Fig.4 2D island morphologies results under diあerent deposition rates T=600 K

2 埋層深度對量子點生長的影響

2.1 襯底表面的應變能計算

在GaAs襯底上生長InAs量子點的過程中,為了得到應變弛豫的圖形襯底,通常會采用埋藏量子點的方法,即在GaAs襯底中預先埋藏一層InAs量子點,而后在此襯底上進行自組織外延生長.以上的模擬在考慮影響量子點生長因素時,應變弛豫襯底上的應變能都取定值,也就是說襯底是一樣的.然而,事實上在不同的埋層深度下得到的襯底表面的應變能是不同的,大量研究也已證實應變能對于量子點的生長有著重要影響作用[10,13,17],因此需要考慮埋藏量子點的深度對應變弛豫襯底上應變能的影響,進而研究對量子點生長產生的影響.為了考慮埋層深度的影響,需要計算出在不同深度下襯底表面的應變能的分布情況及大小.如果將埋藏在襯底中的量子點看成夾雜物,即將襯底看作夾雜問題中的半無限大基底,則此問題便可以作為夾雜問題來處理.而這類各向異性非均勻夾雜問題可用格林函數法得到精確封閉形式解[18-20].

在將InAs量子點埋藏在GaAs襯底時會產生一個特征應變,該特征應變可由晶格失配度決定,InAs與GaAs的晶格失配度約為7%,故在InAs/GaAs系統中特征應變取為0.07[21-22].對于各向異性非均勻夾雜問題,在夾雜物內部存在特征應變,內外部材料參數分別為和Cijlm.如果等于Cijlm,則非均勻夾雜問題就退化成為均勻夾雜問題.

首先,用D表示基體區域,在基體區域V中有一個特征應變.假設是在無限遠處施加的應力,如果材料是均質的,即沒有非均勻性,則在基體區域中應變場所對應的本構關系為=.因此,一個非均勻問題,即在區域V內外部的本構關系可以根據經典的Eshelby非均勻夾雜關系得到[23-24]

式中,σij和γij為產生的應力和應變場.

引入等效特征應變概念,式(5a),(5b)能夠被等價成

令

則式(6)可以被寫成

對于如上所述的特征應變問題,區域V中的應力平衡方程可表示為[25-26]

將式(8a)代入式(9)可以得到

可以清楚看到,等式(10)的右邊部分代表等效體積力,可寫為

式(11)即為新的特征應變與等效體積力之間的關系式.

通過格林函數法可以得到夾雜區域V中任意點x(x,z)處的特征應變,在點X(X,Z)處產生的位移式中,被積函數是等效體積力與線源格林函數的乘積,x;X)為線源格林函數[19],Cijlm為GaAs材料的彈性常數.

進一步,可以得到相應的應變、應力以及應變能:

通過計算得到應變能的分布形式與之前模擬中提出的應變能分布形式近似相同,是完全對稱的近似于正四棱錐形狀,只是在不同埋層深度下的應變能大小不同.表1給出了在不同埋層深度下應變能的最大值Estr2和最小值Estr1.從表中可以清晰地看出,隨著深度的增大,應變能的大小逐漸減小,而應變能的變化幅度也隨之減小.為了簡化計算量,本模擬中只改變在不同埋層深度下對應的應變能的最大值和最小值,且應變能分布形式均選取正四棱錐形.

表1 不同深度下的應變能計算結果Table 1 Corresponding strain energy results under diあerent depths

2.2 模擬結果分析

在運用KMC方法模擬在GaAs襯底上生長InAs量子點的過程中,選取的模擬參數為模擬區域大小為100×100網格,整個區域又分成10×10的小區域,每個小區域上都是應變能呈正四棱錐分布的漸變的應變弛豫圖形襯底形式.生長溫度T=600 K;生長時間t=1.5 s;沉積速率F=0.01 mL/s;表面覆蓋率為15%,即模擬中沉積原子數為1 500個;埋層深度分別為100,300,500,700,900和1 000 nm.相對應的應變能大小由表1給出.

圖5為埋層深度在100,300,500,700,900和1 000 nm情況下,生長結束時的襯底上2D島的形態.從圖中可以明顯看到,在埋層深度較深時,2D島的分布比較混亂,島的尺寸也較小,甚至在達到1 000 nm時幾乎看不到成形的島,也就是說此時量子點幾乎不生長.而隨著深度變淺,島形狀開始變得有規律,開始出現較大的島,模擬中發現埋層深度在100 nm或更淺時,有些小島會逐漸被較大的島吸引,小島逐漸消失,從而形成更大的島.從模擬結果中還可以看出,深度在100,300和500 nm時都出現島成核現象,在700 nm深度下有少量原子在聚集,但原子分布依舊比較混亂,而在900 nm或1 000 nm時并無島形成,而2D島的成核過程將直接影響到3維量子點的生長,這說明埋藏量子點的深度越深越不利于量子點的生長與制備.

上述結果可從計算得到的應變能角度來解釋.量子點埋藏越深,其在襯底表面所產生的應變能越小且應變能的變化幅度越小,如埋層深度達到100 nm時,襯底表面的應變能大小是1 000 nm深度下應變能的200～300倍,變化幅度約為25 000倍.由式(2)可知,應變能越大躍遷概率越大,沉積原子越活躍也越容易在襯底表面遷移.由于在生長過程中量子點更偏向于向著應變能較大的位置生長,故應變能的變化幅度越大越有利于量子點的生長.而埋層深度在1 000 nm時襯底表面應變能的變化幾乎可以忽略,因此在此深度下沒有出現2D島成核情況.

圖5 不同深度下2D島的形態Fig.5 2D island morphologies under diあerent depths

圖6 不同大小島的數量分布Fig.6 Island quantities distributions with diあerent island sizes

雖然從模擬結果中看,埋藏位置越淺越好,但也不能過淺.首先,在深度100 nm條件下已經出現部分較小的島被較大的島吸引的情形,從而出現小部分離散的原子,這并不利于生長大小分布均勻的量子點.圖6顯示不同大小島的數量分布.從圖6中可看出,相對于300 nm情形下,埋層深度在100 nm時島的大小分布更不均勻,出現較多的大于理想尺寸的島(島的尺寸用島所包含的原子數計量,理想的分布均勻的島含有15個原子).為了衡量2D島分布的混亂程度,本工作計算了島的尺寸分布的標準差,標準差越小說明島分布越均勻.計算結果發現,在100 nm情形下標準差值為4.533,而300 nm時為3.307,顯然在100 nm時島的分布比300 nm時更混亂.另外,如果埋藏的量子點位置太淺,埋藏在襯底中的量子點有可能會直接影響到沉積原子與襯底之間的結合能.當埋層深度太小時,其在襯底表面的應變能會非常大,在埋層深度為50 nm時的大小約為100 nm時的103倍.從躍遷概率式(2)可知,躍遷概率隨著應變能大小呈指數級增長,此時的躍遷概率將會非常大,原子會非常活躍,從實際生長制備角度考慮,生長速度太快不容易控制應變能.

3 結束語

利用KMC方法對GaAs襯底上生長InAs量子點的過程進行了模擬,并考慮了溫度、沉積速率和埋層深度對量子點生長前期單原子層階段的影響.模擬結果顯示,溫度對量子點早期2D成核過程有顯著影響,通過控制溫度可以形成均勻有序排列的2D島,能得到最有利于早期單原子層生長的溫度為600 K左右.另外,通過控制沉積速率也能影響2D島的聚集,在溫度較低時可通過適當降低沉積速率來得到更加均勻分布的島;而在較高溫度時,可以選取較高的沉積速率來提高制備效率.埋藏量子點的深度將改變襯底表面應變能的大小,從而影響2D島的成核過程,結果表明當埋藏過深時量子點形成的圖形襯底上的應變能變小,使島成核趨勢不明顯,甚至不會成核.模擬結果顯示,在埋層深度介于100 nm和500 nm之間時更有利于量子點的生長,而在900 nm時沉積原子將不再聚集.這些模擬結果為制備均勻有序的量子點以及優化生產工藝提供了依據.

[1]田強,楊錫震.半導體低維體系簡介[J].物理實驗,2001,21(7):3-6.

[2]魏育新,陳蕊麗.InAs/GaAs量子點激光器增益特性[J].河北大學學報(自然科學版),2016,36(3):232-236.

[3]劉珂,馬文全,黃建亮,等.含有AlGaAs插入層的InAs/GaAs三色量子點紅外探測器[J].物理學報,2016,65(10):329-333.

[4]李樹瑋,小池一步,矢野滿明.垂直堆垛InAs量子點材料的分子束外延生長[J].稀有金屬,2004,28(1):207-209.

[5]SHCHUKIN V A,BIMBERG D.Spontaneous ordering of nanostructures on crystal surfaces[J].Reviews of Modern Physics,1999,71:1125-1171.

[6]LI X L,WANG C X,YANG G W.Thermodynamic theory of growth of nanostructures[J].Progress in Materials Science,2014,64:121-199.

[7]SCH¨OLL E,BOSE S.Kinetic Monte Carlo simulation of the nucleation stage of the self-organized growth of quantum dots[J].Solid-State Electron,1998,42:1587-1591.

[8]THANH V L,YAM V.Superlattices of self-assembled Ge/Si(0 0 1)quantum dots[J].Applied Surface Science,2003,212/213(2):296-304.

[9]ISHIKAWA T,NISHIMURA T,KOHMOTO S,et al.Site-controlled InAs single quantum-dot structures on GaAs surfaces patterned by in situ electron-beam lithography[J].Applied Physics Letters,2000,76(2):167-169.

[10]楊紅波,俞重遠,劉玉敏,等.影響半導體量子點生長因素的分析[J].人工晶體學報,2004,33(6):1018-1021.

[11]胡冬枝,趙登濤,蔣偉榮,等.Si(0 0 1)斜切襯底上Ge量子點的固相外延生長[J].半導體學報,2002,23(6):604-608.

[12]WANG Q,HUANG Y Q,REN X M.The inf l uence of growth parameters on the formation on InAs/GaAs by MOCVD[J].The International Society for Optical Engineering,2012(3):277-298.

[13]PAN E,SUN M,CHUNG P W,et al.Three-dimensional kinetic Monte Carlo simulation of prepatterned quantum-dot island growth[J].Applied Physics Letters,2007,91(19):193110.

[14]EVANS J W,THIEL P A,BARTELT M C.Morphological evolution during epitaxial thin f i lm growth:formation of 2D islands and 3D mounds[J].Surface Science Reports,2006,61(1/2):1-128.

[15]宋禹忻.動力學蒙特卡羅法仿真量子點生長過程的統計分析與參數優化[D].北京:北京郵電大學,2007.

[16]楊園靜,涂潔磊,李雷,等.生長溫度對InAs/GaAs量子點太陽電池的影響研究[J].人工晶體學報,2014,43(2):461-464.

[17]NURMINEN L,KURONEN A,KASKI K.Kinetic Monte Carlo simulation of nucleation on patterned substrates[J].Phys Rev B,2001,63(3):03054.

[18]SUN L G,XU K Y,PAN E.Irregular inhomogeneities in an anisotropic piezoelectric plane[J].Journal of Applied Mechanics,2012,79:021014.

[19]CHEN Q D,XU K Y,PAN E.Inclusion of arbitrary polygon with a graded eigenstrain in an anisotropic piezoelectric half plane[J].International Journal of Solids and Structures,2014,51(1):53-62.

[20]PAN E.Eshelby problem of polygonal inclusions in anisotropic piezoelectric full-and halfplanes[J].J Mech Phys Solids,2004,52:567-589.

[21]何菊生,張萌,肖祁陵.半導體外延層晶格失配度的計算[J].南昌大學學報(理科版),2006,30(1):63-67.

[22]呂建國,王成彪,于翔,等.薄膜的晶格失配應力分析[J].科技創新導報,2008,20:14-15.

[23]ESHELBY J D.Elastic inclusions and inhomogeneities[J].Progress Solid Mechanics,1961,2:87-140.

[24]MURA T.Micromechanics of defects in solids[M].2nd ed.Dordrecht:Kluwer Academic Publishers,1987.

[25]ZOU W N,HE Q C,ZHENG Q S.General solution for eshelby’s problem of 2D arbitrarily shaped piezoelectric inclusions[J].International Journal of Solids and Structures,2011,48:2681-2694.

[26]PAN E.Eshelby problem of polygonal inclusions in anisotropic piezoelectric bimaterials[J].Proceedings of the Royal Society A,2004,460:537-560.