探討聚亞醯胺與金屬層之介面破壞行為

國 立 交 通 大 學

機械工程學系

碩士論文

探討聚亞醯胺與金屬層之介面破壞行為

Investigating the interfacial failure between

polyimide and metal layer

研 究 生：王泰元

探討聚亞醯胺與金屬層之介面破壞行為

Investigating the interfacial failure between polyimide and metal layer

研 究 生：王泰元 Student: Tai-Yuan Wang 指導教授：蔡佳霖 Advisor: Jia-Lin Tsai

國 立 交 通 大 學

機 械 工 程 學 系

碩 士 論 文

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

June 2010

探討聚亞醯胺與金屬層

之介面破壞行為

研究生：王泰元

指導教授：蔡佳霖 博士

國立交通大學機械工程學系碩士班

摘 要

本研究著重在調查具有不同原子間作用力時 polyimide(PI)聚合物 與膠黏接合面金屬層之間的介面間破壞現象。藉由分子動力模擬，施 加負載於膠黏接合面(三明治)結構─將 PI 夾於上下兩層薄金屬板 間，進行拉伸測試，並討論三種不同程度的凡得瓦(van der Waals, vdW)

作用力下PI 與金屬板的破壞現象。結果指出當分子間的 vdW 作用力

較大時，PI 變得較為脆性，其降伏強度以及楊氏係數皆上升，而在三

明治介面間連結的部份，當 PI 表現出較脆的性質時，破壞過程(從裂

縫成長開始到結束)相當的短，幾乎為瞬間的破壞；而另外一方面， 當凡得瓦力較弱時，則表現出破壞過程較長的延性破壞。此外，在不

Investigating the interfacial failure between polyimide

and metal layer

Student : Tai-Yuan Wang Advisor : Dr. Jia-Lin Tsai

Department of Mechanical Engineering National Chiao Tung University

Abstract

This research aims to investigate the interfacial failure between the polyimide polymer and the metal layer associated with different atomistic interaction. Three different degrees of vdW interactions between the polyimide chains were considered in characterizing the mechanical behaviors of the bulk polyimide polymer as well as the interfacial bonding strength. The failure between the polyimide and metal sheet was evaluated by performing tensile loading on the designed sandwich structures with polyimide embedded with the thin metal sheet using molecular dynamic simulation. Results indicated that when the vdW interaction between the molecular chains is getting stronger, the polyimide is becoming brittle resulting in the increase of the yielding stress as well as the modulus. For the interfacial bonding, it is found that when the polyimide is brittle, the failure process (failure begin to failure end) is quite short and almost close to sudden failure. On the other hand, when the extent of vdW interaction is reducing, the associated failure mechanism is dominated by the ductile failure in which the failure process is longer. In addition, the moisture is taken into

the enhancement is more significant as the temperature is lower. In the future, the investigation will be focusing on the influence of moisture on the interfacial failure between the polyimide polymer and the metal sheet.

誌 謝

在這份畢業論文完成時，代表著兩年內的學習有了階段性的成 果，也意味著我在邏輯判斷上有所成長。在指導教授 蔡佳霖博士的 督導鞭策下，使得我走在研究的道路上能不偏離正確的方向，逐步的 朝向真理邁進。在此由衷的感謝老師平時的指導教誨，他認真的求學 精神，也帶領著實驗室的夥伴們進步。另外也要感謝世華學長，在兩 年內的諸多幫忙，運用其豐富的知識帶領我們解決問題；也謝謝與我 一起口試的孟哲，在很多事情上常常提醒我，不致於逾期而未處理； 以及碩一的學弟妹們，認真的處理實驗室的各項事務，使得我們不用 費心。最後也感謝我的父母親，讓我能夠在求學的過程中，沒有經濟 上的憂慮，單純的做研究，希望在未來的博士生涯中，也能夠有讓他 們驕傲的表現。

目錄

中文摘要 i 英文摘要 iii 誌謝 v 目錄 vi 圖表目錄 viii 第一章、緒論...1 1.1 研究動機與文獻回顧 ...1 1.2 研究方法...5 第二章、分子動力學基本概念...6 2.1 運動方程式與積分演算法 ...6 2.2 週期性邊界與截斷半徑 ...7 2.3 系綜(Ensemble)與溫度修正方法 ...8 2.4 原子尺度應力計算 ...9 第三章、PI 與膠合面結構 ...10 3.1 勢能函數之選用以及模型建構 ...10 3.2 PI 之機械性質...15 3.3 膠合面結構之拉伸測試 ...21

4.1 水分子模型 TIP4P/2005 ...30 4.2 PI 之水氣效應...32 第五章、總結...36 第六章、未來工作...37 參考文獻...41 附表...46 附圖...56

圖表目錄

表1 分子動力學之系綜[11]...46 表2 原子代號與質量電荷表[13,20]...47 表3 延展鍵結勢能參數表[13] ...47 表4 角度勢能參數表[13] ...48 表5 扭轉勢能參數表[13] ...48 表6 高分子凡得瓦勢能LJ 參數表[13] ...49 表7 水分子模型與金屬板凡得瓦勢能 LJ 參數表[12,43]...50 表8 各方向楊氏係數比較表...51 表9 不同凡得瓦力之PI 密度比較表 ...51 表10 不同凡得瓦力之PI 機械性質比較表 ...52 表11 不同介面強度之膠合面破壞強度比較表...52 表12PI 密度比較表 ...52 表13 乾燥PI 之勁度矩陣與楊氏係數表 ...53 表14 含水PI 之勁度矩陣與楊氏係數表 ...54 表15 水氣效應於PI 之楊氏係數 ...55

圖2.2.1 截斷半徑示意圖 ...56 圖3.1.1PI 單體 ...56 圖3.1.2 鍵結能量示意圖 ...57 圖3.1.3 碳原子之凡得瓦勢能與凡得瓦力...57 圖3.2.1PI 初始結構 ...58 圖3.2.2PI 平衡溫度能量圖 ...59 圖3.2.3PI 之負載 ...59 圖3.2.4 三方向應力與模擬時間圖...60 圖3.2.5PI 之應力應變圖(1ΕPI) ...60 圖3.2.6PI 之應力應變圖(5ΕPI) ...61 圖3.2.7PI 之應力應變圖(10ΕPI)...61 圖3.2.8 不同凡得瓦力之應力應變比較圖...62 圖3.3.1Z 方向密度分佈計算方式示意圖 ...62 圖3.3.2 秩序參數向量示意圖：(A)PI 上的代表向量 (B) 參考座標Z 軸與代表向量所夾的角度...63 圖3.3.3 膠合面結構密度Z 方向分布圖 ...63 圖3.3.4 膠合面結構秩序參數Z 方向分布圖 ...64 圖3.3.5 應變定義示意圖...64

圖3.3.7 金屬所造成之凡得瓦力計算...65 圖3.3.8 膠合面結構之應力計算示意圖 ...66 圖3.3.9 上下層金屬板之凡得瓦吸斥力與應變之關係圖 ...66 圖3.3.10 上下層金屬板之凡得瓦力與應變之關係圖 ...67 圖3.3.11 平移後之上下層金屬板吸斥力與應變關係圖 ...67 圖3.3.12 正常凡得瓦力之應力應變關係圖 ...68 圖3.3.13 正常凡得瓦力之平均應力與 VIRIAL STRESS比較圖...68 圖3.3.14 計算VIRIAL STRESS時 Z 方向有效長度示意圖...69 圖3.3.15 膠合面結構之拉伸試驗表面破壞圖...69 圖3.3.16 介面強度破壞過程與應變關係圖(5ΕPI) ...70 圖3.3.17 介面強度破壞過程與應變關係圖(10ΕPI) ...70 圖4.1.1TIP4P/2005 水分子模型 ...71 圖4.2.1 含水PI 初始結構 ...71 圖4.2.2 含水 PI 平衡結構...72 圖4.2.3 實驗之PI 水氣效應與溫度關係圖[4] ...72 圖4.2.4 本篇模擬之PI 水氣效應與溫度關係圖 ...73

第一章、緒論

1.1 研究動機與文獻回顧

奈 米 複 合 材 料(nanocomposite) 具 有 極 好 的 比 強 度 (specific strength)、耐高溫、耐腐蝕等特性，在近十年內已備受重視，其中以 高分子為核心、膠黏劑的三明治、膠合面結構的使用量因廣泛運用在 不同的領域(汽車、航太等)而急速的成長，其主因即是因為相較於其 他的結構，三明治具有相當高的撓曲勁度(flexural stiffenss)與重量之 比率、較低的側向形變、以及較好的抗彎性。大多數的三明治面結構 具有兩個表面層，而材料與厚度根據其首要目標：抵抗平面上(in-plane) 以及側向的負載來決定。在假設核心材料為單一均勻時，通常將三明 治視為一對中平面(mid-plane)對稱的結構[1]。這類的三明治結構其中 一種運用，即是做為海上航行之遊艇、漁船之建材，但在長期處於潮 濕環境下、承受海水侵蝕等惡劣環境下，其強度相較於在乾燥環境 中，理應有所不同；儘管一般應用在車輛、列車以及船隻上的材料外 部均有防水處理，但是在擦撞破損、外表掉漆等情況出現後，水氣依 然會對其產生影響。Kolat等人[2]探討了在船舶相關用途上，不同的

系統是由強化玻璃纖維(Glass Reinforced Plastics, GRP)構成的表面 層，分別搭配多種不同的核心材料，藉以討論不同的核心材料下，所 受到海水侵蝕，對於其破壞強度的影響。首先將材料放置在50 ℃， 含有5 % 氯化鈉(NaCl)的蒸氣中，長達72小時，模擬受到水氣的侵 蝕，再利用島津萬能材料試驗機進行Mode-I CSB(Cracked Sandwich Beam)測試。試驗的結果發現：木質(wooden)核心以及夾板(plywood) 核 心 之 斷 裂 韌 性(fracture toughness) 是 下 降 的 ， 而 以 core-mat 、 polyurethane為核心的三明治則略為上升，作者量測了以上材料之斷 裂韌性，但是沒有深入探討這些材料所受到海水侵蝕之機制，而這也 是本篇所感興趣的─不同的親疏水性材料，受到水氣在化學或物理反 應的影響所導致的變化。除了三明治、膠合面(adhesive joint)結構外， 當水氣攝入於複合材料後，對於高分子基材(matrix)與強化纖維間連 結強度的影響，亦是學者廣泛研究的重點。Selzer以及Friedrich[3]， 調查了水氣對於以碳纖維強化之高分子複合材料所造成之機械性 質、破壞行為的影響。其中發現以環氧樹脂(Epoxy)為基礎之複合材 料 ， 僅 在 負 載 方 向 垂 直 於 強 化 纖 維 方 向 時 ， 楊 氏 係 數(Young’s modulus)、拉伸強度(tensile strength)以及破壞前之伸長量有明顯下 降，其原因為負載方向如與強化纖維平行時，依然由纖維承受主要的 負載，而纖維本身幾乎不受水氣的影響；另外當水氣吸附於高分子基

材內時，不僅造成了高分子軟化的現象使得楊氏係數下降，同時也減 弱了纖維與基材接觸面的連結，在此影響下使得破壞前之伸長量亦減 低。Athijayamani等人[4]針對天然纖維與聚酯(polyester)之複合材料， 控制纖維的重量比例、纖維之長度以及含水氣量進行拉伸、彎曲、以 及衝擊試驗。當複材處於乾燥狀態下，隨著纖維的重量比例、長度之 提升，拉伸、彎曲模數皆有明顯上升；一旦材料被浸置於水中後，提 高纖維的重量比例、長度所造成機械性質上升的現象則不明顯；此一 水氣所產生的影響，亦歸因於大量進入的水分子使得纖維與高分子基 材間黏合效果降低，直接導致機械性質減弱。 在介紹完上述受到水氣影響的複合材料，大致上可將水氣的侵 蝕，理解為連結強度的降低，不論在表面與核心材料之間，或是基材 與強化纖維之間均是如此。那麼關於水氣的攝入對於高分子本身的影 響，是否也是一味的削弱其強度呢？ Bernier以及Kline[5]對聚亞醯胺 (Polyimide, PI)─Vespel 材 料 所 進 行 的 研 究 ， 將 未 經 任 何 處 理 (as-received material) 的 Vespel 分 成 三 組 ， 對 第 一 組 僅 進 行 烘 乾 (dried) ； 對 第 二 組 材 料 首 先 在 溫 度 範 圍 80 K ~ 800 K 內 進 行 DMP(Dynamic mechanical properties)測試，紀錄數據，之後持續放置 在真空狀態、100℃狀態下烘烤50小時，成為乾燥材料並且進行同樣

50℃狀態下50小時後，變為含水氣量約1.37 wt%之含水材料，亦在同 樣溫度範圍內進行DMP並紀錄。第三組材料則維持原樣，不作任何 處理(as-received)。由實驗結果發現，在室溫時，第二組的乾燥材料 密度較未處理時下降0.8%，而成為含水材料後之密度則較乾燥時上升 了0.7%；並且在低溫環境下(低於230 K)，乾燥材料之楊氏係數亦明 顯低於未處理時，將此乾燥材料放置於潮濕環境後，所包含的1.37 wt%的水氣使得楊氏係數回復至接近未處理材料。這一系列的實驗中 指出：在低溫範圍，材料的楊氏係數隨著水氣的吸附而明顯增加。這 意味著在低溫時，水氣並非視為塑化劑(plasticizer)[6]來軟化材料，反 而促使結構更加的堅硬，此一情況作者猜測可能緣於水分子與PI結構 中 的 羰 基 內(carbonyl group) 之 氧 原 子 產 生 強 烈 的 氫 鍵 (hydrogen bond)，但是當溫度提升至270 K以上時，大部分的氫鍵斷裂，使得未 處理材料、烘乾後材料、以及含水氣材料三者之楊氏係數顯得非常相 近。 這個實驗結果十分有趣，不僅告訴我們水氣可能導致反塑化作用 [7]，使得結構更緊密，並且顯示出，水氣的效應是隨著溫度變化的， 同時與氫鍵的生成與破壞有關聯；因此本篇主要模擬的對象即是水氣 效應對於高分子與膠合面(類三明治)結構所產生之影響。

1.2 研究方法

本研究使用分子動力軟體DL_POLY2.15 版[8]，進行模擬奈米尺 度的PI 高分子材料、金屬─PI─金屬之膠合面結構，受到微小量水氣 攝入時，對於物理性質、機械行為的影響。第二章將介紹於分子動力 學中所運用到的基本概念，第三章將運用分子動力學以適當之力場建 構 PI、膠黏接合面的模型，首先探討 PI 本身，在凡得瓦力變化下， 所造成楊氏係數上的變化，並針對膠合面結構的破壞過程進行分析， 並討論隨著較大凡得瓦力所產生的結構脆化；第四章討論水分子模型 的選用，並且模擬2 wt%水氣的攝入對於 PI 本身在楊氏係數的影響， 並且模擬其溫度相依性。未來將試著把水氣導入膠合面結構，模擬水 氣效應對其拉伸測試的影響。

第二章、分子動力學基本概念

2.1 運動方程式與積分演算法

分子動力學之理論由 Irving 與 Kirkwood[9]率先於統計力學範疇 所提出，其主要概念為藉由勢能函數(potential function)主導各分子間 之作用力，並且假設所有分子運動遵守古典力學之牛頓運動定律如 (2.1.1 式)所示。 j i j i N 1 j 2 2 i i r r r r dr (r) dΦ dt r d m F − − − = =

∑

= (2.1.1) 其中Fi為第i 顆原子所受到的力，mi為第 i 顆原子之質量，Φ 為勢能 函數，勢能函數主導了原子所受到的力，對於整個模擬而言是影響結 果的重要角色，關於本篇所選用的各種勢能函數將於 3.1 節詳述。ri 與 rj分別為 i、j 原子之座標。藉由(2.1.1)式得到時間點 t 時原子上之 受力，以及已知的時間點 t 時之原子座標，即可運用 Verlet leap-frog

algorithm[8]。 Verlet leap-frog algorithm 藉由已知時間點 t 時的原子所

受到的力F(t)，將 F(t)除以原子質量得到時間點 t 之加速度，如(2.1.2) 式，選定時間步階Δt，從 t 2 1 t− Δ 時之速度推算出 t 2 1 t+ Δ 時之速度，進

m (t) F Δt Δt) 2 1 (t v Δt) 2 1 (t v + ← − + _(2.1.2) Δt)Δt 2 1 (t v (t) r Δt) (t r + ← + + _(2.1.3)

2.2 週期性邊界與截斷半徑

Rapaport[10]在針對分子動力模擬的著作中提到，對於如何使用 相對小的有限系統來忠實的求出近似無限系統之解，這個問題並沒有 唯一答案，而週期性邊界條件即是其中一種可能，它是分子動力學在 模擬有限的原子數時，為了模擬在真正巨觀尺度下，僅有極小部分的 原子會處在接近表面(surface)的狀態而進行的假設。由主胞室(Primary Cell)週遭環繞的鏡像胞室(Image Cell)，也會與主胞室之原子互相影 響，以此避免表面效應(surface effect)。在模擬系統中，所考慮的非鍵 結作用力(non-bond interaction)，當原子間相距一定距離以上時已近似 於零，為了減低計算所需耗費的時間，在計算非鍵結作用力時通常會 給予一指定的截斷半徑(cutoff radius)，與目標原子之距離大於截斷半 徑之原子則不與考慮，如圖 2.2.1，本篇對於凡得瓦勢能以及靜電勢 能所採用之截斷半徑皆為10Ao 。

2.3 系綜(Ensemble)與溫度修正方法

在 分 子 動 力 模 擬 的 微 觀 系 統 中 ， 模 擬 胞 室 為 孤 立 系 統(isolated system)，胞室內的原子、能量均不能離開此系統，而為了控制系統內 之原子、能量、溫度、系統體積、壓力等要素，分為諸多系綜，如表 1 所列，為系綜之分類以及特性[11]。由 2.1 節所提的演算法計算出的 速度，將會直接藉由動能與溫度之關係(2.3.1)式而反應在溫度上，Ek 即為動能(kinetic energy)， mi為:第 i 顆原子之質量，vi為:第 i 顆原子 之未調整速度， ' i v 為第 i 顆原子之調整後速度，N 為總原子數量， k 為波茲曼常數(Boltzmann constant)，T 為由未調整速度所計算出之 系統溫度，Td為目標溫度。由計算所得之溫度 T 與所系統設定溫度 Td相異，為了使系統維持在Td，所使用的方法即是速度調整(Velocity rescaling)[18]。由(2.3.2)式將調整後的速度作為計算動能之依據，進 而達到穩定之溫度Td。 NkT 2 3 v m 2 1 E n 1 i 2 i i k =

∑

= = (2.3.1) T T v v d i ' i = (2.3.2)

2.4 原子尺度應力計算

在原子規模時，系統所受到的應力可以定義為如(2.3.3)式所示， 稱為 Virial stress[12]σ_ij為系統胞室的應力張量、V為系統胞室之體 積、_mα_{為第α 顆原子之質量，} α i v 為第 α 顆原子之 i 方向速度分量、 α j v 為第α 顆原子之 j 方向速度分量、 αβ i F 為α 與 β 原子間之 i 方向作用力 分量、

r

jαβ: α 與 β 原子間之 j 方向距離分量。此式前項為動能，後項 為分子間之作用力乘與相距距離所得之能量，將兩項相加再除以系統 之模擬胞室體積，可視為將功分布於空間中，即為能量密度(Energy density)之概念[8]。

∑

∑

≠ + − = α β α αβ j αβ i α j α i α ij F r ) 2 1 v v (m V 1 σ _(2.3.3)

第三章、PI 與膠合面結構

3.1 勢能函數之選用以及模型建構

本篇膠合面結構由PI黏著劑(adhesive)與兩層金屬板─被黏附物 (adherend)所組成，其中PI所採用之勢能為Dreiding勢能函數，是於 1990年由Mayo等人[13]所提出，作為分子動力模擬的通用勢能函數； 其簡單的形式為其優點；也因此，學者們在模擬特定的物理行為時， 可以適當的調整其中的參數，以期能夠表現出適當的特性，本篇所使 用之原子代號以及質量列於表2，而PI之化學式如圖3.1.1所示，選用 的鏈長(chain length)為10個重複單體。本篇PI所考慮之勢能函數可以 表示如(3.1.1)式，總能量可分為鍵結能(bond energy) (3.1.2)式與非鍵 結能(non-bonded energy) (3.1.3)式。

Utotal = Ubond + Unon-bonded (3.1.1)

Ubond = Ustretch + Uangle+ U tosion + Uinversion + Unon-bonded (3.1.2)

鍵結能如(3.1.4)式~(3.1.7)式、圖 3.1.2 所示，包含以下數種能量： (1) 延展鍵結勢能 Ubond：由兩顆原子間的距離變化，而產生的能量， 其中kr為勁度常數(stiffness constant)，其參數如表 3 所列。 2 0 r stretch k (r r ) 2 1 U = − _(3.1.4) (2) 角度鍵結勢能 Uangle ：由三顆原子所構成之角度，角度變化時產 生之能量，其中kθ為能量參數，各參數為表4 所列。

( )

( )

(

)

2 0 θ

angle k cos θ cos θ

2 1 U = − _(3.1.5) (3) 扭轉勢能 Utorsion :由四顆相鄰原子所組成的兩個平面，兩平面夾 角變化時所產生之能量，其中Aτ為能量參數，各參數如表5 所列。

(

)

[

1 cos mτ δ

]

A U_torsion = _τ + − (3.1.6) (4) 反向勢能 Uinversion :由四顆相鄰原子所組成，其中一顆原子與另外 三顆所組成之平面所夾之角度變化而產生之能量，其中 Aφ為能量參 數，本篇之所有Aφ皆採用40 Kcal/mol-rad2。

( )

[

− ϕ

]

= A_ϕ1 cos U_inversion (3.1.7)

⎪ ⎪ ⎩ ⎪⎪ ⎨ ⎧ ≥ < ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ − ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ = cutoff ij cutoff ij 6 ij 12 ij vdw r r 0 r r r σ r σ 4ε U _(3.1.8) ij j i ij 0 tic electrosta r q q ) ε(r 4π 1 U ε = _(3.1.9)

非鍵結能的部份，如(3.1.3)式，可以分為凡得瓦勢能(van der Waals potential)與靜電勢能(electrostatic potential)。凡得瓦勢能所採用的是 Lennard-Jones(LJ) potential[14]，如(3.1.8)式，此勢能為原子間距離之 函數，對其微分可得凡得瓦力；分為前項排斥力(repulsive force)與後 項吸引力(attractive force)兩部份,其中 ε 為能量參數(Energy well)，σ 為使能量為零之距離；LJ potential 中的參數僅考慮同一種原子間 [14]，例如 Dreiding 勢能裡，碳與碳間的能量參數為 0.0951 Kcal/mol， σ 為 3.473 Å，圖 3.1.3 即為碳原子之凡得瓦勢能與凡得瓦力與原子間 距離之關係；當需要考慮碳對於其他原子時，如碳與氫之間之凡得瓦 勢能，則必須藉由結合規則來概算。Lorentz–Berthelot combining rule[15]是一種廣為使用的結合規則，包括 Dreiding、Amber[16]等勢 能 皆 採 用 此 規 則 。 利 用(3.1.10) 式 、 (3.1.11) 式 ， Lorentz–Berthelot combining rule 可概略的推算出不同原子間的凡得瓦參數，其方法 為：將兩個不同種類原子 i 與 j 各自之零勢能距離 σii與 σjj相加取平

均作為此兩種分子間之零勢能距離σij；將兩種不同原子i 與 j 的各自 能量參數 εii與 εjj相乘開根號作為此兩種分子間的能量參數 εij。本篇 之水分子裡的氧原子對於其他考慮凡得瓦勢能之原子以及金屬板對 於其他考慮凡得瓦勢能之原子，皆遵循此規則；本篇使用之參數列於 表6、表 7。

2

)

σ

(σ

σ

ii jj ij

+

=

_(3.1.10) jj ii ij

ε

ε

ε

=

_(3.1.11)

而本篇所採用之金屬，其結構為面心立方(Face Centered Cubic, FCC)堆積，如圖 3.1.5。由於表面層之材料強度往往高於核心材料數 倍以上，為了計算上的簡化，本模擬不討論金屬板內的現象，將金屬 固定凍結(freeze)[8]於空間中，僅考慮金屬對於高分子所造成的凡得

瓦力；其中所選用的金屬原子凡得瓦參數參考Adnan 與 Sun[12]之數

值，εmetal為0.1668 Kcal/mol ，而 σmetal為 2.315 Å。於第三章中所討

論之 PI 以及膠合面結構，為了簡化模擬過程而不考慮靜電勢能，僅

顆原子的部分電荷(partial charge)，本篇作法參考諸多學者[17-19]之方 法 藉 著 Discovery tool (Materials Studio 4.0, Accelrys) ， 選 用 CVFF(Consistent-Valence Force Field)[20]此種能依照原子種類以及鍵 結決定部分電荷的半經驗勢能(empirical potential)。關於靜電力的計 算使用Coulomb Sum with Distance Dependent Dielectric[8]，如(3.1.9)

式，其中ε0 ≈8.854187817⋅10-12 (F/m)為真空電容率[21]，因為靜電 力能造成長距離的影響，所以使用庫倫勢能(Coulomb potential)時，輔 以距離函數之電介質(dielectric)ε(r_ij)，可以加速庫倫勢能的收斂，其 中 ij * ij) ε r ε(r = ，ε*為常數，可依預定的收斂速度做調整，本篇為1.0。

3.2 PI 之機械性質

此小節將針對 PI 模型本身之楊氏係數進行模擬，首先將六條各

自具有十個重複單體(monomer)，以 Dreiding 勢能建構之 PI 鏈放置於

邊長為 200 Å 之立方體模擬胞室中，如圖 3.2.1，於三方向皆考慮週

期性邊界，時間步階(time step)選用 0.5 fs(femtosecond, fs)，第一步驟

於800 K 的環境下，採用 NPT 系綜，壓力設定為 200 MPa、持續 500 ps (picosecond, ps)，藉著壓力將結構壓縮至緊密狀態，接著第二步驟 同樣在800 K，採用 NPT 系綜，壓力設定為 0.1 MPa (一大氣壓)、持 續500 ps，讓結構在經歷第一步驟之壓縮後，有自由鬆弛的空間，而 後第三步驟則在NPT、0.1 MPa 狀態下，降溫到 600 K，持續 500 ps， 第四步驟則將結構於NST 系綜、0.0 MPa 下降溫到 300 K，最終第五 步驟保持在NST、0.0 MPa 下，降溫至絕對零度 0 K，圖 3.2.2 顯示 在平衡過程中，結構的溫度以及能量變化，最後在0 K 達到穩定，此 時的密度約為1.295 g/cm3，尚在合理範圍內[22-24]。在 Wei[25]運用

MD 探討奈米碳管(carbon nanotube, CNT)埋入 Polyethylene (PE)之結

構附著性的文章中提到，複合材料相較之楊氏係數於PE 本身有提升

質，而在較高溫的情況下其強度則大幅減弱。作者認為當 CNT 埋入 PE 時，藉著凡得瓦力，CNT 吸附了大量的 PE 在表面(遠超過其體積 比率所應該容納的數量)，產生了強化的作用，複合材料不論在小應 變或大應變時，CNT 的應變量均小於 1%，說明了鍵結勢能的貢獻相 當少；而溫度提升所導致的削弱現象，則是因為熱膨脹發生時，使得 CNT 表面的分子密度大幅下降，連帶使得凡得瓦勢能所能影響的原 子數下降，最終導致強度下降。前文提到本研究中將 PI 的結構在 0 K 穩定平衡，即是為了達到略去溫度影響之目的，單純討論 0 K 下結構 之應力應變曲線圖。Wei 所討論的 PE-CNT 複合材料相較於純 PE 本 身，在 50 K 時楊氏係數提升為 3.2 倍，並且點出了凡得瓦力對於機 械性質的影響；本模擬所討論的PI 將在 3.3 節做為膠合面結構的膠黏 劑材料並在 4.2 節考慮靜電力並置入水分子，討論水氣所產生的影 響；在不增加模擬的複雜度的前提下，本篇論文試著改變純 PI 本身 的凡得瓦能量參數 εPI(包括碳、氫、氧、氮這四種原子所各自具有的 能量參數εC、εH、εO、εN

，

以及藉著結合公式所算出的其他項)分別放 大為 5εPI 以及 10εPI，來討論凡得瓦力對於應力應變關係的影響；而 PI 與金屬板之間的介面黏著強度 εinterface 因為藉著前述之結合法則可

表示成：ε_interface= ε_metalε_PI 之式子(其中 εmetal 金屬之勢能參數維持原

5εinterface以及 10 εinterface，將在 3.3 節中討論不同的 εPI下所造成的變 化。在得到穩定的 0 K 無預應力的純 PI 結構後，接著進行楊氏係數 的探討所受到不同凡得瓦力的影響，至於應力率(stress rate)或應變率 (strain rate)則不在本篇討論的範圍。施加的負載情況如圖 3.2.3 所示， 本模擬所運用的第一個方法為在NST 系綜下，施加壓力於材料 Z 方 向，並觀察在此壓力下Z 方向的應變─其定義如(3.2.1)式，將長度變 化量(拉伸後之長度減去原長)除以原長。圖 3.2.4 顯示當模擬時間為 50 ps 時，在指定的 Z 方向已經達到穩定的目標應力 20 MPa，而 X 與 Y 方向則無受力。將目標應力逐步上升(每 50 ps 提高拉伸應力)，並 記錄下每次拉伸的應變量即可得到應力應變關係圖 3.2.5，在 5% 應 變之前，目標應力每次提升 20 MPa，曲線的上升呈現幾乎線性的狀 態，但此斜率在 5% 應變之後小幅度下降，且當 240 MPa 時可以看 到有較大的應變量，可視為降伏的現象發生，為了能觀察此部份的變 化，降伏後每次僅提升10 MPa。 0 0 L L L Strain = − _(3.2.1) 將 5 %應變前的區域視為線性區，並以前五筆(應力值以及相對

另外亦對X、Y 方向利用 20 ~ 100 MPa 的應力應變關係計算而得的 X、Y 方向楊氏係數列於表 8 中，可以發現我們雖然假設此高分子材 料為等向性，但是各方向的楊氏係數仍有不少的誤差，方法一在此先 針對Z 方向來討論。而第二個計算其楊氏係數之方法，從前述的五個 平衡步驟得到無應力PI 平衡結構後，給予材料 X 方向 1 % 均勻應變 (uniform strain)(ε11 =0.01、ε22 = ε33 = 0)，並於 NVT 系綜下平衡；藉由 平衡後三方向的系統應力值進行 C 矩陣之計算；將這些施加 X 方向 應變後所產生之三方向應力σ1、σ2、與σ3，藉由應力應變關係[27]將 這些值帶入(3.2.2)式中可求得勁度矩陣(stiffness matrix)中之 C11 、 C21、與 C31。再分別依照此流程施加應變於Y、Z 方向，即可得到完 整的勁度矩陣。 ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎣ ⎡ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎣ ⎡ = ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎣ ⎡ 33 22 11 33 32 31 23 22 21 13 12 11 33 22 11 ε ε ε C C C C C C C C C σ σ σ (3.2.2) 藉由上述流程所得到之C 矩陣，運算反矩陣 S，如(3.2.3)式進而 求得 E1、E2以及 E3；為了滿足本篇所假設：無結晶型(amorphous)的 PI 為等向性材料(isotropic material)[28, 29]，因此將三方向楊氏係數之 總和取平均視為本篇模擬之PI 之楊氏係數。

⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ − − − − − − = 3 2 23 1 13 3 32 2 1 12 3 31 2 21 1 E 1 E ν E ν E ν E 1 E ν E ν E ν E 1 S (3.2.3) 藉由方法二求出的楊氏係數約為在3.89 GPa，與方法一的三方向平均 楊氏係數3.84 GPa 相當接近，代表在本模擬中，藉由使用 NST 系綜 (方法一)下，對模型進行拉伸，利用所得到結構的應變量與所控制的 目標應力來計算之三方向平均楊氏係數，吻合於利用 NVT(方法二) 系綜，先給予結構指定應變量再平衡以求得該狀態下之各方向應力， 分別實施於三方向而求得之平均E 值。另外關於材料的降伏強度推算 方法，在此僅對 Z 方向拉伸應力應變關係做討論；一般在對陶瓷材料 約使用0.05 % offset，金屬材料為 0.2 % offset，而在這裡使用一般針 對高分子材料的2 % offset method[30-32]，如圖 3.2.5 即是將約 5%應 變前的斜率近似線，橫軸平移2 %應變後與應力應變曲線之交點，約 225 MPa 視為材料的降伏強度，較文獻值[26](127 MPa)高出許多，因 為此處降伏強度的準確度並非本模擬重點，所以仍將此值紀錄為1εPI 的 PI 結構之降伏強度。從前述文章中提到這個階段已經完成了以兩 種方法計算 1εPI 的 PI 結構之楊氏係數，並且估算出其降伏強度，接

時之密度列於表9，並將藉著每次於 Z 方向提升 20 MPa 應力之 5εPI、

10εPI各自的應力應變圖繪於圖 3.2.6、圖 3.2.7。可以發現僅放大凡得

瓦力就可以讓密度有大幅上升的趨勢，反映著結構的在較高的凡得瓦

力下變得緊實；繼續看到圖3.2.6 中，5εPI的應力應變曲線在拉伸過程

中，2 %應變以前大致維持線性，僅在 440 MPa 應力值時有較多的應

變，同樣利用 2 % offset method，可以推估其降伏強度約為 690 MPa；

而圖3.2.7 中 10εPI則呈現近乎完美的斜直線，從5εPI 的降伏強度相較 於1εPI來看，估計10εPI之降伏強度應該相當高，在此先不進行計算。 統計不同 εPI結構之應力應變曲線的線性區前五個應力值以及相對應 的應變資料，來計算並且平均而得的 Z 方向楊氏係數列於表 10，可 以看到楊氏係數隨著凡得瓦力參數的調整放大，也有相對應的放大趨 勢，這個結果也顯示出在此種應力應變關係計算下，非鍵結的凡得瓦 力即是決定高分子楊氏係數的主要因子。在探討完不同凡得瓦作用 力，對於PI 本身所造成的機械性質影響後，在下節將進行以 PI 為膠 黏劑材料的膠合面單軸拉伸試驗模擬，同樣的，PI 也分別於這三種凡 得 瓦 力 狀 態(1εPI、5εPI、10εPI ) 下 ， 與 金 屬 之 間 考 慮 介 面 強 度 PI metal erface int

ε

ε

ε

=

_{(εmetal 不變)，並且比較這些凡得瓦力狀態下對破} 壞過程所產生的影響。

3.3 膠合面結構之拉伸測試

在探討水氣效應對於膠合面結構之拉伸強度的影響之前，在此節 中先對乾燥結構的拉伸過程進行分析。將六條以 Dreiding 勢能建構 好之PI 放置於模擬胞室中，X、Y 方向考慮週期性邊界，將兩層金屬 板分別置於Z 方向邊界內之兩端，藉著與高分子間凡得瓦力之作用， 當上下兩層金屬板朝彼此互相接近時，可將置於其中的高分子壓縮為 膠合面材料之核心部份─膠黏劑，得到膠合面結構後，為了避免核心 部份之高分子皆為平行(層狀)排列，在 NPT 系綜下，壓力設定為 0.0 MPa，分為四個步驟，在 800 K 下持續 500 ps，讓高分子鏈可以充分 運動，接著分別持續500 ps 降溫至 600 K、300 K，最後到 0 K，經由 這四個程序，得到在0K 溫度下的結構，在建構膠合面結構時，上下 金屬板在Z 方向逼近黏著劑材料，有可能壓縮或藉著凡得瓦力拉伸高 分子，因此在進行拉伸測試前，會先針對金屬板的位置進行微調，使 得整體結構Z 方向的 Virial Stress 近似於零，藉以得到預應力極小之 膠合面初始結構。如圖 3.3.1，其 X、Y 方向邊長約 34 Å，Z 方向邊 長約 44 Å。在此溫度下，對於平衡好之結構進行物理性質的分析， 分別為密度以及秩序參數(order parameter)。密度的計算如圖 3.3.1，

Å) ，計算該區間內之高分子原子質量總和(不考慮金屬的質量)再除 以區間體積。而秩序參數之定義為(3.3.1)式[33]，如圖 3.3.2 (a)所示， 將較為堅硬(幾乎維持平面)的胺基與苯環選為代表 PI 之向量，而藉著 這個向量與向量 Z 軸間的夾角變化 (圖 3.3.2 (b))，可以知道 PI 在空 間中排列的情形，如平行Z 軸則秩序參數為 1，如垂直於 Z 軸則秩序 參數為 -0.5，而秩序參數越接近 0 則代表排列越不規則；而與密度相 同沿著Z 方向分別在各個區間內計算。

(

3 cos φ 1

)

2 1 (f) parameter order ₌ 2 ₋ (3.3.1) 由圖3.3.3 可發現在此三種不同強度的 εPI下，膠合面結構之密度，在 靠近金屬板之接合面附近，皆因為受到凡得瓦吸引力使得結構相當緊 密、有較高的密度，但是差距並不明顯，同樣的從秩序參數來看，圖 3.3.4 也展示出，在接合面處，高分子被吸附而緊貼在金屬面上與 Z 軸垂直，使得秩序參數接近-0.5，而遠離金屬板影響之中間高分子部 份，其秩序參數則較無規則，呈現亂數分布的情況。將 1εPI、5εPI、 10εPI各自的每個區間密度總和取平均分別為1.24 g/cm3 、1.35 g/cm3、 1.37 g/cm3，皆在前述表 8 所列之合理範圍內，與純 PI 本身的密度增 大現象一致；而在三個 εPI狀態下，它們在接合面處之密度與秩序參

做了調整(介面強度 εinterface是當εPI放大後再藉著結合法則所求出，並 不做直接更動)，所以並不會在與金屬之接合面處隨著 εPI、εinterface的 變化而有明顯之趨勢，但是在 Adnan 與 Sun[12]所模擬之 PE 與金屬 構成的膠合面中，直接調整介面間的凡得瓦能量參數，而不改變 PE 間的能量參數，使得靠近金屬接合面之PE 密度將隨著參數的放大而 有明顯成長。 討論完拉伸前初始結構之密度分布以及在空間中的排列狀態，接 著將進行膠合面結構的拉伸測試，本篇模擬參考文獻[12]藉著金屬板 與高分子間的凡得瓦勢能 εInterface，用以代表介面強度(黏著劑之效 力)，透過改變金屬板之位置，於 NVT 系綜下以時間步階 0.5 fs，首 先將每一萬步(step)內之第一步將上下兩層金屬板於 Z 方向拉伸至 0.1 %之應變，在之後的步數內進行平衡，使高分子能夠在這些步數內調 整位置並穩定，並在第一萬步紀錄 Z 方向以(2.3.3)式所計算之系統

Virial stress 以及平面間應力(interplanar stress)。金屬板之應變定義如 (3.2.1)式、圖 3.3.5，將未拉伸時最上層以及最下層金屬板之間的距離

視為材料 Z 方向之原長 L0，而藉著改變金屬板的位置進行拉伸，計

算Z 方向邊長 L 與原長 L0之差值，除以原長L0，即可得到Z 方向之

Virial Stress 並不適合處理非週期性的結構，並且對於結構內部的應力 較為模糊不清，另一方面 Machova[35]也討論了溫度在 0 K，原子尺 度下，當原子間作用力的範圍內有非均質(inhomogeneous)的應變狀 況，以(2.3.3)式所計算之系統 Virial stress，並不能恰當的描述其應力， 應該以平面間應力(interplanar stress)來計算則較適合，因此本篇模擬 的膠合面結構裡，即是將所有金屬原子對於每個在截斷半徑內之高分 子原子皆計算彼此間凡得瓦力的 Z 方向分量，如圖 3.3.7；此凡得瓦 力分量即可視為垂直穿過系統胞室截面積之作用力，將此作用力總和 除以截面積，可以得到金屬對於高分子產生之應力，如圖3.3.8 所示。 圖 3.3.9、圖 3.3.10 為此拉伸測試之過程中，上下層金屬板對於高分 子所造成的凡得瓦力與應變關係圖，其中圖 3.3.9 將凡得瓦力的吸 引、排斥力分別討論。由圖 3.3.9 可以發現，上下層所分別造成的吸 斥力有差異，但圖 3.3.10 中合力的部份在低應變(0 % ~3 %)的部份卻 相當一致，這個現象可能是由於上下層金屬各自截斷半徑內之高分子 原子數量不同，因此上下層金屬所造成的凡得瓦吸斥力不同；為了清 楚的觀察應變提高時，上下層金屬吸斥力隨著應變的提升所反應出的 變化，在此將上下層金屬吸斥力進行平移，使上層與下層之吸引力起 點相同、同樣的排斥力的起點也相同，如圖 3.3.11。由圖 3.3.10 可發 現，當開始測試時，上下金屬板在0 % ~ 3 %應變之所造成的凡得瓦

合力相同，隨著拉伸的進行，在3 %~8 %應變內凡得瓦力持續上升， 但上下層金屬板所造成的凡得瓦力已有差距，並且此差距逐漸增大， 此現象可以視為裂縫的產生，並且隨著應變增加逐漸成長，直到約8 %應變時，下層金屬板之凡得瓦力開始驟降，即是下層金屬板已經遠 離了高分子，破壞已經在下層發生；同樣的這個現象也可以從圖3.3.11 看到，上下層吸斥力在3 %應變後開始出現差距，在 8 %應變前上下 層金屬板吸斥力皆維持緩慢降低的狀態，在8 %應變後，下層吸引力 驟降，而應可視為下層已經在此時破壞；接著上層吸引力逐漸回復， 代表在3 %~8 %應變間，隨著金屬板在 Z 方向的移動，高分子核心材 料因受到上下層金屬板各自的凡得瓦力，而無法較貼近上下金屬板的 表面，當拉伸進行到 8%應變後，下層金屬板的吸引力已經無法與上 層競爭，所以高分子變得較貼近上層金屬板表面，上層金屬板之吸引 力以及排斥力皆上升。將上下層金屬板的凡得瓦力除以膠合面結構的 截面積，即可得到應力應變關係圖，如圖3.3.12。與圖 3.3.10 的趨勢 相同，當拉伸進行到 3 %應變時，上下層金屬所造成的應力開始分 歧，劃下第一條虛線，視即破壞的開始，接著拉伸到 8%應變時，下 層金屬所造成的應力達到高值(約 265 MPa)，隨即開始明顯下降，在 此劃下第二條虛線，代表破壞的完成，並且將此時破壞層的極限應力

面間應力值取平均，與系統的 Virial Stress 比較如圖 3.3.13，其中 Virial Stress 的有效體積會隨著拉伸測試進行時而改變，因為在 X、Y 方向 膠合面的邊長是固定的，所以當金屬板在Z 方向移動時，Z 方向的有 效長度取為最外層上下金屬板間的距離，如圖 3.3.14；可以看到在 8 % 應變前平均的平面間應力與Virial Stress 的值相當吻合，而 8 %應變 後則略有差異，顯示出Virial Stress 僅能考慮整個結構，所以與平均 後的平面間應力較為接近，但是平面間應力的計算則可以更進一步區 分出上下層破壞時，該層的拉伸強度，因此本篇將利用平面間應力所 算出的極限強度(ultimate strength) 265MPa(下層)，視為此膠合面結構

的破壞強度，最終的破壞情況如圖3.3.15，可以看出下層金屬板已經

與高分子脫離(破裂)。繼續將討論不同凡得瓦強度時，對於破壞過程

中所造成的影響，將原本具有正常凡得瓦介面強度(1εPI)之平衡結構，

藉著結合法則分別調整為 5εPI以及 10εPI，再於 NPT 系綜下時間步階

0.5fs，升溫至 300 K、0.0 MPa，進行 1000 ps 的擾動，接著降溫至 0K、 0.0MPa，持續 500 ps 的平衡，並且亦對 Z 方向 Virial Stress 微調，以

得到 5εPI以及 10εPI之初始結構，將此二結構亦進行和正常介面強度

狀態下完全相同的拉伸測試，以探討介面強度的強弱對於破壞過程的 影響。前文在正常的介面強度下，討論了破壞過程，包含破壞的起始， 從上下金屬板各自造成的應力有差異算起，以及破壞的結束─破壞的

該層平面間應力開始驟降時，所以兩個虛線間的應變量，代表正常膠 合面結構在破壞過程中約拉伸了5 %的應變，圖 3.3.16 為 5εPI時的破 壞過程，可以發現在破壞過程中應變量降低到約2.5 %應變，而由圖 3.3.17 中更可以看到 10εPI時的破壞過程僅拉伸有1.5 %應變；此一現 象與一般材料力學[36]中，延性材料承受拉伸下在破壞前能夠有較多 的伸長量，而脆性材料則不易伸長即破壞的現象相似，代表著當調整 凡得瓦力時，原子間的吸引力放大，使得整個結構趨向脆性，由圖 3.3.16、圖 3.3.17 中可以看到，經過放大凡得瓦力之膠合面結構，其 極限強度亦受影響而變大，如表11 所列，在正常凡得瓦強度 1εPI下， 極限強度約為265 MPa，而 5εPI時極限強度則約為 640 MPa，另外在 10εPI時則為980 MPa。 本篇模擬接合面的方法參考Adnan 以及 Sun[12]運用分子動力學 討論高分子黏著層的厚度對於接合面強度的影響，在奈米尺度時，金 屬層與「有效凡得瓦區域(Effective van der Waals Zone)」(即截斷半徑 內)是造成高分子黏著層較為緊密之原因，所以較薄之黏著層，其強 度較高，因高分子黏著層完全處於凡得瓦有效區域內，較為緊密；反 之如黏著層之厚度較大，則僅有部分落在凡得瓦有效區域內，使得接 合面的強度下降，這個結果雖然與一般廣泛認知的較大尺度之接合面

時，凡得瓦力所能影響到的範圍甚小，故不能直接將凡得瓦有效區域 當做大尺度時厚度與強度之關聯主因。另外他們改變金屬層與PE 黏 著層間的凡得瓦力，以代表介面間連結強度之強弱，當連結強度越好 時，破壞僅發生在高分子黏著層裡(cohesive failure)，而在正常或者較 弱的連結強度時，破壞則可能發生在接合處(adhesive failure)或者黏著 劑內。本篇研究則是將這個接合面的拉伸破壞行為，藉著平面間應力 所畫出的應力應變曲線，進一步分析出破壞發生於上層或下層，經由 分別建立多個模型，我們發現破壞並非固定的發生於某層，而是跟拉 伸前的初始結構有關，在 3.2 節中，先討論了 PI 本身在放大 εPI下的 機械性質的變化，接著在此小節討論以這些PI 為黏著劑的接合面(其 介面強度依然由結合公式推算)；所以本研究表現了當黏著劑本身機 械強度的不同，對於在拉伸測試下所產生的變化：較堅硬、脆性的高 分子結構(5εPI、10εPI)，有較短的「裂縫成長期」(圖 3.3.16、圖 3.3.17 虛線間的應變量)，而具有延展性的結構此成長期則較長，另外在多 次的模擬中，上下層都有發生破壞的紀錄，甚至當金屬板移動的距離 過大，也可能同時發生破壞，而藉著分析平面間應力，能了解到上下 層各自的應力應變狀態，可以知悉破壞層的極限強度，相較於使用 Virial stress 僅能計算整個系統的應力而言，應該是較為合適的。。這 個章節探討了乾燥的PI 與以 PI 為黏著層的黏著膠合面，分析了膠合

面Z 方向單軸拉伸時的破壞過程，以及調整凡得瓦力後，對於兩者在 機械性質上的影響；在下個章節將介紹關於水分子模型的建立以及水

第四章、水分子建構以及 PI 之水氣效應

4.1 水分子模型 TIP4P/2005

如何將水的物理性質，逼真的模擬出來，一直是分子動力學領域 中熱門的題目[37-39]。在生科領域中，也有學者探討DNA與RNA如 何在水中運動，並且形成雙股螺旋結構[40]，以下將針對如何選定水 分子模型做說明。SPC/E (Simple Point Charge/Extended) [37, 41]水分 子模型是1990年前後所提出，與能夠極化(水分子氫氧分解)的模型相 比，能有效降低計算時間而常被採用，此水分子為由三顆原子所組成 的共平面結構，兩顆氫原子各自帶有+0.4238 e (e為基本電荷，約為 1.6021892×10-19 Coulomb) 的 電 荷 ， 並 與 氧 原 子 構 成 之 夾 角 θ 為 109.47∘，與氧原子相距為1.0 Å，氧原子則帶有-0.8476 e之電荷。在 模擬純水結構中，僅考慮氧原子與其他氧原子間的凡得瓦勢能，參數 為0.1553 Kcal/mol、σ為3.473 Å(LJ-potential)；而TIP4P(Transferable Intermolecular Potential with 4 Points)[42]則是以四點來表現分子間勢 能，並且可以為了重現實驗結果而調整參數。2005年由Abascal以及 Vega[43]提出TIP4P/2005水分子勢能函數，並且在他們的研究中指 出，在描述密度、三態的變化溫度等物理性質上，TIP4P/2005是優於

SPC/E、TIP4P以及其他的水分子模型的。如圖4.1.1，此水分子為共 平面結構，兩顆氫原子各自帶有+0.5564 e (e為基本電荷，約為 1.6021892×10-19 Coulomb) 的 電 荷 ， 並 與 氧 原 子 構 成 之 夾 角 θ 為 104.52∘，分別與氧原子相距L1為0.9572 Å，一顆虛構、無質量的原 子被放置在θ之分角線上，與氧原子相距L2為0.1546 Å，並帶有-1.1128 e之電荷。與SPC/E模型相同的是，整個水分子結構中，僅有氧原子考 慮 與 其 他 分 子 間 的 凡 得 瓦 勢 能 ：0.1852 Kcal/mol、 σ 為 3.1589 Å (LJ-potential)，而靜電勢能(electrostatic potential)只由虛構原子與氫原 子來表現。值得注意的一點是，上述兩個水分子模型中，均沒有獨立 描述氫鍵的項，而是藉著適當的凡得瓦勢能與靜電勢能，表現出符合 現實中水的物理性質，並且整體水分子視為一不會分解的剛體(rigid body)。Docherty[44]等人在分子動力學領域著眼於甲烷(methane) 在 水中所可能產生之水化物，藉著考慮Lorentz–Berthelot combining rule 來表達甲烷與水分子間能量(the methane-water LJ interaction energy ε)，使用TIP4P/2005之水分子模型，僅需將ε增大7%就能忠實的呈現 甲烷在水中之過量化學能(excess chemical energy)；而Goujon[45]等人 亦藉著TIP4P/2005模擬了表面改質之石墨板(graphite)接枝許多具負 電荷的phenyl-NTA分子，探討水分子與這個改質過的石墨片間的作

水分子使得周圍結構更緊密。基於以上的討論，顯示此水分子模型的 在重現水的性質有不錯的表現，屬於目前較為準確的模型，因此本篇 模擬將採用此模型，討論水氣對於高分子本身以及膠合面在機械性質 上所產生的影響。

4.2

PI 之水氣效應

本節將藉由分子動力學，探討 2 wt%的水氣對於 PI 本身結構所 造成之影響。水分子與 PI 原子間，考慮了凡得瓦力(水中的氧原子對 於其他的PI 原子們)與靜電力的作用；首先分別建立兩個完全相當的 PI 模型初始模型氣在其一任意置入了 2 wt%的水分子，如圖 4.2.1 所 示，TIP4P/2005 模型將氫鍵以凡得瓦力、靜電力來模擬，如果像 3.2 節中的乾燥 PI 不考慮靜電力，則可能無法完全的模擬出水的特性， 因此本節中所建立的乾燥、含水氣的兩個 PI 模型，均考慮靜電力， 而所選用的靜電力勢能以及PI 原子上的部份電荷如 3.1 節中所說明。 再進行平衡的流程，第一步驟於800 K 的環境下，採用 NPT 系綜， 壓力設定為200 MPa、持續 500 ps (picosecond, ps)，藉著壓力將結構 壓縮至緊密狀態，接著第二步驟同樣在800 K，採用 NPT 系綜，壓力 設定為0.1 MPa (一大氣壓)、持續 500 ps，讓結構自由鬆弛，而後第 三步驟則在NPT、10 MPa 狀態下，降溫到 600 K，持續 500 ps，第四

已經分別得到室溫、一大氣壓下，乾燥與含水 PI 這兩個結構，其中 含水 PI 的結構如圖 4.2.2 所示，水分子不規則的分布在 PI 中。其密 度的比較如表12 所列，在本篇模擬中，含水 PI 之密度較乾燥的上升 了約1.6 %，與多於文獻中密度提昇 0.7 %，但是趨勢相同，而乾燥 的PI 密度為 1.23 g/cm3亦在合理的密度範圍內。接著使用 3.2 節中的 方法二分別對在0 K、100 K、以及 300 K 之溫度下平衡(stress-free) 的乾濕 PI 結構，進行楊氏係數的計算。而因為此節中，具有 100 K、 300 K 的結構，能量與應力在這些溫度下會有小幅度之振盪，所以當 得到應力平衡結構後，對於在每個溫度下分別給予材料X 方向 1%正 向應變(ε11 =0.01、ε22 = ε33 = 0)，並於 NVT 系綜下平衡；由於在平衡 過程中系統應力會收斂但不會呈現定值，所以要決定系統之應力時， 先對結構進行持續長達1500 ps 之平衡，並採用 1000 ps 至 1500 ps 模 擬時間內之應力總和平均視為系統應力，當作穩定的應力值進行 C 矩陣之計算；進行總和平均的目的是為了避免在變化的應力下，使用 瞬間的應力值計算所造成的誤差。將這些施加 X 方向應變後所產生 之三方向應力σ1、σ2、與σ3，將這些值帶入(3.2.2)式中可求得勁度矩 陣中之C11 、C21、與 C31。再分別依照此流程施加應變於Y、Z 方向， 即可得到完整的勁度矩陣，運算反矩陣 S，最終求得 E1、E2以及E3；

細數值紀錄於表13、表 14。藉著在 0 K、100 K、300 K 所記錄之楊 氏係數，與文獻[5]在不同溫度下所得到之楊氏係數比較整理如表 15 所示。由圖 4.2.3、圖 4.2.4 以及表 15 中可發現，在本篇模擬中，不 論在哪個溫度下，置入水氣皆使得 PI 的楊氏係數上升，這一現象與 普遍將水氣視為塑化劑，吸收將使結構軟化、強度下降之觀念相悖； 本篇模擬之水氣效應在0K 時甚為明顯，而在 100 K 時略微減弱，而 300 K 時影響則相當小，整體趨勢與文獻相符；但較為細部的現象， 仍尚有差異。如文獻中在室溫(300K)時，水氣效應相當微弱，對於楊 氏係數之增幅僅有3 % ; 當溫度低於氫鍵斷裂的溫度(約 230 K)，也 就是水氣效應發揮作用的關鍵點時，含水氣 PI 的楊氏係數會大幅度 的提升，在100 K 時對於含水氣材料之楊氏係數提高約 22 %；而本 篇模擬在300K 時即提高了約 14 %，於 100K 時則約提高 15 %，在 0 K 時則提高 23 %，未能明顯表現出此影響水氣效應的關鍵溫度。另 外一方面在乾燥PI 本身的楊氏係數，隨著溫度上升有大幅度的下降， 顯示出目前的模擬勢能在表現 PI 之楊氏係數方面尚有不足，依照化 學結構的不同，PI 的楊氏係數約落在 2.0 ~ 3.52 GPa 間[46-48]，本模 擬室溫時所計算出之楊氏係數大約符合此低值，但這也是未來在模擬 不同的高分子材料時，期望能改進的目標：尋求更理想的勢能函數， 進 而 運 用 MD 計 算 常 溫 的 高 分 子 機 械 性 質 。 本 章 節 選 擇 了 以

TIP4P/2005 水分子模型，此剛體水模型將水的氫鍵以凡得瓦勢能與靜 電力來考慮，因此在進行與此水模型相關的研究時，皆不能略去電荷 的貢獻，所以將乾燥 PI 選用 CVFF 之部分電荷，使之能與水分子藉 著靜電勢能交互作用。由於各方向皆考慮週期性邊界，所以此方法無 法進行如水的擴散性(Diffusivity)、蒸發量等模擬；在不同溫度下，2 wt%水氣對於 PI 在楊氏係數上都造成強化效果：低溫時(0 K、100 K) 有顯著的效果，反之在室溫(300 K)則較無差異，這個現象支持前述凡 得瓦力的論述[25]─將帶有較強的凡得瓦力的水分子置入 PI 結構 中，使楊氏係數上升，並且顯示出，在0 K 時，沒有溫度的影響，靜 電力的考慮即造成乾燥 PI 楊氏係數的些微上升；因為靜電力與凡得 瓦力皆為距離函數，在 Virial stress 的計算下，較高的溫度時各原子 的熱振盪距離也較大，使得這個強化的現象變得不明顯，恰好可以表 現出水氣效應的溫度依賴性。由於要探討膠合面吸附水氣後對於破壞 強度產生之影響時，尚有一些困難點存在，包括如何決定金屬板之電 荷以及選用適當的靜電勢能等，因此將這個部分放在未來工作。

第五章、結論

本研究運用分子動力模擬，藉著施加單軸拉伸測試，調查了具有 不同原子間作用力時，PI 聚合物之機械性質與以 PI 為膠黏劑之膠黏 接合面之間的介面間破壞現象。另外並討論水氣對於 PI 高分子本身 的影響，分別在不同溫度下，比較其楊氏係數的變化。所得的結果可 以歸納下列幾點： 1. PI 本身密度與機械性質(包括楊氏係數以及降伏應力)隨著凡 得瓦力的放大而上升，整體結構變得緊密堅硬。 2. PI 膠黏劑/金屬接合面在凡得瓦力的放大下，膠黏劑本身脆 化，使得拉伸的極限強度提升，但破壞過程應變量減少。 3. 含水氣的 PI 在低溫時具有明顯的反塑化作用，楊氏係數大幅 提升，而隨著溫度上升此強化現象減弱。

第六章、未來工作

在本文第三章的部份，完成了乾燥的PI 與以 PI 為膠黏劑材料的 膠合面分析，其中由兩種計算楊氏係數的方法中，顯示了結構在等向 性的表現上尚有改進的空間，未來在討論高分子材料時，將可利用更 大的模型、更恰當的平衡過程以建構較為均質的模型。接著於第四章 中選定了TIP4P/2005 水分子模型，討論加入了 PI 中，對其密度與楊 氏係數隨著溫度所產生的變化；未來將延續此水分子模型的運用，將 之加入於膠合面內，模擬攝入水氣時，對於膠合面介面破壞的影響。 而這個模擬有一些需要考慮的部分，首先是在乾燥的情況下，膠合面 的拉伸僅需考慮金屬板對於 PI 間的凡得瓦力，但是在置入水氣後， 為了重現氫鍵的特性，靜電力必須列入考量，因此如何決定金屬板上 所帶的電荷以及搭配適合的靜電勢能，是需要解決的問題。而水氣的 分布情況，對於此拉伸行為也會產生影響，如果水氣分散於 PI 中， 則在拉伸金屬時，較無直接的作用力，(PI 同時受到金屬層與置於其 中的水氣影響)；若水氣集中吸附於介面，則可能對單邊的破壞降伏 強度有明顯的影響。本文選擇了 PI 此種高分子，其中具有的官能基

當這類疏水性材料受到浸泡而使水氣進入結構，這時水氣則可能發揮 塑化劑的作用，使得材料變軟。所以未來在不同種類的高分子材料 上，水氣效應對於不同(數量)的親水基、以及親水基的有無等情況， 所表現出來的物理、機械性質的變化，亦是本模擬的主要目標，並且 與實驗數據比對。而最終目標，也是水氣最耐人尋味的性質，隨著水 氣含量的逐步增加，對於機械強度、密度影響的效果會經過一折返點 而 呈 現 相 反 的 趨 勢 ─ 塑 化(plasticization) [49] 與 反 塑 化 作 用 (antiplasticization)。Chamarthy 等人[50]以微晶纖維素(microcrystalline cellulose ,MCC )為研究對象，指出不論是製造過程中或者儲藏時讓 MCC 暴露於潮濕環境下，將會影響 MCC 之機械性質。而一般視為 塑化劑之水氣，在低含量時反而會有反塑化劑(antiplasticizers)之效 果。以微觀角度來探討反塑化作用，MCC 之滲透性(permeability, Dp) 起先隨著水氣含量的上升而下降，直到5% RH 後，才隨著 RH 的 提高而上升；獨立的將MCC 對高度疏水性(hydrophobic)之辛烷(octane) 也進行滲透性與RH 之關係圖，亦得到和水氣的滲透性完全相同的趨 勢。水氣與辛烷之 Dp減少，意味著低含量的水氣，增加了高分子的 堅硬度。當水氣剛進入乾燥的高分子時，水分子本身(氫鍵的搭建)阻 礙了自己以及辛烷的移動；直到 RH 高於 5%，才發揮有助於滲透的 效果。從巨觀角度來探討反塑化作用，藉由三點彎曲方法(The triple

point beam bending method )對材料進行在不同 RH 下，楊氏係數以及 拉伸強度的測試，亦可發現在低RH 時，水氣對於楊氏係數之強化效 果顯著，此效果在約30% RH 時開始驟降。值得注意的是，在討論不 同性質時，其反塑化與塑化作用的轉換臨界點(turnabout point)是不同 的，作者解釋為：如討論微觀角度、分子的行為，則轉換的門檻較低 (在較低 RH 即發生轉換)；若討論巨觀層面、機械行為，反應整個高 分子的內聚能(cohesive energy)之性質，則此轉換門檻則較高。 另 外 DLUBEK 等 人 [51] 利 用 正 電 子 湮 滅 時 間 光 譜 (positron

annihilation lifetime spectroscopy, PALS)[52]量測尼龍 6(nylon 6)之區 域自由體積(local free volume)隨著 RH 上升的變化，大致可以將實驗

過程分為兩個階段。第一階段為RH 小於 50%時，自由體積較未含水 時的略為降低，可以視為已經存在的自由體積被水分子填滿所致；而 第二階段為RH 大於 50%時，自由體積又恢復並且隨著 RH 增加大量 超越了未含水氣時，這個現象則為隨著過多水氣的攝入，弱化了原本 已經存在之水分子與胺基(amide group)間的氫鍵，並導致塑化作用的 發生。 前述兩篇文獻討論了水氣的塑化與反塑化作用，在不同性質上的 影響，未來本研究將期望能模擬出多種高分子材料(不僅止於 PI)，建

參考文獻

[1] J. Vinson, "On the Mechanical Behavior of Advanced Composite Material Structures," Major Accomplishments in Composite Materials and Sandwich Structures, pp. 431-439, 2009.

[2] K. Kolat, G. Neser, and C. Ozes, "The effect of sea water exposure on the interfacial fracture of some sandwich systems in marine use," Composite Structures, vol. 78, pp. 11-17, 2007.

[3] R. Selzer and K. Friedrich, "Mechanical properties and failure behaviour of carbon fibre-reinforced polymer composites under the influence of moisture," Composites Part A, vol. 28, pp. 595-604, 1997.

[4] A. Athijayamani, M. Thiruchitrambalam, U. Natarajan, and B. Pazhanivel, "Effect of moisture absorption on the mechanical properties of randomly oriented natural fibers/polyester hybrid composite," Materials Science and Engineering: A, vol. 517, pp. 344-353, 2009.

[5] G. A. Bernier and D. E. Kline, "Dynamic mechanical behavior of a polyimide," Journal of Applied Polymer Science, vol. 12, pp. 593-604, 1968.

[6] H. K. Glen and A. L. Jennifer, "Rheological Property Evolution in Concentrated Cement-Polyelectrolyte Suspensions," Journal of the American Ceramic Society, vol. 85, pp. 2989-2994, 2002.

[7] J. K. Sears and J. R. Darby, The technology of plasticizers. New York: John Wiley & Sons, 1982.

[8] W. Smith, C. W. Yong, and P. M. Rodger, "DL_POLY: application to molecular simulation," Molecular Simulation, vol. 28, pp. 385-471, 2002.

[9] J. H. Irving and J. G. Kirkwood, "The Statistical Mechanical Theory of Transport Processes. IV. The Equations of Hydrodynamics," The Journal of Chemical Physics, vol. 18, pp. 817-829, 1950.

Dinola, and J. R. Haak, "Molecular dynamics with coupling to an external bath," The Journal of Chemical Physics, vol. 81, pp. 3684-3690, 1984.

[12] A. Adnan and C. T. Sun, "Effect of Adhesive Thickness on Joint Strength: A Molecular Dynamics Perspective," The Journal of Adhesion, vol. 84, pp. 401-420, 2008.

[13] S. L. Mayo, B. D. Olafson, and W. A. Goddard, "DREIDING: A generic force field for molecular simulations," The Journal of Physical Chemistry, vol. 94, pp. 8897-8909, 1990.

[14] J. E. Lennard-Jones, "Cohesion," Proceedings of the Physical Society, vol. 43, pp. 461-482, 1931.

[15] M. P. Allen and D. J. Tildesley, Computer simulation of liquids. Oxford: Oxford university press, 1989.

[16] J. Wang, R. M. Wolf, J. W. Caldwell, P. A. Kollman, and D. A. Case, "Development and testing of a general amber force field," Journal of computational chemistry, vol. 25, pp. 1157-1174, 2004. [17] Y. Cheng, Q. Meng, W. Wang, and J. Wang, "Structure–function

relationship of three neurotoxins from the venom of Naja kaouthia: a comparison between the NMR-derived structure of NT2 with its homologues, NT1 and NT3," Biochimica et Biophysica Acta (BBA)/Protein Structure and Molecular Enzymology, vol. 1594, pp. 353-363, 2002.

[18] J. K. Seydel and M. Wiese, Drug-membrane interactions: analysis, drug distribution, modeling. weinheim: Vch Verlagsgesellschaft Mbh, 2002.

[19] T. C. Clancy, S. J. V. Frankland, J. A. Hinkley, and T. S. Gates, "Molecular modeling for calculation of mechanical properties of epoxies with moisture ingress," Polymer, vol. 50, pp. 2736-2742, 2009.

[20] U. Dinur and A. T. Hagler, "New approaches to empirical force fields," Reviews in computational chemistry II, p. 99, 1991.

[21] J. D. Jackson, Classical Electrodynamics, 3rd edition. New York: John Willey & Sons, 1999.

[22] G. Rabilloud, High-performance polymers: chemistry and applications. Paris: Editions Technip, 2000.

[23] G. M. Odegard, T. C. Clancy, and T. S. Gates, "Modeling of the mechanical properties of nanoparticle/polymer composites," Polymer, vol. 46, pp. 553-562, 2005.

[24] T. C. Clancy, "Multi-scale modeling of polyimides," Polymer, vol. 45, pp. 7001-7010, 2004.

[25] C. Wei, "Adhesion and reinforcement in carbon nanotube polymer composite," Applied Physics Letters, vol. 88, p. 093108, 2006. [26] D. M. Delozier, R. A. Orwoll, J. F. Cahoon, J. S. Ladislaw, J. G.

Smith, and J. W. Connell, "Polyimide nanocomposites prepared from high-temperature, reduced charge organoclays," Polymer, vol. 44, pp. 2231-2241, 2003.

[27] I. M. Daniel and O. Ishai, Engineering mechanics of composite materials. New York: Oxford University Press, USA, 1994.

[28] V. J. Vasudevan and J. E. McGrath, "Atomistic modeling of amorphous aromatic polybenzoxazoles," Macromolecules, vol. 29, pp. 637-645, 1996.

[29] A. Adnan, C. T. Sun, and H. Mahfuz, "A molecular dynamics simulation study to investigate the effect of filler size on elastic properties of polymer nanocomposites," Composites Science and Technology, vol. 67, pp. 348-356, 2007.

[30] I. M. Ward and J. Sweeney, An introduction to the mechanical properties of solid polymers. Chiehester, England: John Wiley & Sons Inc, 2004.

[31] B. Hartmann, "Tensile Yield in Crystalline Polymers," Marcel Dekker, Inc., Handbook of Polymer Science and Technology., vol. 2, pp. 101-125, 1989.

[32] H. Mahfuz, A. Adnan, V. K. Rangari, M. M. Hasan, S. Jeelani, W. J. Wright, and S. J. DeTeresa, "Enhancement of strength and stiffness of Nylon 6 filaments through carbon nanotubes reinforcement," Applied Physics Letters, vol. 88, p. 083119, 2006.

[33] M. Rubinstein and R. H. Colby, Polymer physics. New York: Oxford University Press, USA, 2003.

[34] D. H. Tsai, "The virial theorem and stress calculation in molecular dynamics," The Journal of Chemical Physics, vol. 70, pp. 1375-1382, 1979.

[35] A. Machova, "Stress calculations on the atomistic level," Modelling and simulation in materials science and engineering, vol. 9, pp. 327-337, 2001.

term in effective pair potentials," The Journal of Physical Chemistry, vol. 91, pp. 6269-6271, 1987.

[38] K. Takemura and A. Kitao, "Effects of water model and simulation box size on protein diffusional motions," The Journal of Physical Chemistry B, vol. 111, pp. 11870-11872, 2007.

[39] I. C. Lin, A. P. Seitsonen, M. D. Coutinho-Neto, I. Tavernelli, and U. Rothlisberger, "Importance of van der Waals Interactions in Liquid Water," Journal of physical chemistry. B, 2009.

[40] G. L. Seibel, U. C. Singh, and P. A. Kollman, "A molecular dynamics simulation of double-helical B-DNA including counterions and water," Proceedings of the National Academy of Sciences, vol. 82, p. 6537, 1985.

[41] I. M. Svishchev and P. G. Kusalik, "Structure in liquid water: A study of spatial distribution functions," The Journal of Chemical Physics, vol. 99, pp. 3049-3058, 1993.

[42] W. L. Jorgensen and J. D. Madura, "Temperature and size dependence for Monte Carlo simulations of TIP4P water," Molecular Physics, vol. 56, pp. 1381-1392, 1985.

[43] J. L. F. Abascal and C. Vega, "A general purpose model for the condensed phases of water: TIP4P/2005," The Journal of Chemical Physics, vol. 123, p. 234505, 2005.

[44] H. Docherty, A. Galindo, C. Vega, and E. Sanz, "A potential model for methane in water describing correctly the solubility of the gas and the properties of the methane hydrate," The Journal of Chemical Physics, vol. 125, p. 074510, 2006.

[45] F. Goujon, C. Bonal, B. Limoges, and P. Malfreyt, "Molecular Dynamics Description of Grafted Monolayers: Effect of the Surface Coverage," The Journal of Physical Chemistry B, vol. 112, pp. 14221-14229, 2008.

[46] W.-J. Chou, C.-C. Wang, and C.-Y. Chen, "Characteristics of polyimide-based nanocomposites containing plasma-modified multi-walled carbon nanotubes," Composites Science and Technology, vol. 68, pp. 2208-2213, 2008.

[47] D. M. Delozier, R. A. Orwoll, J. F. Cahoon, J. S. Ladislaw, and J. G. Smith, "Polyimide nanocomposites prepared from high-temperature, reduced charge organoclays," Polymer, vol. 44, pp. 2231-2241, 2003.

"Preparation and morphological, electrical, and mechanical properties of polyimide-grafted MWCNT/polyimide composite," Journal of Polymer Science Part A: Polymer Chemistry, vol. 45, pp. 3349-3358, 2007.

[49] A. Bos, I. G. M. Punt, M. Wessling, and H. Strathmann, "CO2-induced plasticization phenomena in glassy polymers," Journal of Membrane Science, vol. 155, pp. 67-78, 1999.

[50] T. Carvajal, S. P. Chamarthy, and R. Pinal, "Moisture Induced Antiplasticization of Microcrystalline Cellulose," 2006, pp. 236-236.

[51] G. Dlubek, F. Redmann, and R. Krause-Rehberg,

"Humidity-induced plasticization and antiplasticization of polyamide 6: A positron lifetime study of the local free volume," Journal of Applied Polymer Science, vol. 84, pp. 244-255, 2002. [52] D. Cangialosi, H. Schut, A. van Veen, and S. J. Picken, "Positron

annihilation lifetime spectroscopy for measuring free volume during physical aging of polycarbonate," Macromolecules, vol. 36, pp. 142-147, 2003.

附表

表 1 分子動力學之系綜[11] 系綜 控制條件 系綜描述 NVE 原子數目N 系統體積V 系統能量E 原子數目、系統體積為定值，遵守能量 守衡，使模擬進行中之能量維持定值。 NVT 原子數目N 系統體積V 系統溫度T 原子數目、系統體積為定值，進行溫度 調節，使模擬進行中之溫度接近初始設 定值。 NPT 原子數目N 系統壓力 P 系統能量T 原子數目為定值，進行溫度與壓力調 節，藉由調整系統體積與原子位置，使 模擬進行中之溫度，三方向壓力之平均 接近初始設定值。 NST 原子數目N 系統各方向應力S 系統溫度T 原子數目為定值，進行溫度與壓力調 節，藉由調整系統體積與原子位置，使 模擬之溫度，各方向應力接近初始設定 值。

表 2 原子代號與質量電荷表[13, 20] 原子種類 原子代號與狀態描述 原子質量(g/mol) 電荷(e) C_2 :與氧以雙鍵連結 0.380 碳(Carbon) C_R :構成苯環之碳 12.011 0.000 HW :水分子中的氫 1.008 0.5564 氫(hydrogen) H_ :接於碳上之氫 1.008 0.100 OW:水分子中的氧 15.9994 0.000 O_2 :與碳以雙鍵連結 -0.380 氧(oxygen) O_3 :與碳以單鍵連結 15.9994 -0.056 氮(Nitrogen) N_3 :與碳以單鍵連結 14.0067 -0.110 金屬 Al3 :為 FCC 排列金屬 26.9815 0.000 表3 延展鍵結勢能參數表[13] Kr (Kcal/mol- Å2) r0(Å) N_3─C_R 700 1.362 N_3─C_2 700 1.392 C_R─C_2 1050 1.360 C_R─C_R 1050 1.390 C_R─O_3 700 1.350 C_2─O_2 1400 1.220

表4 角度勢能參數表[13] Kθ (Kcal/mol) θ0 (drgree) XX─N_3─XX 109.000862 106.70 XX─ C_R─XX XX─ C_2─XX 133.333333 120.00 XX─O_2─XX XX─O_3─XX 106.697940 104.51 表5 扭轉勢能參數表[13]

Aτ(Kcal/mol-rad2) δ (degree) m (constant) XX─N_3─C_R─XX

XX─N_3─C_2─XX 0.125 -180 6

XX─C_R─C_R─XX

XX─C_R─C_2─XX 1.25 180 2