低維孔隙結構對熱電優值影響之研究

國 立 交 通 大 學

機 械 工 程 學 系

碩士論文

低維孔隙結構對熱電優值影響之研究

EFFECTS OF LOW-DIMENSIONAL POROUS STRUCTURES

ON THE THERMOELECTRIC FIGURE-OF-MERIT

研 究 生： 李 式 堯

指導教授： 陳 俊 勳

低維孔隙結構對熱電優值影響之研究

Effects of Low-Dimensional Porous Structures

on the Thermoelectric Figure-of-Merit

研 究 生：李式堯 Student：Shih-Yao Lee

指導教授：陳俊勳 Advisors：Chiun-Hsun Chen

曲新生 Hsin-Sen Chu

國 立 交 通 大 學 機械工程學系

碩士論文

A Thesis

Submitted to Department of Mechanical Engineering College of Engineering

National Chiao Tung University in Partial Fulfillment of the Requirements

for the Degree of Master

in

Mechanical Engineering

July 2008

獻給

雙親

低維孔隙結構對熱電優值影響之研究

研究生：李式堯 指導教授：陳俊勳 曲新生 國立交通大學機械工程學系 摘 要 本論文旨在探討奈米線線寬、孔隙密度及孔隙間距高寬比對多孔 矽奈米線徑向熱電優值之影響。輻射傳輸方程式將被用來同時描述電 子和聲子於介質中的運動行為，並佐以鬆弛時間近似的假設來簡化能 量載子之間複雜的碰撞機制。當載子運動至材料內任一界面時，採用 彈性鏡面式散射直接對其動量造成改變來反映出界面阻力將成為尺寸 效應下傳輸係數變化的關鍵因素。研究結果顯示，於散射相依的區間 中，多孔奈米線的熱電優值將會隨著孔隙間距高寬比增加而出現週期 性震盪遞減的行為，此時載子的紐森數、孔隙密度以及孔隙間距高寬 比三者被用來綜合判斷孔隙界面是否有效地被散射所利用。當奈米線 線寬為 10 奈米，孔隙間距高寬比介於 0.72 及 1.38 之間時，運動於孔 隙間的能量載子彼此破壞其入射時的狀態而呈現出較窄的能量分佈， 這使得熱電優值在此區間不降反升；當孔隙間距高寬比介於 1.38 及 2.13 之間時，孔隙間的載子能夠持續地保持入射時的本徵態來提高其 返回原材料入射邊界的機率，這會降低材料能量傳遞的能力而使熱電 優值再度下滑。於本文研究中所發現到多孔奈米線的徑向熱電優值將

Effects of Low-Dimensional Porous Structures

on the Thermoelectric Figure-of-Merit

Student：Shih-Yao Lee Advisors：Chiun-Hsun Chen Hsin-Sen Chu

Department of Mechanical Engineering National Chiao Tung University

ABSTRACT

The purpose of this study is to investigate the effects of wire thickness 、 pore density and interpore aspect ratio on the radial thermoelectric figure-of-merit in porous nanowire systems, which incorporate the description of electrons and phonons transport behaviors into radiative equations due to their ballistic nature. Constant relaxation time approximation and elastic specular scattering on the boundaries are utilized to simplify the complicated collision mechanisms for carrier-carrier and carrier-boundary interactions. Although the nanoporous matrix is a connected system where perfect confinement cannot take place, we show how the difference in the magnitude of carrier-boundary scattering give the possibility to obtain relative high thermoelectric figure-of-merit. These results suggest that carefully choosing a set of parameters composed of Knudsen number、pore density and interpore aspect ratio will be very crucial to optimize thermoelectric processes in low-dimensional porous structures, and it could be advantageous for

誌 謝

於此階段完成之際，承蒙恩師 曲教授新生引領下進入極具前瞻性 的研究領域而拓展了我的視野，老師膽大心細、積極且熱情的個人特 質，對我在求知精神及待人處事上兼備言教和身教的意義，在此願致 上最誠摯的謝忱與敬意。 感謝陳教授俊勳長期對師門有形無形的照顧；感謝陳教授朝光在風 雨交加的颱風天前來給予學生許多寶貴的建議和討論，以及黃教授吉 川對論文的指正都使我獲益匪淺。 感謝時明、建新、旭昉、孝彥、木勝、鵬樟、建評、清益、家瑜、 孟熹、子淵、智堡、文耀、彥成學長及純怡學姊在這段日子裡的愛護； 感謝同儕榮祥、祥哲、建宏、炳坤、迪翰、耀文的幫助，和我多年的 老友志煌、允程、子鈞、迪傑、建智、富凱、文聖、定罡、穎志、育 奇，謝謝你們的支持。 特別感謝曾經啟迪我思想的許教授淙慶、王教授星豪、李教授正 中、陳教授特良以及我最親愛的家人，因為有您，也才會有現在的我。 至此，要感謝的人有太多了，請原諒我無法一一表達對您的謝意， 您對我的關懷均將永銘於心。 最後，我要再次感謝恩師無比的耐心、愛心和持續對我的鼓勵， 謝謝您給我扭轉生活態度的動力，以及即將帶來的改變。

目 錄

摘要

...

i

誌謝

...

iii

目錄

...

iv

表目錄

...

vii

圖目錄

...

viii

符號說明

...

xi

一

緒論

...

1

1.1

前言

...1

1.2

熱電應用及原理簡介

...2

1.2.1 發電機及致冷器之應用概況

...2

1.2.2 熱電效應及熱電優值侷限

...3

1.2.3 傳輸係數之微觀探究

...5

1.3

文獻回顧

...7

1.4

本文探討主題

...10

二

理論分析

...

20

2.1

統御方程

...20

2.1.1 波茲曼方程式

...20

2.2

熱電優值分析

...21

2.2.1 熱導率分析

...21

2.3.1 聲子傳輸方程式

...29

2.3.2 邊界條件

...30

2.3.3 熱導率計算

...31

2.4

奈米線徑向二維導電系統

...31

2.4.1 電子運動方程式

...31

2.4.2 邊界條件

...32

2.4.3 電導率計算

...33

2.4.4 塞貝克係數計算

...33

三

數值方法

...

39

3.1

有限差分法

...39

3.2

S

N

近似法

...40

3.3

奈米線徑向熱電系統

...41

3.3.1 二維輻射傳輸差分方程式

...41

3.4

收斂條件

...43

3.5

數值方法驗證

...44

3.5.1 格點測試

...44

3.5.2 結果比較

...45

四

結果與討論

...

54

4.1

奈米線結構之熱電優值分析...54

4.1.1 熱導率之尺寸效應

...54

4.1.2 功率因子之尺寸效應

...55

4.2

單孔奈米線結構之熱電優值分析...55

4.2.1 線寬對傳輸係數之影響

...56

4.4

孔隙式奈米線結構之綜合分析...60

五

結論與建議

...

94

5.1

綜合結論

...94

5.2

建議研究與未來展望

...95

表目錄

表 2- 1 不同空間尺度下之傳輸理論分類 ...34

表 2- 2 能量載子特徵 ...35

表 3- 1 SN近似法之不同方向數與其相應權重關係...46

圖目錄

圖 1- 1 晶片熱通量變化趨勢 ...11 圖 1- 2 各式致冷器性能比較 ...12 圖 1- 3 放射性同位素熱電發電機結構 ...13 圖 1- 4 熱電致冷晶片裝置 ...14 圖 1- 5 塞貝克效應示意 ...15 圖 1- 6 珀爾帖效應示意 ...16 圖 1- 7 湯姆森效應示意 ...17 圖 1- 8 熱電優值發展走勢 ...18 圖 1- 9 熱電文獻發表篇數統計 ...19 圖 2- 1 奈米線徑向系統示意 ...36 圖 2- 2 多孔奈米線徑向系統示意 ...37 圖 2- 3 鏡面反射示意 ...38 圖 3- 1 有限差分法示意 ...48 圖 3- 2 數值計算流程 ...49 圖 3- 3 熱導率格點測試 ...50 圖 3- 4 功率因子格點測試 ...51 圖 3- 5 熱導率驗證 ...52 圖 3- 6 電導率及塞貝克係數驗證 ...53 圖 4- 1 奈米線徑向二維系統示意 ...61 圖 4- 2 二維系統之無因次化溫度分佈(L=10nm)...62 圖 4- 3 功率因子及熱導率隨奈米線線寬變化 ...63 圖 4- 4 無因次化熱電優值隨奈米線線寬變化 ...64 圖 4- 5 單孔奈米線徑向二維系統示意 ...65

...67 圖 4- 8 單孔系統之無因次化溫度分佈(L=1μm, LP=50nm ) ...68 圖 4- 9 單孔系統之功率因子及熱導率隨孔隙尺寸變化 ...69 圖 4- 10 單孔系統之無因次化熱電優值及孔隙密度隨孔隙尺寸變化 ...70 圖 4- 11 單孔系統之無因次化熱電優值隨孔隙密度變化 ...71 圖 4- 12 雙孔奈米線徑向二維系統示意 ...72 圖 4- 13 雙孔系統之功率因子及熱導率隨孔隙尺寸變化 ...73 圖 4- 14 雙孔系統之無因次化熱電優值及孔隙密度隨孔隙尺寸變化 ...74 圖 4- 15 雙孔系統之無因次化熱電優值及孔隙間距高寬比隨孔隙密 度變化...75 圖 4- 16 雙孔系統之無因次化溫度分佈(L=10nm , LP=1.5nm ) ...76 圖 4- 17 雙孔系統之無因次化溫度分佈(L=10nm , LP=2.1nm ) ...77 圖 4- 18 雙孔系統之無因次化溫度分佈(L=10nm , LP=2.7nm ) ...78 圖 4- 19 雙孔系統之無因次化溫度分佈(L=10nm , LP=2.9nm ) ...79 圖 4- 20 雙孔系統之無因次化溫度分佈(L=10nm , LP=3.2nm ) ...80 圖 4- 21 雙孔系統之無因次化溫度分佈(L=10nm , LP=3.4nm ) ...81 圖 4- 22 雙孔系統之無因次化溫度分佈(L=10nm , LP=3.6nm ) ...82 圖 4- 23 雙孔系統之無因次化溫度分佈(L=10nm , LP=3.7nm ) ...83 圖 4- 24 四孔奈米線徑向二維系統示意 ...84 圖 4- 25 四孔系統之功率因子及熱導率隨孔隙尺寸變化 ...85 圖 4- 26 四孔系統之無因次化熱電優值及孔隙密度隨孔隙尺寸變化 ...86

圖 4- 29 四孔系統之無因次化溫度分佈(L=10nm , LP=2.9nm ) ...89 圖 4- 30 四孔系統之無因次化溫度分佈(L=10nm , LP=3.4nm ) ...90 圖 4- 31 四孔系統之無因次化溫度分佈(L=10nm , LP=3.7nm ) ...91 圖 4- 32 多孔系統之無因次化熱電優值隨孔隙密度變化交叉比對 ..92 圖 4- 33 多孔系統之無因次化熱電優值隨孔隙間距高寬比變化交叉 比對...93

符號說明

字母

a 加速度,

[ ]

m/s A 截面積,

[ ]

m Ar 孔隙間距高寬比, L /LP W

[ ]

-c 比熱,

[

J/kg-K

]

V c 體積比熱, _{⎡J/m -kg-K}3 _⎤ ⎣ ⎦ C 熱容量,

[ ]

J/K D 狀態密度, ⎡_⎣1/cm3⎤_⎦ e 基本電荷, 1.602x10-19

[ ]

C E 電場,

[

V/m

]

f 分佈函數, _⎡1/cm3_⎤ ⎣ ⎦ h 普朗克常數, 6.626x10-34

[ ]

J-s = 狄拉克常數, h/2π

[ ]

J-s I 輻射強度, ⎡_⎣W/m -sr2 ⎤_⎦ J 電通量, ⎡_⎣A/m2⎤_⎦ k 熱導率,

[

W/m-K

]

B k 波茲曼常數, 1.381x10-23

[ ]

J/K n K 紐森數, Λ/L

[ ]

-l 材料長度,

[ ]

cm

* m 有效質量,

[ ]

kg n 數量密度, ⎡_⎣1/cm3⎤_⎦ P 功率因子, ⎡_⎣W/m-K2⎤_⎦ Pd 孔隙密度, L /L2P 2

[ ]

-q 熱通量, _⎡W/m2_⎤ ⎣ ⎦ Q 熱流率,

[ ]

W r 空間座標,

[ ]

m R 電阻,

[ ]

Ω S 塞貝克係數,

[

V/K

]

t 時間,

[ ]

s T 溫度,

[ ]

K v 速度,

[ ]

m/s V 電位,

[ ]

V Z 熱電優值,

[ ]

1/K ZT 無因次化熱電優值,

[ ]

-β _{湯姆森係數,}

_[

_V/K

_]

Π 珀爾帖係數,

[ ]

V ε 單位能量,

[ ]

J F ε 費米能量,

[ ]

eV Carnot η 卡諾效率,

[ ]

-η _{方向餘弦, sin cosθ φ}

_{[ ]}

rad

D

θ _{德拜溫度,}

_{[ ]}

K

ν 頻率,

[ ]

Hz ρ _電阻率,

_[

Ω-cm

_]

σ 電導率,

[

S/cm

]

τ 鬆弛時間,

[ ]

s Λ 平均自由徑,

[ ]

nm ω 角頻率, 2πν

[ ]

Hz Ω 立體角,

[ ]

sr

標號

e 電子 ph 聲子 C 低溫 H 高溫 collision 碰撞 o 平衡狀態 < > _平均

一 緒 論

1.1 前言

2007 年，前美國副總統 高爾(Al Gore)及聯合國跨政府氣候變遷 小組(IPCC)獲頒諾貝爾和平獎，諾貝爾委員會於頌詞中表揚其長期致 力於建立及散播人為氣候變遷的知識，並奠定採取必要對策的基礎， 可見得氣候變遷問題已受到極大的關注。目前所知，氣候異常的原因 以六大關鍵性溫室氣體中的二氧化碳影響最甚，其內碳氧雙鍵的共振 頻率正好落於紅外線區，能有效地吸收太陽輻射熱而促進地球暖化。 輻射熱能若是被物質吸收後，根據熱力學第二定律，其運動是屬 於較不具有特定方向性的低品質能量，因此，如何將熱能轉化成高品 質的電能來提供生活所需，除了直接影響各國的經濟發展外，更能減 緩暖化對生態上的衝擊。 熱電裝置亦即能夠直接進行熱電轉換，其內無移動件、環境友善、 反應迅速、無方向性限制及安靜等諸多相對優勢，在發電的應用端， 利用核能熱能直接照射熱電模組產生電力的方式，早已實行於 1960 年代登陸月球的阿波羅太空船中，之後如伽利略號等太空任務皆可見 其蹤跡；在致冷的應用上，需要快速或精準的溫控領域都有發揮的空 間，如醫療器材、實驗設備及計算機工業等。根據英代爾45 奈米製程 技術下產出的微處理器，其內單位晶片上的電晶體數量及生熱密度分 別可達近109顆及每平方公分40 瓦以上，如圖 1- 1所示，因此迅速地 分散或是移除局部高熱將有助於降低熱阻對電訊傳遞的耗損，甚至是 由熱應力所造成的材料破壞。由此可知，於致冷、發熱及環境保護層

1.2 熱電應用及原理簡介

熱電裝置是根據熱電原理發展而成的能量轉換工具，然而其能量 轉換效率僅約為傳統壓縮機的 1/3 [1]，因此造成了廣泛應用上的困難， 如圖 1- 2所示。在本節中，將首先提及熱電裝置於發電和致冷上的應 用後，從基本原理的角度切入瞭解熱電性能難以提昇的原因，再進一 步探討能量載子對於構成材料熱電性能的物理性質有何影響，期望能 由宏觀及微觀的趨近找出提升熱電性能的著眼點，以及研究方向上的 建議。

1.2.1 發電機及致冷器之應用概況

熱電於應用上可以區分為以環境提供的高低溫差來產生電流的發 電機，或是提供電流給熱電模組來製造溫差的致冷器兩大類。實際應 用分別簡述如下： 熱電發電廣泛被使用的場合包含有溫度量測系統、交通運輸工具 的廢熱回收系統及太空船的電力系統等[2,3]。以應用於太空船的電力 系統為例，其利用放射性同位素(Radioisotope)進行核分裂反應時所放 射的輻射熱來提供熱電轉換所需要的熱能來源，此種產生電力的方式 早已被美國太空總署(NASA)應用於超過 20 艘航宇飛行的太空船中， 如圖 1- 3所示。根據 1989 年發射升空的伽利略號(Galileo)太空船，其 內使用兩組放射性同位素熱電發電機(GPHS-RTG)，並以鈽(238Pu)和矽 鍺(Silicon-Germanium, Si-Ge)分別作為核燃料及熱電材料，單組熱電發

電子裝備冷卻系統等。以電子裝備的冷卻為例，根據英代爾提供的新 型散熱設計中，將熱電致冷器置於散熱片與微處理器之間來加速散熱 以 提 供 更 低 的 晶 片 表 面 溫 度 ， 若 是 在 有 致 冷 器 的 情 況 下(ZT=1 , COP=2)，提供其 60 W來移除晶片上 120 W的熱能，相較於沒有致冷 器之下能夠使得晶片表面溫度下降達 12 之多，而散熱片與環境間 的熱阻( o_C SA θ )將由 0.18 /W降低至 0.08 /W[4]；抑或將熱電材料埋入 晶片用以降低熱點(Hot spot)對微處理器運算效能及材料上的影響，如 o_C o_C 圖 1- 4所示，都是近年來被廣泛討論的話題[4-8]。

1.2.2 熱電效應及熱電優值侷限

1822 年，德國科學家塞貝克(Thomas Johann Seebeck)為熱電研究 揭開序幕後，同法國錶匠珀爾帖(Jean Charles Athanase Peltier)及英國 物理學家湯姆森(William Thomson, Lord Kelvin)的發現奠定了熱電效 應 的 基 礎 。 就 目 前 所 知 ， 熱 電 效 應 可 分 為 ： 塞 貝 克 效 應(Seebeck effect)、珀爾帖效應(Peltier effect)和湯姆森效應(Thomson effect)。

塞貝克效應： A和B兩種導體串聯形成封閉迴路，若是材料相接觸的兩端點具有 溫度差，則可於封閉迴路內量測到電位差及其所伴生的電流，如圖1- 5所示，且相異接觸端點上的溫度差和電位差成正比的關係： V S T Δ = − Δ (1- 1) 其 中 為 塞 貝 克 係 數 (Seebeck coefficient) 或 稱 之 為 熱 電 勢 (Thermoelectric power)，起因於塞貝克效應是描述熱發電的物理現 象，同時也是熱電偶(Thermocouple)的基本原理。 S

電流通過材料相異接觸點時，分別會出現吸熱或放熱的反應，如 圖 1- 6所示，且吸收或釋放的熱量正比於電流大小： Q= ΠJ (1- 2) 其中 為珀爾帖係數(Peltier coefficient)，Q為熱流率， 為電流。珀爾 帖效應是描述電流能夠製造高低溫差的物理景象，即熱電致冷器 (Thermoelectric cooler)的基本原理。 Π J 湯姆森效應： 導體內同時存在溫度梯度及電流的情形下，其內不同位置皆會有 有吸熱或放熱的行為稱之。當電子行進的方向與溫度遞減的方向一致 時，則電子於運動過程中將不斷地釋放能量來與環境達成熱平衡，此 時材料會向外界釋放其多餘的熱能；反之，若電子流動與溫度遞減的 方向相反，則其於傳輸過程中將向材料吸收能量來增加自身的位能而 使得材料對外界進行吸熱反應，如圖 1- 7所示。湯姆森並找出三種熱 電係數之間的關連性，如下所示： ST = Π (1- 3) dS T dT β = (1- 4) 式(1- 3)及式(1- 4)合併稱為凱爾文關係(Kelvin relation)，其中β為湯姆 森係數(Thomson coefficient)。 上述之熱電效應皆是以導體為研究對象時所發現到的物理定律， 且皆為可逆並具有相關性(Correlation)，根據熱電裝置的需要，定義出 熱電優值(Thermoelectric figure-of-merit, Z)[2]為其性能指標，ZT 為無

conductivity)、 為塞貝克係數、S T 為材料平均溫度。由式(1- 5)可知， 提高材料的電導率、塞貝克係數或是降低熱導率皆能增加ZT，然而以 導體為材料時，根據維德曼-夫蘭茲定律(Wiedemann-Franz law)[9]，熱 導率與電導率成正比關係： k LT σ = (1- 6) 其中L為勞侖茲常數(Lorenz number)，T為絕對溫度。從上式中發現到 熱導率與電導率將同時上升或下降，而電導率和塞貝克係數卻往往呈 現相反的變化趨勢，此三個參數相互耦合的結果對於 ZT 在高導電、 低導熱及高塞貝克係數的要求上，造成了嚴重的限制。

1.2.3 傳輸係數之微觀探究

影響熱電性能的傳輸係數，由下而上(Bottom-up)的面向將能夠對 其有更具體的認識。以電導率來說，決定導電的載子僅為電子，而且 位於原子核束縛合力愈小的能級區間中愈有機會受到外界刺激而產生 群體移動，甚至有脫離功函數(Work function)的可能。從固態物理發 展出的能帶理論(Band theory)得知，處在費米能階(Fermi level)附近的 價電子對導電最有貢獻，因此藉由摻雜(Doping)來提高費米能階，或 者改變材料溫度等增加導帶電子數量的技術，都能直接提升電導率。 在金屬中，導帶與價帶部分重疊的結果，電子並不被束縛於特定的原 子核作用力範圍，這使得導體與絕緣體在電導率的差異上，足足橫跨 有 25 個數量級(Order of magnitude)之多[10]。 電子受力行進期間，能量及動量的變化與散射的機制有關，散射 若是以輻射的說法，包含有反射(Reflection)、折射(Refraction)及繞射

-晶界及電子-晶體表面邊界等諸多影響因素。當材料溫度遠小於德拜 溫度時(T θD)，電子主要受到晶格缺陷等結構性阻礙而多發生彈性碰 撞，則電導率與溫度的關係並不顯著；當材料溫度接近於德拜溫度時 (T ~θ_D)，聲子與電子的交互作用嚴重影響了電子的平均自由徑(Mean free path,Λ )， 此 時 電 導 率 與 溫 度 成 五 次 方 反 比 的 急 遽 變 化 關 係 ( 5)；當材料溫度遠大於德拜溫度( T σ _∝ − D T θ )，電子的散射幾乎是以 聲子為主，電導率將與溫度成一次方反比的關係( )[9,10]。因 此，室溫下材料的導電能力深受聲子的密度影響，而聲子(Phonon)為 晶 格 振 動 能 的 量 子 單 位 ， 其 允 許 的 量 子 態 能 量 為 1 T σ _∝ − 1 ( ), 0,1, 2,3,... 2 n h n n ε = ν + = ，不同的量子態將對應到相異的振動頻率，反 過來說，材料內的高溫區會具有高能且高密度的聲子來進行能量傳 送，由於聲子在傳遞的本質上即是晶格間相互影響的熱波行為，晶格 內原子間位能的非諧和本性(Anharmonic nature)導致聲子間的作用多 屬於非彈性散射，亦即碰撞前後要發生能量上的變化，此種非彈性碰 撞可分為正常過程(Normal process)及變向過程(Umklapp process)。

於正常過程中，二聲子碰撞後合而為一或是一分為二皆不對動量 造成影響，在這種狀況下，雖然聲子能量獲得了重新分配，但是在材 料內傳遞時並無對其動量造成阻礙，這將會得到熱導率趨近於無限大 的結果；而變向過程的模式認為聲子碰撞前後在能量及動量上都要進 行交換，此種解釋與實驗瞭解到有限的熱導率是一致的，然而影響熱 導率的重要參數尚有材料比熱(Specific heat, )。當材料溫度小於德拜c

趨近於杜龍-柏帝定律(Dulong-Petit law)預測的常數( )，此時聲子 平均自由徑的降低將會直接造成導熱率下降。除了聲子之間的非彈性 碰撞過程，如同電子一般，聲子與材料內缺陷及邊界的作用多為不造 成能量上衰減的彈性碰撞[9,10]。 3 c= R 前述聲子在熱導率上的貢獻，同樣發生於電子。在導體中，自由 電子(Free electron)同時具備導電以及導熱的能力，但是在非導體中， 雖然並無自由電子來進行傳熱，但聲子卻實質地幫助了熱能的傳遞， 這使得熱導率在材料上的變化範圍僅達到五個數量級左右。 另外，當導體兩端具有溫度差時，熱端能量只需要達到一定數量 級( )，即足以提供電子由價帶躍遷至導帶所需要的能量，此時 熱端的電子會自發性地往低能量的冷端流動而造成材料內電荷分佈不 均的現象，這種類似霍爾效應(Hall effect)的結果將導致材料內電勢的 產生來阻擋電子持續地由熱端向冷端擴散，一旦電子所受之庫倫力與 熱擴散力達成平衡時，材料內電勢不再增加，塞貝克係數正好用來衡 量二者的平衡位置。 ~ (O k TB )

1.3 文獻回顧

自從19 世紀初期發現熱電效應以來，熱電在應用上的發展相當緩 慢，起因於金屬中的自由電子同時支配傳熱及傳電而無法改善ZT參數 耦合的問題，當代的困境於百年後獲得重要舒緩。1950 年代半導體工 業發展之初，利用半導體本身在導熱及導電載子上的差異，配合摻雜 異質的合金技術來增強導電降低導熱，以當時最好的熱電材料碲化鉍 合金(Bi2Te3 alloy)為例，就算其室溫下的ZT成功地被提升至 1 左右 (~10%

η

)，仍舊無法與傳統壓縮機的性能相比[12]，如圖 1- 8所示。

如圖 1- 9所示，在此期間，零星的研究仍然著重於如何獲得更低的熱 導率上，然而所採取的策略卻已經大不相同。例如Yao[13]以分子束磊 晶(Molecular beam epitaxy)的技術，沈積出砷化鋁/砷化鎵(AlAs/GaAs) 薄膜交互堆疊的超晶格(Superlattices)結構來研究材料特徵尺寸方向的 熱傳情形。其研究指出，較高的堆疊週期數能夠降低材料的傳熱能力， 並且隨著堆疊週期數不斷地上升，熱導率將會接近於鋁鎵砷三元合金 (Al0.5Ga0.5As alloy)所能提供的低導熱係數，事實上，在此之前以人工 的(Man-made)方式製備而成的超晶格結構已經開始受到了學界的注 意。1982 年Ren和Dow即對超晶格結構進行研究，其基於卡拉威模式 (Callway’s model)下，考慮三聲子的變向散射(Umklapp scattering)後認 為相異的薄膜堆疊週期數的確造成熱傳係數下降，至此以後，利用調 整結構的方式來改變材料熱物理性質上的研究從 1990 年代起逐漸蓬 勃發展。 1993 年，Hicks 和 Dresselhaus[14,15]探討當時最好的熱電材料碲 化鉍合金，若是將其製作成量子井結構，則當薄膜厚度由1 下降至 時，ZT 可從 2.5 上升至 5；若是製作成量子線結構，則當線寬為 時，ZT 更可高達 14。此種低維度結構之所以能夠突破 1960 年 代無法有效提升熱電優值的瓶頸，在於低維結構能改變塊材(Bulk)的 能態密度(Density of state)來提高功率因子(Power factor)同時降低熱導 率。在降低導熱的研究上，Majumdar[16]探討鑽石薄膜特徵尺寸方向 的傳熱時認為固態材料在微尺度中的熱傳行為並無法利用宏觀下發現 的傅立葉定律來描述，其依循能夠描述局部非平衡態區域內載子運動 nm 0.5nm 0.5nm

之多，此種低維材料的物理性質隨著材料尺寸變化的現象被稱之為古 典尺寸效應(Classical size effect)。基於聲子輻射的概念下，聲子在界 面處的散射將會嚴重地影響低維材料的傳熱能力，而界面熱阻的理論 模式早於 1950 年代末期被提出用來分析當熱能流經兩塊材之間時， Little[17]考慮彈性鏡面反射的模式來描述聲子穿越界面所遭遇到的阻 力。Swartz 和 Pohl[18]於 1989 年再根據同一個主題進行探討後認為聲 子行經界面時，應該發生不具有能量損耗的彈性漫射較為合理。以上 兩種界面熱阻的基礎模式成為了 1990 年代以後一系列低維材料的熱 傳研究中所不可忽略的重要關鍵[19-27]，因此，低維傳熱的主題開始 由二維結構往更低維度前進[28-30]，並探討於材料內添加異質結構來 增加散射面積的效應，例如將矽奈米線週期排列於鍺基材中，研究矽 鍺奈米線徑向傳熱時發現當固定矽奈米線線寬，隨著鍺的莫爾分率下 降並無法獲得導熱增強的效果，表示宏觀中由材料本身的物性來決定 其傳輸係數的大小在介觀(Mesoscope)下已經完全不同，這是由於載子 在介觀中的運動傾向於彈道式(Ballistic)傳輸，因此界面對聲子運動造 成的阻力將會直接影響材料熱物理性質的變化。導熱率的尺寸效應同 樣受到實驗上的證實，2004 年 Song 和 Chen[31]量測多孔單晶矽薄膜 時發現到當薄膜長寬各為數毫米及數百微米，微米級的孔隙間距已經 使得導熱率出現強烈的尺寸效應，其於隨後的研究更進一步提出若是 以添加孔隙的方式來提升材料的 ZT，則孔隙尺寸座落於聲子及電子二 者平均自由徑之間較佳[32]。 在功率因子的研究方面，Hicks 和 Dresselhaus[14,15]同時對低維 度線材及面材進行研究後指出 ZT 將隨著材料尺寸縮小而急遽的上

雜度。Sofo 和 Mahan[33]認為超晶格結構中兩層薄膜間應該考慮電子 的穿遂效應，其應用克若尼-潘尼模型(Kroniq-Penny model)來描述固體 中電子的能量隨著動量變化的色散關係(Dispersion relation)，此作用是 來自於晶體內原子排列的週期性將使得電子運動受到以布里淵區 (Brillouin zone)為週期的有限位能屏障影響，再配合鬆弛時間假設後求 解波茲曼方程發現到 ZT 隨著薄膜厚度及薄膜間距變化會有極值出 現。Boido 等[34-38]同樣從量子井、量子線延伸至複合材料的研究， 由於複合結構具有穩定其內線材以及降低製作成本等優勢，可能有助 於增強熱電材料於應用端的發展。 近年在實驗上的研究，如 2001 年 Venkatasubramanian 等[39]製作 碲化鉍/碲化銻(Bi2Te3/Sb2Te3)超晶格，量測其室溫下 ZT 可達到 2.4；

Harman 等[40]將鉍摻雜入 N 型鉛硒碲/碲化鉛(N-type PbSeTe/PbTe)形

成量子點超晶格結構，其最佳 ZT 值更可高達 3.5。

1.4 本文探討主題

綜合過去研究得知，提高熱電優值的概念有 選擇傳熱和傳電機 制差異大的材料 材料物性要具有容易被調整的潛力(3 以摻雜的技 術提高電荷載子的濃度增強導電 利用低維度結構產生較差的能帶 重疊性來獲得不連續的狀態密度或者較不對稱的載子能量分佈，皆能 提供較高的塞貝克係數(5)增加材料界面的方式來破壞晶格間的交互 作用降低導熱。因此，本文選擇極具發展潛力的本徵半導體矽為材料 [41,42]，以增加界面的方式探討低維奈米線結構內含孔隙的情況下， (1) (2) ) (4)

Year of Announcement 1996 1998 2000 2002 2004 2006 Ch ip He at F lu x ( W /c m 2 ) 0 10 20 30 40 50 60 McKinley i & pSeries Freeway Sun UltraSPARC Pentium 3 Pentium 4 _zSeries Intel Power 4 (1.3 GHz)

Air Cooling Limit for proceesors

Temperature Ratio, T_hot /T_cold Coef ficient of P e rformance, COP 0 1 2 3 4 5 6 1.0 1.1 1.2 1.3 1.4 1.5 1.6 1.7 1.81.92.0 Carnot Cycle Stirling Refrigerators Household Refrigerators and Air-Conditioners Thermoelectric Refrigerators Stirling Cryocoolers 2 1 0.5 7 4 10 m ZT 圖1- 2 各式致冷器性能比較[19]

Heat Source Support

Cooling Tubes Aluminum Outer Shell Assembly

General Purpose Heat Source (GPHS)

Active Cooling System (ACS) Manifold Multi-Foil Insulation Silicon-Germanium (Si-Ge) Unicouple Pressure Relief Device Midspan Heat Source Support Gas Management Assembly RTG Mounting Flange 圖1- 3 放射性同位素熱電發電機結構[2]

l L

L

Ac

b

Insulating Layer/Silicon Substrate

Membrane Hot Junction Cold Junction Cold Area Metalization P-Type Peltier Thermoelement N-Type Peltier Thermoelement Bonding Pads 圖1- 4 熱電致冷晶片裝置[5]

A

B

Ι

Heat

input

Electrical

power

output

圖 1- 5 塞貝克效應示意[2]

A

B

Ι

Heat

absorbed

Electrical

power

input

圖 1- 6 珀爾帖效應示意[2]

A

B

T

T+ΔT

Ι

Heat

released

Electrical

power

input

Heat

absorbed

圖 1- 7 湯姆森效應示意[2]

Year 1950 1960 1970 1980 1990 2000 2010 Dimens ionles s T h ermoelect ric Figure-of -M e rit , ZT 0 1 2 3 4 HIGH TEMPERATURE ROOM TEMPERATURE MnTe ZnSb SiGe PbTe TAGS Skutterudites (CIT-JPL) LAST (MSU) PbSeTe/PbTe QDSL (MIT-LL)

BiTe Spin-cast BiTe (Japan) BiTe/SbTe SL (RTI) US DAPRA/ONR Program Start US DoD Goals (2004) 圖1- 8 熱電優值發展走勢[2]

Year 1950 1960 1970 1980 1990 2000 Pe rc e n t ( % ) 0.00 0.01 0.02 0.03 0.04 0.05

Source: Web of Science Keyword: ‘thermoelectric’

二

理論分析

當材料特徵尺寸下降至微奈米級時，能量載子在局部區域內不一 定滿足平衡態分佈，此時傳輸方程的對象將由局部平衡態性質，轉變 成對載子分佈的探討。以下首先介紹本研究依循之統御方程，進而說 明傳輸係數的計算模式。

2.1 統御方程

於不同空間尺度下描述相同物理量之變化時，Chen[43]以材料特 徵尺寸(L)、粒子平均自由徑(Λ)及其相干波長(λC)來界定出相應的傳 輸方程，如表 2- 1所示。當LΛ，粒子在運動期間不斷地發生碰撞， 而使其傳遞類似以擴散的方式進行，基於此宏觀角度下，如牛頓剪應 力定律、傅立葉定律及歐姆定律等皆可因應而生；當Λ > >L λ_C，粒子 間的碰撞距離可能大於其與界面的距離，此時界面效應足以改變材料 的物性；當LλC，粒子的相位訊息將影響其動力學性質。以下便針 對廣泛使用於介觀中的粒子性波茲曼方程作介紹。

2.1.1 波茲曼方程式

根據呂維爾定理(Liouville’s theorem)，粒子狀態若於傳輸時不發 生任何變化，則其相空間分佈函數(Phase space distribution function, f )

守恆於同一條流線上，可表示為：

及加速度。 若考慮粒子於運動過程中發生碰撞，則其分佈函數將於傳遞時發 生變化，表示如下： ( ) r v c o l l i s i o n f v f a f t t f δ δ ∂ _{+ ⋅ ∇ + ⋅ ∇ =} ∂ (2- 2)

上 式 即 為 波 茲 曼 傳 輸 方 程 式(Boltzmann transport equation) 。 其 中

( )c o l l i s i o n f t δ δ 為粒子發生碰撞時分佈函數的時變項，此粒子碰撞變化 項極為複雜，而使得欲直接求得波茲曼方程之解析解幾乎是不可能 的。在氣體分子動力學領域中，Bhatnagar、Gross 和 Krook 於 1954 年 提出的近似法(BGK approximation)是以鬆弛時間(Relaxation time)的概 念來簡化粒子於真實情況下碰撞的複雜度，除此之外，Fokker 和 Planck 提出的近似法(Fokker-Planck approximation)也普遍使用於該領域中 [44]。 由於波茲曼方程式並不受限於局部平衡態的假設，因此可被用來 描述局部非平衡區域內的傳輸現象。本研究即是依循粒子性波茲曼方 程，探討局部非平衡系統中能量載子的傳輸性對熱電優值的影響。

2.2 熱電優值分析

熱電優值是綜合傳熱、傳電及熱電三種物理現象所得到用來衡量 材料熱電轉換能力的優劣指標，構成熱電優值的性質分別為熱導率、 電導率和塞貝克係數，其中電導率與塞貝克係數又可合併稱為功率因 子，各傳輸係數的內涵將分述如下。

2.2.1 熱導率分析

熱導率。根據熱能在全頻電磁波之波長區間的定義為0.1 ~ 100 mμ [11]， 能夠在此波長區間中幫助傳熱發生的熱載子，可以是電子、聲子、光 子及分子四種[43]，各別特徵如表 2- 2所示。在固態和室溫下，並無 分子的對流熱傳且隱含光子輻射於電子與聲子之中後，一旦高溫施加 於材料一端時，Kaganov等[45]認為比熱較小的電子會將熱能先行吸收 後再傳遞給聲子，被稱為兩階段模式(Two-step model)，此時熱能將從 高溫端往材料低溫端傳送，這種逐漸影響的熱波行為，在利用短脈衝 雷射加熱(Short-pulse laser heating)金屬的研究中佔有極其重要的影響 力，隨著系統達到穩態時，傅立葉定律(Fourier’s law)則可用來描述宏 觀下的熱傳行為。 因此，傅立葉定律是在穩態及滿足局部熱平衡條件下才成立的， 從宏觀的角度來說，熱導率是由材料的本性來決定，以下將基於古典 動力學理論，在局部熱平衡的假設下對熱導率的本質進行初步瞭解。 根據熱通量的定義為單位時間下通過單位面積的熱流量，若以朝 向x方向的熱通量qx為例，每單位體積中的熱載子數目n具有x方向的 平均速度<vx >和能量ε，則x方向的熱通量可表示為： x x q = <n v >ε (2- 3) 其中熱載子的能量ε可以用材料的熱容量(Heat capacity, C)置換後得 到： x x q = <n v > ΔC T (2- 4) 由於溫度為x的函數，且利用鬆弛時間τ 近似來表示溫度變化之距離， 上式改寫為：

能量的1 3，即 = 2 <vx> 2 1 3< >v 。式(2- 5)整理後可得： 2 1 3 x dT q nC v dx τ = < > (2- 6)

再以體積比熱(Volumetric specific heat, )及平均自由徑(Mean free

path, )取代後得到： V c Λ 1 3 = < > Λ x V dT q c v dx (2- 7) 上式為熱通量的展開形式，雖然前述的推導過程仍然是在古典動 力學理論之下，卻已經足以發現一些有趣的性質來作為宏觀與微觀的 連接橋樑。其與傅立葉定律比較後可將熱導率表示： 1 3 = _V < > Λ k c v (2- 8) 其中cV 為材料體積比熱，< >v 為熱載子的平均群速， 為熱載子的平 均自由徑。在純金屬中，傳熱主要是藉由電子來進行的，然而於合金 或是摻雜半導體(Doped semiconductor)中，電子及聲子皆為熱載子，此 時熱導率表示如下： Λ 1 1 3 3 = + = e < > Λ + ph < > Λ electron phonon V e e V ph ph k k k c v c v (2- 9) 由 1.2.3 節分析可知，構成熱導率的三個參數彼此並非完全獨立的， 若從平均自由徑的面向來看，熱載子於傳遞時所受到的散射將會直接 影響其平均碰撞距離，而馬賽厄斯法則(Mattheisen’s rule)正可用來簡 化複雜的散射行為，其使用的是阻力串聯的概念，可表示為： 1 1 , 1.2.... i i i τ =

∑

τ = (2- 10) 其中τ 為粒子由非平衡至平衡態所經過的時間，即鬆弛時間。從上式 可以瞭解到，馬賽厄斯法則的重要假設為各種阻力彼此是獨立的，因

可轉化成平均自由徑的形式，而各種散射機制的影響可以從平均自由 徑與散射截面積(Scattering cross section)的關連，延伸至散射截面積與 粒子尺寸參數(Particle size parameter)的關係來一探究竟。

以上是在滿足局部平衡態的情況下對熱導率所進行的探討，若是 當材料尺寸降低至粒子平均自由徑的數量級時，局部非平衡態的描述 將更具有實質上的意義。從微觀的角度切入，以本文研究之 4A 族本 徵半導體(Intrinsic semiconductor)矽為例，聲子為其主要熱傳載子，則 式(2- 3)中朝向x方向熱通量qx將改寫為： ( ) x x q =

∫

v =ωf Dω ω ωd (2- 11) 其中每單位體積下的熱載子數目 ，已經由處於能量區間n [ ,ω ω+dω]的 載子狀態數目D( )ω 及其單位狀態的機率密度函數 f_ω兩者取而代之。其

中 為狄拉克常數(Dirac constant)，= ω為角頻率(Angular frequency)，而

分佈函數 f_ω必須從波茲曼方程式推衍得到。 由於聲子傳播速度約為定值(_{~ (10 )}2 m O s )，因此可忽略 fω隨著速度 變化的影響(∇vfω)，並考慮穩態的情形，則波茲曼方程式變化成： ( ) r c o l l is i o n f v f t ω ω δ δ ⋅ ∇ = (2- 12) 再將粒子間的碰撞以鬆弛時間近似之後，朝向x方向變化的分佈函數 展開為： o x d f f f v d x ω ω ω τ − = (2- 13) 假設分佈函數 f_ω隨著x的改變接近於平衡態變化，而 o f_ω 又是溫度的函

o o x d f d T f f v d T d x ω ω = ω −τ (2- 15) 再將其帶回x方向熱通量的表示式，則式(2- 11)改寫如下： 2 _{( )} o x x df dT q v D dx dT ω (2- 16) d τ ω ω ω = −

∫

= 上式即為傅立葉定律的微觀形式，熱導率為： 2 _{( )} o x df k v D d dT ω τ ω ω ω =

∫

= (2- 17) 由於聲子為玻色子(Boson)，因此平衡態分佈函數滿足玻色-愛因斯坦 分佈(Bose-Einstein distribution)： 1 exp( ) 1 = − ₋ o B f k T ε μ (2- 18) 其 中ε 為 能 量 ，μ 為 化 學 勢(Chemical potential)， 為 波 茲 曼 常 數 (Boltzmann constant)，T為材料局部溫度。在一般情況下， B k μ可利用費 米能量εF來取代[9,10]。 經過以上在多重尺度下對熱導率的分析後，相同的邏輯將應用於 下節對電導率及塞貝克係數的探討。

2.2.2 功率因子分析

當施加不隨時變的電位差於材料兩端時，材料內部電子會由低電 位往高電位移動，此時歐姆定律(Ohm’s law)正好用來描述宏觀穩態 下，電流與電位差之間的關係，表示為： V RI Δ = (2- 19) 其中R為電阻(Electric resistance)，ΔV 為材料兩端電位差，I為材料內 的電流。若材料自身的電阻較高，則外加電位差所能驅動的電荷流量

上式可見，沿著電子行進方向的材料長度l增加或者縮減材料截面積 A，都會使得材料內電阻上升，而R與長度參數 l A之間的關係稱之為 電阻率(Electric resistivity, ρ)。再帶回式(2- 19)中得到： l V I A ρ Δ = (2- 21) 由電場E的定義為單位長度所能產生的電位差，則可進一步改寫為： I E A ρ = (2- 22) 此時，形同於傅立葉定律的通量形式已呼之欲出。定義電導率σ為電 阻率之倒數，上式將可變化成： J =σE (2- 23) 其中 為電通量(Electric flux)，即每單位時間下通過單位面積的電量，J E為外加電場，σ 為材料的電導率。由此可見，σ 在宏觀下是材料的 物質本性。 根據電通量的定義為單位時間通過單位面積的電荷量，若以朝向 x方向的電通量J_x為例，每單位體積中的電荷載子數目 具有n x方向的 平均速度<v_x >和各別攜帶的電量 ，則e x方向的電通量可表示為： x x J = <n v >e (2- 24) 電荷載子所擁有的方向性運動速度以加速度置換後得到： x x J =na eτ (2- 25) 其中τ 為粒子發生碰撞的平均時間，即鬆弛時間。將電荷載子受力運 動與電場的關係帶入整理後可得： 2_τ (2- 26)

2 * ne m τ σ = (2- 27) 由上可知，較小的有效質量、較高的載子密度及較長的散射距離都能 有效提高電導率。若是從散射距離切入來討論電子的散射機制，大致 可分為晶格缺陷(Defect)、雜質(Impurity)或是與聲子間交互作用的影 響，再利用上節提及的馬賽厄斯法則，則電子傳輸過程時所遭遇到的 阻力表示如下： 1 1 1 o p (2- 28) h τ τ= +τ 其中τ_o為晶體不完美造成的結構性散射項，τ_ph是由晶格震動所提供。 根據能帶理論，當電子運動於完美晶體內則能夠持續地處於某個本徵 態(Eigenstate)而不受到任何影響，一旦材料溫度升高，晶格振動將破 壞週期性位能場，則電子平均自由徑的大小就會開始從晶體結構性缺 陷轉向由與聲子作用的劇烈程度來決定。 以上是基於局部平衡態假設下對電導率的描述，同樣地，由下而 上的探討將從分佈函數開始。考慮朝向x方向運動的電荷載子，則其 電通量Jx可以表示為： ( ) x x J =

∫

v ef D_ε ε εd (2- 29) 式(2- 24)中之每單位體積下的電荷載子數目 ，已經由能量區間n [ ,ε ε+dε]中載子的狀態數目D( )ε 及其單位狀態的機率密度函數 f_ε兩者 取而代之，而分佈函數 f_ε必須從波茲曼方程式推衍得到。另外，電荷 載子可以為帶負電的電子抑或帶正電的電洞，往後將持續以電子來簡 化說明。 穩態時，電子波茲曼方程式表示如下：

分佈函數隨著速度變化的加速度項。將電子的加速度置換成與電場的 關係，並以鬆弛時間近似來描述分佈函數的碰撞時變項，朝向x方向 變化的分佈函數展開為： ∂ ₊ ∂ ₌ − ∂ ∂ o x x f e E f f f v x m v ε ε ε τ ε (2- 31) 假設分佈函數 f_ε隨著x的改變接近於平衡態變化，而 o fε 又是費米能階 的函數，上式左項可近似為： o F x x o x x f e E f f d e E f v v x m v d x m v ε ε ε ε ε ε ε ε ∂ ∂ ∂ ∂ ∂ + = + ∂ ∂ ∂ ∂ ∂ (2- 32) 式(2- 31)經過整理後，則分佈函數 f_ε變成： o o F x f d f f v e E d x ε ε ε ε τ ε ∂ ⎛ ⎞ = + _⎜ + _{⎟ ∂} ⎝ ⎠ (2- 33) 再將上式帶回x方向電通量的表示式，整理後則式(2- 29)改寫如下： 2 2 _{( )} 3 o x f e dV J v D dx ε d τ ε ε ε ∂ = ∂

∫

(2- 34) 上式即為電通量的微觀形式，其中V 為材料端電位。與式(2- 23)相比 較後則電導率為： 2 2 _{( )} 3 o f e v ε D d σ τ ε ε ∂ = ∂

∫

ε (2- 35) 由於電子為費米子(Fermion)，因此平衡態分佈函數滿足費米-狄拉克分 佈(Fermi-Dirac distribution)： 1 exp( ) 1 = ₋ + o B f k T ε μ (2- 36) 其中ε為能量，μ為化學勢， 為波茲曼常數，T為材料局部溫度。在 一般情況下， B k μ可利用費米能量ε 來取代[9,10]。

子的分佈梯度可由費米能階及溫度展開後表示為： o o o F F f f d f dT x dx T dx ε ε ε ε ε ε ε ε ∂ ∂ − ∂ = − − ∂ ∂ ∂ (2- 37) 代回式(2- 32)整理後再帶入式(2- 34)，則來自於溫度差所產生的電流 為： 2_{( )} 3 o F F x f d e dT J v eE D dx T dx ε ε ε ε ( )d τ ε ε ε ∂ − − = + + ∂

∫

(2- 38) 上式可再改寫成如下形式： 11( ) 12( x dV dT J L L dx dx = − + − ) (2- 39) 其中電場已由電位梯度∇V 取代， 為電導率， 為溫度梯度與電流 之間的關係係數，又被稱為耦合係數(Coupling coefficient)。當材料冷 端電荷密度達到飽和時，於此形同開路電路(Open circuit)的情況下定 義出塞貝克係數： 11 L L₁₂ 2 12 2 11 ( ) ( ) / 1 / _{( )} ∂ − − _∂ = = = − ∂ ∂

∫

∫

o F o f v D L dV dx S f dT dx L eT v D d ε ε d τ ε ε ε ε ε τ ε ε ε (2- 40) 上式即為塞貝克係數的微觀形式，其所代表的物理意義為材料端電位 達穩態時，每個熱電子擁有的平均動能大小。

2.3 矽奈米線徑向二維導熱系統

本文研究系統之二維傳熱統御方程式、邊界條件及熱導率的計算 方法，將分述如下。

2.3.1 聲子傳輸方程式

1993 年，Majumdar[16]於研究鑽石薄膜內聲子傳熱行為時，將聲

別於光子傳輸時所受到的散射、吸收及放射效應，而改以近似平衡態 的鬆弛時間處理，並且忽略頻率的影響，即灰體假設。 EPRT 使得原先處理粒子分佈 f_ω的問題，轉變成對能量I_ph的探 討，轉換如下： ( , ) ( ) ph I =

∫

vθ φ ω= f D_ω ω ωd (2- 41) 其中I_ph為聲子輻射強度(Intensity)，=ω為聲子的能量單位，v( , )θ φ 為聲 子的運動速度，而D( )ω 為其狀態密度。 二維穩態時，EPRT 將會變成如下形式： o p h p h p h p h x y I I I I v v x y τ ∂ ∂ − + = ∂ ∂ (2- 42) 其中 o ph I 為平衡態的聲子強度，τ 為鬆弛時間。若將平均自由徑引入上 式並展開平衡態聲子強度後，則上式可改寫為： 4 1 4 Ω − ∂ ∂ + = ∂ ∂ Λ

∫

p h p h p h p h I d I I I x y π π μ η (2- 43) 其中 為曲面所對應之單位球心角，或稱之為單位立體角，該項所代 表的物理意義為局部平衡態成立時，此區域內朝向各個單位立體角的 輻 射 強 度 皆 必 須 相 同 。 dΩ μ 和η為 相 應 於 x 軸 及 軸 的 方 向 餘 弦 。 ( y

cos , sin cos

μ = θ η= θ φ)

2.3.2 邊界條件

( ) 0 , Iph 0 III Y y ∂ = = ∂ (2- 46) ( ) , Iph 0 IV Y L y ∂ = = ∂ (2- 47) 當奈米線內存有孔隙時，則給予孔隙邊界絕熱的假設。

2.3.3 熱導率計算

於局部非平衡區域內，若從溫度的角度來探討材料熱導率的尺寸 效應並不恰當，而應從熱通量的角度下來進行較為合理，則系統所涵 蓋空間中的任何一點在x方向的熱通量表示如下： ( , ) ( , ) x ph q x y =

∫

μI x y dΩ (2- 48) 上式描述了空間中的輻射強度Iph( , )x y 在x方向的淨能量通量即為x方 向的熱通量，其中μ為相應於x軸的方向餘弦。 在二維穩態系統中，x方向的熱通量將隨著 軸位置而變，因此y x 方向的平均熱通量可以表示為： 0 1 ( , ) L x x q q x y dy (2- 49) L =

_∫

則熱導率的計算方式如下： ( ) x H C q k T T L = − / (2- 50)

2.4

矽奈米線徑向二維導電系統

本文研究系統之二維傳電統御方程式、邊界條件及電導率和塞貝 克係數的計算方法，如下分述之。

2.4.1 電子運動方程式

處理粒子分佈 f_ε的問題，轉變成對電量I_e的探討，轉換如下： ( , ) ( ) e d I =

∫

v θ φ ef D_ε ε dε (2- 51) e I 被稱為電子輻射強度，其所代表的意義為單位時間單位面積單位能 量區間中，往單位立體角傳遞的電量，因此電子的輻射強度有別於以 能量為單位的聲子輻射強度。其中vd( , )θ φ 為電子的漂移速度(Drift velocity)，D( )ε 為電子的狀態密度，而 為其電量單位。 e 在本系統條件下，忽略電子於二維空間中的加速度效應，並以鬆 弛時間近似的概念來簡化電子複雜的碰撞機制，則二維穩態電子輻射 傳輸方程式表示如下： o e e e x y e I I I I v v x y τ ∂ ₊ ∂ ₌ − ∂ ∂ (2- 52) 其中 o e I 為平衡態的電子強度，τ 為鬆弛時間。若將平均自由徑引入上 式並展開平衡態電子輻射強度後，則上式可改寫為： 4 1 4 Ω − ∂ ₊ ∂ ₌ ∂ ∂ Λ

∫

e e e e I d I I I x y π π μ η (2- 53) 其中 為曲面所對應之單位球心角，或稱之為單位立體角，該項所代 表的物理意義為局部平衡態成立時，此區域內朝向各個單位立體角的 輻射強度皆必須相同。 dΩ μ和η為相應於x軸及y軸的方向餘弦。

2.4.2 邊界條件

二維穩態系統中，四個邊界條件將會直接地影響電子輻射強度於 二維平面上的分佈情形，本研究之邊界條件給予如下：

( ) , Ie 0 IV Y L y ∂ = = ∂ (2- 57) 其中x方向的邊界條件隱含了電場對電子的作用，y方向則使用不穿透 材料表面的彈性鏡面式反射假設，若是當電子撞擊奈米線內孔隙邊界 時，則以相同的散射模式來簡化處理。

2.4.3 電導率計算

根據電導率的定義為單位電場作用下通過控制表面的電流密度， 或稱之為電通量，則系統所涵蓋空間中任何一點往x方向的電通量將 表示成為： ( , ) ( , ) x e J x y =

∫

μI x y dΩ (2- 58) 上式描述了空間中的電子輻射強度I x ye( , )在x方向的淨電量通量即為 x方向的電通量，其中μ為相應於x軸的方向餘弦。 於二維穩態系統中，電通量同時為x和 的函數，而y x方向的平均 電通量可表示如下： 0 1 ( , ) L x x J J x y dy (2- 59) L =

_∫

則電導率的計算方式為： x J E σ = (2- 60)

2.4.4 塞貝克係數計算

塞貝克係數為描述材料內電勢與載子熱擴散動能相同時的能量平 衡位置，因此處於平衡位置以上的載子才能有效地於相同溫差下替材 料創造更高的端電位差，基於此概念時，則塞貝克係數的計算方式可

表2- 1 不同空間尺度下之傳輸理論分類[43]

Electrons Phonons Photons Fluids

Wave Quantum mechanics Quantum mechanics Maxwell EM theory Super fluidity Quantum Boltzmann equation Coherence theory Ballistic

transport Ray tracing

Free molecular flow

Boltzmann equation Radiative equation Boltzmann equation Ohm’s law Fourier’s law Newton’s shear stress Diffusion approximation Regimes Part ic le ( )c D<O A ( )_p D<O A ~ ( )c D O A ~ ( )_p D O A ( ) D<O Λ ~ ( ) D O Λ ( ) D>O Λ Coherence length, Phase-breaking length, Mean free path, o 100 A 10 Ao μm~km o 100 ~ 1000 A 100 A ~ 1kmo c A p ≥ Λ A Λ

表2- 2 能量載子特徵[43]

Free Electrons Phonons Photons Molecules

Source Propagation Media Statistics Frequency Velocity Freed from nucleic bonding Lattice vibration Electron and

atom motion Atoms In vacuum or media In media In vacuum or media In vacuum or media

Fermi-Dirac Bose-Einstein Bose-Einstein Boltzmann

0-infinite Debye cutoff 0-infinite 0-infinite

Metal Ceramic Si Y X 圖2- 1 奈米線徑向系統示意

Metal Ceramic Si Y X 圖 2- 2 多孔奈米線徑向系統示意

三

數值方法

隨著計算機工業的進步，對於隨機論性的(Stochastic)，或是無法 及不易求得解析解(Analytic solution)之確定論性的(Deterministic)問 題，皆獲得了新的解決辦法。隨機論性的模式需要從機率的本性進行 極大量的亂數運算來達到統計上的結果，如蒙地卡羅法(Monte Carlo method)；而確定論性的模式如分子動力學法(Molecular dynamics)或是 本研究依循的數值方法，必須先將系統變化規則的統御方程予以離散 成可進行疊代的模式之後，再配合初始和邊界條件來完全描述系統。 由於數值處理的過程中涉及方程式離散化、數值法穩定性以及計算機 有效位數的問題，因此藉由數值模擬得到的近似解必然包含截尾誤差 (Truncation errors)及捨位誤差(Round-off errors)兩種。另外，收斂條件 的設定也相當重要。

於本文研究系統中，統御方程的對流項採用有限差分法(Finite

difference method)，而積分項是以 SN近似法來處理，分別說明如下。

3.1 有限差分法

有限差分法應用於偏微分方程時，即是先將方程式欲求得的對象

以多項式(Polynomials)的方式表示，如泰勒級數展開 (Taylor series

expansion)後，刪除影響較小的高階項次並對時間及空間選擇保留的階 次，再將各微分項以前向差分(Forward difference)、後向差分(Backward difference)或中央差分(Central difference)處理。

1 i i f f df dx x + − = Δ +O

( )

Δx (3- 1) 其中O(Δx)為截尾誤差，將會隨著x方向的格點數增加而下降。顯然 地，前向差分的截尾誤差影響力能夠再次降階，如： 2 4 1 3 2 i i i f f f df dx x + + − + − = Δ +

(

)

(3- 2) 2 ( ) O Δx 後向差分： 當函數 f 僅為x的函數時，則 f 隨著x變化的一階後向差分方程可 表示為： 1 i i f f df dx x − − = Δ +O

( )

Δx (3- 3) 若截尾誤差為

(

_{( )}2

)

O Δx ，則 f 在x方向上的斜率可表示為： 1 2 3 4 2 i i i f f f df dx x − − − − = Δ +

(

)

2 ( ) O Δx (3- 4) 中央差分： 當函數 f 僅為x的函數時，則 f 隨著x變化的一階中央差分方程可 表示為： 1 1 2 i i f f df dx x + − − = Δ +

(

)

2 ( ) O Δx (3- 5) 上述三種基本的差分規則，考慮邊界條件對系統狀態變化方向上 的影響之後，再選擇合適的差分方法來進行狀態演化。前向、後向及 中央差分法如圖3- 1示意。

3.2 S

N

近似法

和 Fiveland 的研究，才首度將 D.O.M.應用於輻射熱傳領域[11]。 在輻射熱傳問題中，由於散射項(Scattering term)為積分形式，因 此輻射熱傳方程為一積微分方程式，SN近似法是以 個方向數來涵蓋 輻射強度於全立體角 N 4π 的分佈，每個方向都有其相對應的權重函數 (Quadrature weights, )，而三維輻射強度的大小則可由各方向的強度 乘上其權重函數後，再相加成得到： w 4 1 N i i i Id w π = Ω =

∑

I (3- 6)

∫

其中下標 代表輻射強度離散的數目，不同方向數所對應的方向餘弦 (Direction cosine)如表 3- 1所示。 i

3.3 奈米線徑向熱電系統

本文以輻射方程式同時描述電子及聲子於奈米線徑向熱電系統中 的運動行為，並利用數值計算的方式進行研究，因此首先須將輻射方 程予以離散成為差分方程後，再配合格點數設置來使計算機能夠執行 疊帶的動作。數值計算流程與計算使用的重要物性參數分別如圖 3- 2 和表 3- 2所示。

3.3.1

二維輻射傳輸差分方程式

輻射傳輸方程中，積分項採用 D.O.M.之 S4近似法來模擬三維空 間的輻射強度於二維平面上的變化情形，此時傳輸方程為四條聯立方 程組，再配合有限差分法處理後，由邊界條件來決定輻射傳遞方向。 第一象限(cosθ > 0,cos > 0φ ) x 方向及 y 方向皆以前向差分處理：

1 1( 1, 1) 2 2( 1, 1) 3 3( 1, 1) 4 4( 1, 1 4 w I x y w I x y w I x y w I x y ) π + + + + + + + + + + + = Λ 1( 1, 1) I x+ y+ − Λ (3- 7) 整理後可得第一象限輻射傳輸差分方程式： 1 1

cos sin cos 1

( 1, 1) 4 w I x y x y θ θ φ π ⎡ ⎤ + + _⎢ + − + _⎥ Δ Δ Λ ⎣ Λ⎦ 1 1

cos sin cos

( , 1) ( 1, ) I x y I x y x y θ θ φ = + + Δ Δ + 2 2( 1, 1) 3 3( 1, 1) 4 4( 1, 1 4 w I x y w I x y w I x y ) π + + + + + + + + + Λ (3- 8) 第二象限(cosθ < 0,cos > 0φ ) x 方向利用後向差分，y 方向以前向差分處理： 2( , 1) 2( 1, 1) 2( 1, 1) 2( 1, )

cos I x y I x y sin cos I x y I x y

x y θ + − − + + θ φ − + − − Δ Δ 1 1( 1, 1) 2 2( 1, 1) 3 3( 1, 1) 4 4( 1, 1 4 w I x y w I x y w I x y w I x y π − + + − + + − + + − + = Λ ) 2( 1, 1) I x− y+ − Λ (3- 9) 整理後可得第二象限輻射傳輸差分方程式： 3 2

cos sin cos 1

( 1, 1) 4 w I x y x y θ θ φ π ⎡ ⎤ − + _⎢− + − + _⎥ Δ Δ Λ ⎣ Λ⎦ 2 2

cos sin cos

( , 1) ( 1, ) I x y I x y x y θ θ φ = − + + − Δ Δ 1 1( 1, 1) 3 3( 1, 1) 4 4( 1, 1 4 w I x y w I x y w I x y π − + + − + + − + + Λ ) _{(3- 10)} 第三象限(cosθ < 0,cos < 0φ )

1 1( 1, 1) 2 2( 1, 1) 3 3( 1, 1) 4 4( 1, 1 4 w I x y w I x y w I x y w I x y ) π − − + − − + − − + − − = Λ 3( 1, 1 I x− y− ) (3- 11) − Λ 整理後可得第三象限輻射傳輸差分方程式： 4 3

cos sin cos 1

( 1, 1) 4 w I x y x y θ θ φ π ⎡ ⎤ − − _⎢− − − + _⎥ Δ Δ Λ ⎣ Λ⎦ 3 3

cos sin cos

( , 1) ( 1, ) I x y I x y x y θ θ φ = − − − − Δ Δ 1 1( 1, 1) 2 2( 1, 1) 4 4( 1, 1 4 w I x y w I x y w I x y ) π − − + − − + − − + Λ (3- 12) 第四象限(cosθ > 0,cos < 0φ ) x 方向利用前向差分，y 方向以後向差分處理： 4( 1, 1) 4( , 1) 4( 1, ) 4( 1, 1

cos I x y I x y sin cos I x y I x y

x y θ + − − − + θ φ + − + − Δ Δ ) 1 1( 1, 1) 2 2( 1, 1) 3 3( 1, 1) 4 4( 1, 1 4 w I x y w I x y w I x y w I x y π + − + + − + + − + + − = Λ ) 4( 1, 1 I x+ y− − Λ ) _(3-13) 整理後可得第四象限輻射傳輸差分方程式： 2 4

cos sin cos 1

( 1, 1) 4 w I x y x y θ θ φ π ⎡ ⎤ + − _⎢ − − + _⎥ Δ Δ Λ ⎣ Λ⎦ 4 4

cos sin cos

( , 1) ( 1, ) I x y I x y x y θ θ φ = − − + Δ Δ 1 1( 1, 1) 2 2( 1, 1) 3 3( 1, 1 4 w I x y w I x y w I x y π + − + + − + + − + Λ ) _{(3- 14)}

3.4 收斂條件

上節可知，二維輻射系統是由四條差分方程聯立之後解得空間中

全傳播於二維平面時，稱之為一次完整疊代。而強度終止計算前的兩 個疊代平面，必須滿足以下的相對誤差設定始達到收斂要求： 1 5 10 + − − ≡ < n n i i error n I I I ε i (3- 15) 由於本文所探討的的二維穩態輻射問題，於傳輸方程式中並無時 間項，因此在未滿足收斂條件的空間輻射分佈並不具有任何物理意 義，僅當各方向強度滿足收斂條件後，才開始一連串求得傳輸係數的 動作。

3.5

數值方法驗證

進行主題研究時，探討的方式可分為直接採樣的實驗法，或是基 於已知理論來描述系統動態的模擬法。以模擬法來說，又可分為直接 對方程式求得解析解或是以近似的方式逼近準確值，其中計算機在取 得近似解的能力上，對於學界及業界均有極大的貢獻。 本文研究亦即以計算機進行數值演算的方式來模擬熱電的傳輸行 為，然而此種研究方式產生的誤差，除了所依循科學理論的假設之外， 尚有前述之數值模擬面的不準度，在這種情況下，數值模擬將有必要 藉由其它的模式來佐證數值結果不隨著時間、地域及研究方法等有所 差異。

3.5.1

格點測試

數值模擬技術中，首先是將整個系統所涵蓋的區域切割成各式網

點將獲得較佳的系統動態解析度，然而格點過小也可能導致由計算機 本身有效位數的限制所產生的捨位誤差，這種狀況較不常見。在本文 系統中，採用正方形網格點來模擬輻射於二維空間的穩態分佈情形， 其中選擇三種不同格點數101x101、501x501 及 1001x1001 進行格點數 測試後發現，格點影響落於工程允許誤差 10％之內，因此文中將採取 利於計算的 101x101 格點數進行往後的研究，如圖 3- 3~4 所示。

3.5.2 結果比較

藉由數值計算來模擬系統的變化，除了必須進行前述之格點測試 外，與其它已發表於高影響力期刊的研究法相互比較，將有助於檢驗 研究模式的完備性，但不表示獲得相異的結果必定是遠離真實的。本 文以矽奈米線進行結果比較，熱電優值中包含三個物理性質，在熱導 率的部分分別與氣液固法(Vapor-liquid-solid method)成長出的矽奈米 線[28]和Li等[48]的模擬結果比較，而功率因子則是由小尺寸放大至近 三維結構時與塊材比對[47]，如圖 3- 5~6 所示。

表3- 1 SN近似法之不同方向數與其相應權重關係[11] Ordinates Weights S₂ (symmetric) S₂ (nonsymmetric) S₄ 0.5773503 0.5773503 0.5773503 0.5000000 0.5000000 0.2958759 0.2958759 0.9082483 Order of Approximation

ξ

η

μ

w

1.5707963 0.7071068 1.5707963 0.2958759 0.9082483 0.2958759 0.9082483 0.2958759 0.2958759 0.5235987 0.5235987 0.5235987 S₆ 0.1838670 0.1838670 0.1838670 0.6950514 0.6950514 0.9656013 0.1838670 0.6950514 0.9656013 0.1838670 0.6950514 0.1838670 0.9656013 0.6950514 0.1838670 0.6950514 0.1838670 0.1838670 0.1609517 0.3626469 0.1609517 0.3626469 0.3626469 0.1609517

表3- 2 數值計算參數[20,43,46,47] Values_300K Specific Heat (Volumetric) 107 Parameters 0.19 m_o Electrons Phonons Units

Drift Velocity_sat. Effective Mass_<100> Mobility

Mean Free Path Group Velocity

cm/s

kg

2

cm /V-s

o

A

m/s

3

J/m -K

1450 10 2604 9.3x105 1804

Backward Exact Forw_ard Central

ƒ

χ

圖3- 1 有限差分法示意

Start Initialization Discrete radiative transfer equation 1 5 1 10 n n n I I I + − − ≤ × End YES NO 1 n n I =I + 圖 3- 2 數值計算流程

Wire Thickness, L (nm) 1 10 100 1000 Thermal Conductiv ity , K (W /m-K) 0 20 40 60 80 100 120 140 101x101 grids 501x501 grids 1001x1001 grids 圖3- 3 熱導率格點測試

Wire Thickness, L (nm) 1 10 100 1000 P o w e r F a ct or, P (W /m -K 2 ) 0.00 0.01 0.02 0.03 0.04 0.05 101x101 grids 501x501 grids 1001x1001 grids 圖 3- 4 功率因子格點測試

BULK Wire Thickness, L (nm) 1 10 100 1000 Th er m a l C o n d u c ti v ity , K ( W /m -K ) 0 20 40 60 80 100 120 140 160 Present results Li et al. (2008) NUMERICAL RESULTS Li et al. (2003) EXPERIMENTAL DATA 圖 3- 5 熱導率驗證