膠體粒子對平面之電泳運動(3/3)

(1)

行政院國家科學委員會專題研究計畫成果報告

膠體粒子對平面之電泳運動(3/3)

計畫類別：個別型計畫計畫編號： NSC93-2214-E-002-003- 執行期間： 93 年 08 月 01 日至 94 年 07 月 31 日執行單位：國立臺灣大學化學工程學系暨研究所計畫主持人：李克強計畫參與人員：池明輝、唐于博、羅仕瀚報告類型：完整報告處理方式：本計畫可公開查詢

中華民國 94 年 10 月 21 日

(2)

摘要本論文探討一膠體粒子垂直於一不帶電平板的電泳運動現象。在此我們解除過去求解電泳現象的諸多限制，首先於第一年裡求解一任意表面電位下的長直圓柱運動；並且在第二年更進一步去探討任意外加電場下的運動情形；最後在最後一年突破以往求解雙球座標單一區間的情形，將原本的計算空間延伸至粒子內部，進一步利用正交配位法及牛頓-拉福生疊代法求解系統的非線性電場及流場方程式。研究結果發現，當電雙層厚度較薄的時候，其電泳速度相對來說比較快，但是當粒子的表面電位越大，極化效應就越顯得重要。可是當外加電場強度逐漸變大了以後，其驅動力將超過粒子表面電位造成的極化效應，使得粒子的運動速度又恢復至正比於電雙層厚度的結果。另外當液滴內外黏度比愈小，液滴內部的拖曳力越小，其電泳速度隨之增加，反之當內外黏度比越大的時候，其電泳速度越接近硬球膠體粒子的結果。此外，發現平板的邊界效應在電雙層厚度很大或是粒子離平板越靠近的時候相當明顯，尤其粒子太接近平板的時候會因為電雙層的變形而使得電泳速度異號，亦即電泳運動受邊界影響甚為重要。關鍵字：電泳現象、雙球座標、電雙層、正交配位法。

(3)

Abstract

The electrophoretic behavior of the spherical, non-conducting particle normal to a plane is investigated in this study. We separate this system into three topics. In the first year, we consider the electrophoresis of a long cylindrical particle with an arbitrary zeta potential normal to the plane. Moreover, we find the electrophoretic behavior under the arbitrary applied field in the second year. Finally, with the aid of orthogonal collocation method and bipolar coordinates, we extend previous analyses in that the flow field inside the spherical particle can not be neglected, and the electric double layer thickness for the particle is arbitrary.

We find that the electrophoretic mobility of the particle is affected by some factors: double layer thickness, zeta potential, the strength of applied field, the viscosity ratio between the fluid inside and outside the liquid drop, and the distance to the solid plane. We find that the magnitude of the scaled electrophoretic mobility increases with the decrease of double layer thickness. However, the polarization effect must be considered when the magnitude of zeta potential is large. Besides, the mobility will increase with the decrease of viscosity ratio. This is because the flow inside the drop enhances the hydrodynamic drag on the liquid drop. In addition, we find that the closer the particle to the plane, the more significant the distortion of electric double layer. The electrophoretic mobility becomes slower.

(4)

目錄中文摘要 I 英文摘要 II 目錄 III 圖表目錄 V 第一章緒論 1 1.1 膠體系統介紹 1 1.2 文獻回顧 5 第二章理論分析 11 2-1 系統描述 11 2-2 主控方程式基本介紹 14 2-3 平衡狀態 17 2-4 擾動狀態 1 9 2-5 電泳速度之計算 3 2 第三章數值方法 3 5 3-1 正交配位法 3 6 3-2 空間映射 40 3-3 兩區聯解問題之處理 44 3-4 牛頓-拉福生疊代法 46 3-5 數值積分 49 第四章圓柱對平板在任意表面電位下之結果與討論 50 第五章球對平板在任意外加電場下之結果與討論 71 第六章液滴對平板結果與討論 78

(5)

第七章結論與成果自評 103 第八章參考文獻 105 附錄 A 座標系統簡介 11 0 附錄 B 流場方程式相關推導 11 5 附錄 C 力積分之詳細推導 1 2 2

(6)

第

一章

序論

1.1 膠體系統介紹

膠體系統是自然界中相當常見的系統，例如：化妝品(cosmetics)、染料 (dyestuffs)、藥物(pharmaceuticals)、乳化物(emulsions)、紡織品(fabrics)…等等，舉凡日常生活中的食衣住行，都是膠體系統應用的範疇，所以膠體科學是一門實用性相當廣泛的科學，值得好好地深入探究。一般膠體系統容易與真溶液(true solution)系統造成混淆，但簡單地說，膠體溶液(colloidal dispersion)與真溶液的基本差異在於真溶液之溶質係“溶解”在溶劑中，形成一均相(homogeneous phase)，稱為溶液；而膠體系統為兩個“不互溶(immiscible)”的相，為非均相(heterogeneous phase)，因此可區分為分散相 (disperse phase)，即是形成膠體粒子之相，與分散介質(dispersion medium)，即是分散粒子所分布的介質。此外，膠體系統的分散相的粒子尺寸遠大於分散介質的分子，這就是膠體系統與真溶液的主要不同，也就造成這兩個系統在物理及化學特性上的差異，所以描述它們的理論及在生活上的應用也就不同。一般而言，膠體溶液分散相的粒徑下限約在1nm左右，若小於此值，則將難以與真溶液加以區分；粒徑上限約為1 mµ ，但沒有一確切的範圍，因為在許多方面，譬如：乳化物、礦物分離過程，所遭遇到的分散粒子皆大於1 mµ 。所以，膠體溶液的分散相範圍約在nm至µm之間，至於分散介質的大小，如同一般分子，是小於分散相的。

(7)

吾人常用親媒性( lyophilic，solvent loving)與疏媒性( lyophobic，solvent hating)來描述膠體界面或其上之官能基與溶液分子之間溼潤(wetting)或媒合 (solvated)的趨勢。若分散介質為水，則使用親水性(hydrophilic) 與疏水性 (hydrophobic)描述之。簡單地說，對於親媒性的膠體系統而言，分散粒子之間傾向“散布”於分散介質中，當這種狀態存在時，因為 Gibbs 自由能(free energy) 的降低而是熱力學穩定(thermodynamically stable)。對於疏媒性的膠體系統而言，若上述的情況出現，則會因為 Gibbs 自由能的增加而不穩定 (thermodynamically unstable)，因此維持分散粒子聚集成團，可以使 Gibbs 自由能產生最小值，呈現穩定狀態。如果想要使疏媒性膠體溶液中的分散粒子均勻分散，則必須對膠體粒子的表面作一些特殊處理，使得分散粒子間存在強的排斥力，達到分散的目的。例如加入電解質使得分散粒子表面代電而相互排斥，或加入高分子物質，使之附著於分散粒子上，以增強分散粒子間的斥力，可避面聚集成團的結果，但需要強調的是，這樣的膠體系統仍然是不穩定的，在長時間放置下，分散粒子仍然會聚集成塊。膠體系統中，表面性質往往決定了這個系統的穩定性以及實用價值，所以界面科學(interface science)實為膠體科學中的重要關鍵。又所謂的“界面”，即是指分散相與分散介質的交界面，其特性就代表了膠體的表面特性，譬如：電荷吸附(adsorption)與電雙層(electric double layer)…等等，這些特性對於整體溶液的物理性質有重大的影響，例如在使用油井鑽探的技術時，某些黏土懸浮液之濃度會受到少量鈣離子或磷離子的加入而改變的界面特性，因而有增濃 (thickening)與稀薄(thinning)等不同作用，所產生巨大的變化，就是膠體界面性質的改變所造成。由此可知界面科學在膠體系統中是不可忽略的一環，必須要掌握膠體的界面特性，才可以妥當的應用膠體系統，發揮其最大功效，提升系

(8)

大部分的物質接觸到一極性介質時，會經由離子化(ionization)、離子吸附及離子解離(ion adsorption and dissolution)等反應機制而產生表面電荷，而這些表面電荷會影響到鄰近溶液中的電荷離子分佈情形，使得在粒子表面與整體溶液間會形成一帶電區域，稱為電雙層(electric double layer)。當粒子表面電雙層受到外力作用時，電雙層對於粒子運動行為產生電力及流力的效應，稱此現象為膠體的電動力學現象(Electrokinetic phenomena)。一般而言，膠體電動力學現象可以分為以下幾類： (1)電泳(electrophoresis)：帶荷電的膠體表面連同其表面附著物，在外部施加電場 的影響下，相對於固定不動液體之移動行為。 (2)電滲透(eletro-osmosis)：液體在外加電場下，相對於不動荷電粒子表面之運動 現象。 (3)沈降電位(sedimentation potential)：帶電粒子相對於不動液體運動時所產生的 電場。粒子的運動可以是重力場或離心力場所造成。 (4)流場電位(streaming potential)：液體沿著不動之帶電粒子運動時，所產生的電 場。液體的運動可以是壓力差所造成。 (5)電導(electric conductivity)：帶電粒子在外加電場作用下所產生電流以導電的現 象。上面的電動現象中，以電泳應用範圍最廣，實用性也較高，並可由此種技術得知膠體表面的帶電荷性質及發展其相關應用。又一般而言，電泳運動具有相當廣泛的用途，譬如生化實驗中的電泳，常用來模擬並分析DNA 或蛋白質等物質的表面現象與帶電情形[1]；又譬如本計畫中的粒子對平板的運動現象，目前已知可以廣為應用於廢水處理、電泳塗佈[2]等相關領域。因此由前述可知，膠體系統的穩定性常與膠體粒子的表面電荷性質息息相

(9)

關，又膠體系統的穩定性可以決定此系統之價值，所以研究及量測膠體之電荷性質是一個重要的課題，其中電泳運動可以提供很多關於膠體統有利用價值的資訊，並且具有可觀之工業應用價值。

(10)

1.2 文獻回顧

膠體粒子電泳長久以來被廣泛地應用在不同的領域，而電泳理論最早是由 Smoluchowski[3]於 1918 年發表膠體電泳速度與表面電位之物理性質關係，唯其假設流體為牛頓流體，且電雙層厚度(係與離子強度與電荷價數相關)相對於粒子半徑極薄，並且不考慮邊界效應(單一膠體粒子在無窮大流場中)的影響下所得到其電泳速度為： E U a E _πµ εζ 4 = (1.1) 而接著Hückel[4]在 1924 年時導出電泳速度在電雙層很厚的時候但在低表面電位及弱外加電場。之後Henry[5]則進一歩延伸了上述電泳理論的適用範圍，他推導出單一膠體粒子在任意電雙層厚度下的電泳速度。其後的研究者多以突破局部之限制，例如假設低表面電位，但電雙層厚度相對於粒徑是極薄或極厚[6,7]的情況，得到簡化系統下的解析解。

最早使用電腦來計算電泳的是Wiersema[8]。至 1978 年，O’Brien 及 White[9] 採用打靶法(shooting method)求得球形膠體粒子於無邊界效應影響下之電泳速度，滿足任意電雙層厚度、任意電位，並考慮電雙層之扭曲變形效應，此一方法必須滿足當外加電場相對於粒子本身帶電所造成的電場很小時，可將外加電場對電位的影響所造成的方程式、兩個濃度方程式和流場方程式線性化，利用此特性，把電泳泳動度的計算簡化，將電泳現象分為兩部分: (1)粒子不動，但受外加電場作用。 (2)粒子以一定速度運動但無外加電場。 接著陸續有人對於單一膠體粒子在稀薄溶液的電泳運動之後[10-11]，也有人對於其它系統中的電泳運動提出其相關的研究。

(11)

在其它的電泳系統中較常見為一群膠體粒子在懸浮溶液中或是膠體粒子在有邊界的範圍內所進行的電泳運動。在探討密集膠體粒子的電泳現象方面，首見是 Levine 和 Neale[12]在 1973 年採用了 Kuwabara[13]的單位細胞模型(Unit cell model)，此一模型可以有效描述一群膠體粒子的泳動。此外 Levine 和 Neale[12] 亦考慮了粒子間的交互作用，及任意電雙厚度，但需假設其粒子表面為低表面電

位。如此而得到相當複雜的解析解。在他之後Kozak 和 Davis[14,15]則推導了一

個更為普遍的理論，可適用於任意表面電位及任意電雙層厚度，但忽略電雙層重疊(Double layer overlapping)效應。而 Ohshima[16]則在 1997 年時利用近似法 (Approximation method)推導出密集球形膠體粒子在懸浮液中之電泳速度，他綜合

上述二組作者之假設，而得到相對誤差在4%以下的電泳速度。

Lee et.al.[17-20]在近幾年中，認為邊界效應對膠體之電泳運動有相當的重要性，乃發展出一套有效率之假性光譜法(pseudo-spectral method)[17-20]，並以孔洞模型(cavity)[21]及 Kuwabara[13]所提出之單一晶格模型(unit cell model)來探討膠體粒子與邊界之間的交互作用，其結果適用在高表面電位、任意電雙層厚度以及考慮了電雙層扭曲的效應，實以大幅提升電泳運動理論應用的範圍。上述之邊界條件理論大多討論球形對襯並且可以一維化(one dimensional)的情形。若考慮一膠體粒子受一平板之邊界效應，Morrison 和 Stukel[22]發表了一篇一球形膠體粒子垂直於一平板之電泳運動，文中提到他們主要之目的為修正 Smoluchowski Eq.(1.1)式在邊界時的應用，故其假設球之半徑遠大於電雙層厚度，且膠體粒子表面之電位夠低，以至可以線性化方程式，使用雙球座標求解析解，其所得之結論為邊界之影響會降低電泳速度，而隨著距離的接近效應越大。之後Keh 和 Lien[23]指出 Morrison 所求之解有錯誤存在，故 Keh 使用了相同的

(12)

平板之電泳、膠體粒子平行於平板之電泳及膠體粒子在一圓柱形孔洞中之電泳行為，文章仍是以無窮小之電雙層厚度及低表面電位，使用Reflection 方法求得解析解，文中指出有邊界之影響有三種基本效應分別是: (1)邊界對於流力部分會造成其速度下降 (2)邊界會改變粒子和外加電場之間的作用力 (3)若邊界帶一電位則會引發電滲透效應。 之後Ennis 和 Anderson[25]做和上述之論文相同幾何系統之電泳分析，可以 計算有限電雙層厚度，但仍限制於低表面電位及電雙層厚度不可以接觸到平板，文中提到有限電雙層厚度更是加強了其降低粒子速度之效應。Feng 和 Wu[26]則是考慮一 Prolate 物體以旋轉方式向可傳導平板前近，但只適合在電雙層很薄以

及低表面電位的情形下。最近Shugai 和 Carnie[27]以 Teubner’s[28]方法求解一個膠體粒子垂直一平板之電泳運動，電雙層的厚度沒有限制但膠體表面電位仍為低電位，主要為修正Henry’s 結果，文中針對了球和平板帶固定表面電位及固定電荷密度來觀察其現象，並証明在某些情況下此方法較Reflection 方法確準。在 2002 年，Sellier[29]提出一套新的計算流程，針對平行或垂直平板的電場作用下，一任意形狀及任意表面電位粒子的電泳現象，但是仍舊限制在電雙層厚度很薄的時候；Tang et al.[29]進一步將表面固定電位的限制解除，探討當膠體粒子表面產生解離反應，即所謂Charge regulation model 時電泳現象的變化，而 Chih et.al.[30] 更廣泛地在任意表面電位和電雙層厚度下，對於球形膠體粒子垂直於平板之電泳行為作進一步地探討。因此，根據 Tang 與 Chih 的基礎，本研究室將這個理論更進一步推展至一長直圓柱粒子對一平板的電泳現象，這也是本次三年計畫的第一年的計畫進度。然而上述文獻相關研究，都集中於弱外加電場的探討。但是在工業應用上，為了增加產率或其他原因，大部分的外加電場的大小往往是不可忽略的，這時候

(13)

採用 O’Brien 及 White[9]線性化的理論基礎將會產生嚴重的誤差。因此在 2002

年由Lin et al.[31]等人把整個系統的問題拓展至任意外加電場的行為，亦即不經

過線性化直接求解問題的非線性主控方程式，而在本計畫的第二年裡，我們也將 Lin et al.[31]與 Chih et al.[30]的結果做一合併，將球對平板的結果推廣到任意外加電場的行為。

以上垂直於平板運動的電泳現象都是假設粒子為光滑不可穿透的硬球粒子，可是截至目前為止尚未出現關於複合粒子或是液體球對平板的電泳現象之研究。液滴在理論分析上與固體粒子相似，但由於內部流場的影響，使得其分析方法更為複雜，故比起一般固體粒子，此領域的文獻較少。關於液體球的電泳最早是由Craxford et al. [32]於 1937 年研究水銀液滴在電場下之運動；之後 Booth[33] 在不考慮鬆弛效應、弱外加電場，以及低表面電位時，對不同液滴內部電荷的分

布研究，得到電泳速度與電雙層厚度、及液滴黏度之間的關係；但是Levich and

Frumkin[34]則是提出完全不同的理論來描述電雙層很薄時之電泳現象，並且其結果與Booth 的結果差距甚大。Levine and O’Brien[35]曾對此二結果作比較，發

現在水銀液滴的電泳現象中因為Booth 並沒有考慮表面偶極的作用，Levich and

Frumkin 的結果較為準確。往後有關於液滴粒子的研究多偏重於水銀液滴，這是因為水銀本身特殊之物理性，使得理論模型的建立上更為理想。例如 Ohshima[36,37]曾探討在稀薄懸浮系統中，水銀液滴表面電位很高的時候，其電泳速度與液滴內外之黏度無關，此稱之為液滴的「固化現象」(10)。

近年來，Baygents and Saville[38]以數值方法對弱外加電場下液滴之電泳行為做了相關的理論研究，其結果包含了非傳導（液滴內部無電解質）與可傳導（液滴內部有電解質）液滴之電泳行為，根據他們的研究指出，當液滴的表面電位增加時，電泳速度將隨之增加且產生一極大值；當電雙層厚度增加時，電泳速度遞

(14)

更為複雜，甚至在某些情形下，例如液滴內、外的電雙層都相當薄的時候，液滴甚至會出現和它所帶的表面電位反向的情形，這是由於其液滴的表面出現了一層緩衝層，其流動的現象也較非傳導液滴更為複雜。此外，Baygents and Saville[39] 也考慮弱電解質溶液中電雙層厚度會受到解離常數影響，針對弱電解質溶液中的

電泳行為加以計算。而 Ohshima 也延續過去[37]探討稀薄懸浮系統的結果，以

Kuwabara[13]之單位晶格模型，在低表面電位，不考慮電雙層重疊的情形下，探討密集帶電水銀液滴之各種電動力學現象[40]，並發現當液滴黏度趨近於無窮大時，其電泳現象可以視之為等同於剛性球體(rigid sphere; hard sphere)，隨後更將此一系統衍生至無鹽溶液(salt free medium)的情形 [41]。Lee et al.[42-44]也根據 O’Brien and White 所提出在弱外加電場下之計算方式，分別探討任意表面電位下的非傳導液滴以及密集水銀液滴的電泳行為，並且指出極化效應以及液滴與外邊界孔隙度的大小將如硬球粒子一般影響電泳速度的結果，並且當內部流體黏度遠高於外部溶液時，其電泳速度的結果將可回歸到硬球的結果；反之當內部流體黏度遠高於外部溶液時，計算出的電泳速度結果將趨近於氣泡(Bubble)。另外於 2005 年 Lee 等人亦將計算模型推廣至非牛頓高分子液滴的電泳情形[45]，強調高分子液滴的 shear thinning 效應造成液滴黏度下降，而電泳速度愈偏離牛頓流體的情況，此現象隨電雙層厚度減少而更加顯著。如果考慮一球形液滴對平面的運動現象，Bart[46]在 1968 年求解一球形液滴受固定外力作用而產生垂直平面運動的速度大小；Wacholder 與 Weihs[47]在 1972 年延續 Bart 的結果直接求解兩顆液滴互相運動靠近以及液滴垂直於平板或自由液面的內外流場解析解；Grashchenkou[48]更進一步探討液體球表面為可滑動面時的結果；而 Chen 與 Keh[49]等人則更進一步探討溫度變化造成液滴對平板的熱泳現象；Kasumi[50]等人則在 2004 年以實驗的方式比較 Maragoni 效應對氣泡熱泳的影響。但是關於電泳現象之理論，過去由於受限於液滴的電泳現象需要同時求解內部與外部流體的運動情形，加上電泳的主控方程式較為複雜，造成計算

(15)

上的困難，因此目前為止尚未有關於此一方面的研究。所以，我們在第三年的計畫中，將球對平板的電泳現象拓展至低表面電位下液滴對平板的系統。

(16)

第二章

理論分析

2-1 系統描述

若考慮系統為一球形帶電非導體粒子，如圖2-1，其半徑為a，在一電價比為 z1:z2的電解質溶液中，z1，z2分別表示為陽離子和陰離子的電價數。為了表示溶液的電中性故陽陰離子之濃度比為n20=n10/α，α為離子價數比( 1 2 z z − )，而此液體球和平板之距離為h，如圖 2-1 所示。在此一系統中我們採用雙球座標或是雙圓柱座標

(

η ,,ξ θ

)

(

η,ξ,x

)

來求解主控方程式，如圖 2-2 ，其定義範圍為 π ≤ ≤ π ≤ θ ≤ π ≤ ξ ≤ ∞ < η ≤ , 0 ,0 2 ,0 x 2 0 [51]。其中η=η₀表示為球表面(η 的定₀ 義為cos−1(h/a))；η=0表示為平板；η>η₀則表示球內部；而ξ=0,π代表對稱軸。其與極座標的關係為 ξ − η η = cos cosh sinh c z ξ − η ξ = cos cosh sin c y 其中c 為焦距。在 z 方向有一外加電場E_Z使得此膠體粒子以一速度U 向上運動。

(17)

Fig.2-1 Geometric configuration of the system in this study

(18)

當圖一的粒子為一圓柱或球形粒子時，我們的假設如下： (1) 粒子外為不可壓縮牛頓電解質流體。 (2) 粒子移動緩慢，故可忽略對流項，流體流動可視為緩流(Creeping Flow)；而流體的物性如黏度、介電常數等為固定值。 (3) 整個系統已達到了擬穩態(quasi-steady state)的狀態，故可忽略方程式中的時間項。 (4) 圓柱粒子為一長直圓柱，亦即圓柱本身的長度遠大於其截面半徑。 (5) 粒子表面均帶有固定電位，並且如果外加電場對離子分佈之影響遠小於粒子表面電位，則可將主控方程式作進一步化簡；反之若外加電場的大小不可忽略，則不能做任何化簡。另外，若粒子本身為一液體球，則系統必須滿足以下的假設： (1) 此液滴內外均為不可壓縮牛頓流體，且液滴本身為非電解質溶液。 (2) 液滴移動緩慢，故可忽略對流項，流體流動可視為緩流(Creeping Flow)；而流體的物性如黏度、介電常數等為固定值。 (3) 液滴維持球形對稱的性質，不受外加電場影響改變其形狀。 (4) 整個系統已達到了擬穩態(quasi-steady state)的狀態，故可忽略方程式中的時間項。 (5)液滴表面帶有固定電位，並且外加電場對離子分佈之影響遠小於液滴表面電荷所造成之反向電場。

(19)

2-2 主控方程式基本介紹 計算一膠體粒子受外加電場所產生的電泳運動須考慮電場、流場及電解質濃度場其隨位置的變化和相互影響。故需電位方程式，流場方程式以及離子守恆式等非線性耦合之電動方程組(Electrokinetic Equations)。其內容詳述如下： (i)電位方程式因為我們假設流體流動緩慢時間項可以忽略，故電位方程式可以用Poisson 方程式來表示

∑

= ε − = ε ρ − = φ ∇ 2 1 j j j 2 z en _(2-1)

其中ρ、ε分別為空間電荷密度(space charge density)及介電常數(permittivity)，而

為電解質之物種j 之電荷數、e 為單一電子帶電量、為電解質之物種j 之濃度。(2-1)式即所謂的高斯定律，表示通過任一封閉曲面的電場通量，與其封閉曲面內的總電荷數成正比。式中的介電常數 j z n_j

ε

可以視為高斯定律在不同介質中的修正項。考慮電場通過物質時，物質本身具有永久電偶極矩（permanent dipole moment）或是分子本身雖為非極性，但在電場的影響下產生感應電偶極矩 (induced dipole moment)。而這些偶極矩通常與外加電場的方向相反，減弱了電場在物質內的效應。 (ii)流場方程式在緩流的假設之下流場可由Navier-Stoker Equation 表示，且此一流體須滿足不可壓縮之牛頓流體，故可表示如下 (2-2) 0 = ∇

⋅

v (2-3) 0 p 2 ₋_∇ ₋_ρ_∇_φ₌ ∇ µ v 其中p 和 µ 分別表示流體之壓力和流體之黏度，而 ρ 之定義和(2-1)式相同；另外為單位體積下粒子所受之電力。 φ ∇ ρ

(20)

在擬穩態的情況下，離子的守恆式可以下式表示 0 = ⋅ ∇ J_j (2-4) 其中Jj代表物種j之密度通量，可表示成如下的型態： v φ n J_j _j _j B j j j _k _T n en z D _⎥+ ⎦ ⎤ ⎢ ⎣ ⎡ ∇ + ∇ − = (2-5) 等號右邊第一項表示為擴散(diffusion)通量的貢獻、第二項表示為傳導(migration) 通量貢獻，而第三項表示為對流(convection)通量的貢獻。合併上述兩方程式，即可得j 物種的濃度守恆方程式為 0 n T k en z D _j B j j j = ⎪⎭ ⎪ ⎬ ⎫ ⎪⎩ ⎪ ⎨ ⎧ − ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ ∇ + ∇ • ∇ n_j φ v j=1,2 (2-6)

其中Dj是物種j的擴散係數(Diffusion coefficient )，kB為Boltzmann常數，T為絕對溫度。在本研究中，如果外加電場的強度所造成的電位與粒子表面電位相當，則我們必須直接求解(2-1)、(2-2)、(2-3)、(2-6) 等式。反之如果吾人假設外加電場的強度，相較於粒子表面所帶電荷形成之電場，是屬於較弱的。因此，可將待解的主要變數寫成平衡狀態(未施加電場時)加上擾動狀態的形式。並根據弱電場的假設，知外加電場對擾動的影響可視為線性，亦即: z E z e( , ) ( , )E ) , (ηξ =φ η ξ +δφηξ φ (2-7) z j e j j( , ) n ( , ) n ( , )E n η ξ = ηξ +δ ηξ (2-8) z )E , ( 0 ) , (η ξ = +δv η ξ v (2-9) z e( , ) p( , )E p ) , ( p η ξ = ηξ +δ η ξ (2-10) 其中「δ 」代表擾動狀態。綜觀上述四式，所有系統變數只與η、ξ有關，此乃考量此物理問題的幾何

(21)

系統乃為軸對稱。另外值得一提的是，我們可以將系統的電位φ視為兩分量的線性相加[17,18]，分別是沒有外加電場影響時的平衡電位φ 以_e ，及由外加電場所造成的電位δφ 亦即φ=φ_e +δφ。綜合以上之敘述我們知道系統為一二維軸對稱系統，和時間無關。而系統待求解之變數有平衡電位(φ )、外加電位(_e δφ)、系統速度(v)、系統壓力(p)和電解質濃度(nj)，這適用於第一年度「圓柱對平板之電泳現象」以及第三年度「液滴對平板之電泳現象」均屬此一範疇。故為了簡化求解之複雜性及快速得到合理之解，我們探討不同情形下的電泳運動時，係根據以上之主控方程式做不同的簡化和推導而得。所以接下來，我們將對電位、濃度分佈、流場方程式做進一步化簡。

(22)

2-3 平衡狀態: 延續前面，我們將待解的系統電位分為平衡狀態及外加擾動狀態的形式，故需分別探討此兩種狀態，以得到整體系統的解。在這裡，我們先討論平衡時系統的電位分佈。根據前一節的描述，我們知道當系統達平衡時，即意指沒有外加電場的施加，若假定膠體粒子只受所帶電位的影響且靜止不動。故我們不需解流場，只需同時求解電場及濃度場。由(2-1)式，可知

∑

= ε − = φ ∇ 2 1 j j j e 2 z en _(2-11) 其中拉普拉斯運算子(Operator of Lapalace )乃是根據雙球座標或雙圓柱座標下的定義，在此定義下的運算子符號可參照Appendix A 的推導。而由(2-6)式，考慮系統乃為平衡，故J_j =0，v=0，故將此兩式合併，可得 0 T k ez n D B j j j = ⎟ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎜ ⎝ ⎛ ∇ + ∇n_j φ_e (2-12) 將(2-12)式作積分，可得 ) T k e z exp( n n B e j 0 j j φ − = (2-13)

而(2-13)式裡的為整體的離子濃度(Bulk phase ion concentration)，而代表。此式即為著名的波茲曼分佈(Boltzmann distribution)。將(2-13)式代入(2-11)式，可得 j n e 0 j n

∑

= φ − ε − = φ ∇ 2 1 j B e j 0 j j e 2 ₎ T k e z exp( en z (2-14) 此即著名的 Poisson-Boltzmann 方程式，從以上的推導可以發現此式為電場及濃

(23)

度場的耦合式。因此只要解得此式，便可了解系統平衡時的情形。然而，在第三年裡考慮液滴表面電位相當的低，則(2-13)式的指數形式可以使用泰勒級數展開並消去高次項後，線性化成下式 ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ φ − = T k e z 1 n n B e j 0 j j (2-15) 則平衡電位方程式可寫成

∑

= ⎟⎟⎠ ⎞ ⎜⎜ ⎝ ⎛ φ − ε − = φ ∇ 2 1 j e j 0 j j e 2 kT e z 1 en z (2-16) 其對應的邊界條件為 ζ = φe , η=η0 (2-17) 0 e = φ , η=0 (2-18) 0 e ₌ ξ ∂ φ ∂ , ξ=0,π (2-19) 0 e = φ , η→∞ (2-20) 其中r 為方向向量(r≡ )，在雙球座標的定義下rrˆ cosh cos ⎟1/2 ⎠ ⎞ ⎜ ⎝ ⎛ + = x c r η ξ 。 (2-17)式代表假設電荷均勻分佈在液滴表面且具有固定電位ζ ；(2-18)式則說明平板上不具電荷；而(2-19)式為滿足系統是軸對稱。另外如果我們考慮一非傳導液滴，最後(2-20)式則闡述了電荷不會累積在液滴內部。

(24)

2-4 擾動狀態: (2-13)或(2-15)式考慮的是平衡狀態，無流動的情況下的離子分佈。若產生外加作用力時，流體流動以及外加電場則有可能對電位分佈造成影響。又因為在本研究中，為了區分高低表面電位的不同，這裡根據年度計畫的不同作一個別之論述。 1. 第一年：圓柱對平板之電泳現象 在這一節我們要探討的是當粒子為一長直圓柱，表面為高電位的情形，在弱外加電場的作用下，此時溶液中的電解質濃度可以用下式來表示： ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ φ +δφ+ − = T k ) g ( e z exp n n B j e j 0 j j , j=1,2 (2-21) 此式和(2-13)式的波茲曼分佈不太相同。當表面電位的大小不可忽略時，電解質濃度除了受平衡電位(φe)的影響也受外加電場的影響(δφ)，gj表示電雙層的變形程度對於電解質濃度的影響。所謂電雙層的變形程度即是因流場的變化而使得電雙層產生變化，如此一來會在膠體附近產生一內電場，而更進一步地影響電解質的濃度分佈。上述效應即是一般所稱的極化效應。由於外加電場仍是一個擾動項，故其產生的外加電位也可視為一擾動電位。極化效應所產生的電位和外加電位是大小差不多，亦可視為一擾動項。而此兩擾動項的大小遠小於平衡電位，因此(2-21)式可以簡化成下式：

(

⎥ ⎦ ⎤ ⎢ ⎣ ⎡ + δφ − ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ φ − = _j B j e B j 0 j j _k _T g e z 1 T k e z exp n n

)

, j=1,2 (2-22) 如此空間電荷密度即可表示成下式：

(

∑

= ⎥⎦ ⎤ ⎢ ⎣ ⎡ + δφ − ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ φ − ≅ ρ 2 1 j B j j e B j 0 j j _k _T g e z 1 T k e z exp en z

)

(2-23) 將上式代入(2-1)式則可得高表面電位之電位方程式為：

(25)

(

)

∑

= ⎥⎦ ⎤ ⎢ ⎣ ⎡ + δφ − ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ φ − ε − = φ ∇ 2 1 j j B j e B j 0 j j 2 _g T k e z 1 T k e z exp en z (2-24) 首先了解平衡電位，由於平衡電位方程式的假設即是無外加電場的存在故其方程式和(2-14)式相同。接下來擾動電位方程式同樣可以由(2-24)式來求得，即是將(2-24)式減去(2-14)式。

(

)

(

j

)

B j 2 1 j B e j 0 j j 2 1 j B e j 0 j j 2 1 j B j j e B j 0 j j 2 g T k e z T k e z exp en z ) T k e z exp( en z g T k e z 1 T k e z exp en z + δφ ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ φ − ε = φ − ε + ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ + δφ − ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ φ − ε − = δφ ∇

∑

= = = (2-25) 在濃度分佈表示式中我們可以看到此處的計算較低表面電位多了二個變數 gj。由於多了此二個變數，因此需多加入二條主控方程式。我們將(2-24)式代入(2-4) 式則可以新的擾動電位方程式： e 1 1 e 1 1 2 D 1 g kT e z g − ∇ φ ∇ = ∇ φ ∇

⋅

v

⋅

(2-26) e 2 2 e 1 2 2 D 1 g kT e z g + α ∇ φ ∇ = ∇ φ ∇

⋅

v

⋅

(2-27) 同樣地在這裡我們也需要求解流場方程式，在高表面電位其流場方程式可以表示如下：(細述於 Appendix B) ⎥ ⎥ ⎦ ⎤ ξ ∂ φ ∂ ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ η ∂ ∂ α + η ∂ ∂ − ⎢ ⎢ ⎣ ⎡ η ∂ φ ∂ ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ ξ ∂ ∂ α + ξ ∂ ∂ ξ − = ψ e 2 e 2 1 e 1 e 2 e 2 1 e 1 B 10 2 2 1 4 g n g n g n g n T k n e z c sin x E (2-28) 其中 ) kT e z exp( n₁_e = − 1 φ_e 為平衡時陽離子的濃度，而 ) kT e z exp( n₂_e = α 1 φ_e 為平衡時陰離子的濃度。

(26)

為了計算上的方便我們需將上述高表面電位之主控方程式做無因次化。我們定義膠體球表面為一固定表面電位，其電位值為ζa。我們即根據此一電位值做無因次化可得

、

，將以上無因次變數分別代入(2-14)、(2-25)、(2-26)、(2-27)及(2-28) 式中即可得到無因次化主控方程式： a j * j =φ /ζ φ * _j ₁₀ j n /n n = * _j _a j g / g = ζ U* ₌U/_ζ_a2_ε/_µa µ ε ζ ψ = ψ* / _a2 a/ i.無因次化平衡電位方程式： )] exp( ) [exp( ) a ( ) 1 ( 1 e r e r r 2 * e 2 * ₋_φ _φ ₋ _αφ _φ φ κ α + − = φ ∇ (2-29) ii.無因次化外加擾動電位方程式：

( )

(

)

[

]

(

( )

)

[

1*e 1* *2e *2

]

2 * * e 2 * e 1 2 * 2 * _n _g _n _g 1 a n n 1 a α + α + κ = δφ α + α + κ − δφ ∇ (2-30) iii.無因次化流動擾動電位方程式： * e * * 1 2 r * 1 * * e * r * 1 2 * _g ₋_φ _∇ _φ _∇ _g ₌_φ _Pe _∇ _φ ∇

⋅

v

⋅

(2-31) * e * * 2 2 r * 2 * * e * r * 2 2 * _g ₊_αφ _∇ _φ _∇ _g ₌_φ _Pe _∇ _φ ∇

⋅

v

⋅

(2-32) iv.無因次化流場方程式 ⎥ ⎥ ⎦ ⎤ ξ ∂ φ ∂ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ η ∂ ∂ α + η ∂ ∂ − ⎢ ⎢ ⎣ ⎡ η ∂ φ ∂ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ ξ ∂ ∂ α + ξ ∂ ∂ ξ α + κ − = ψ * e * 2 * e 2 * 1 * e 1 * e * 2 * e 2 * 1 * e 1 * 2 * 4 g n g n g n g n c sin x ) 1 ( ) a ( E (2-33) 其中Pej =U/aDj

、

n1*e =exp(−φrφe*)

、

n*2e =exp(αφrφ*e)

、

a B 1 r T k e z ζ = φ

。

在高表面電位下膠體電泳的邊界條件大致和之前並無太大的不同，在平衡電位方面我們若假定膠體粒子表面帶一固定電位值ζa，則可表示成： a e =ζ φ at η=η0 (2-34)

(27)

將其無因次化 at η=η 1 * e = φ 0 (2-35) 此外在無窮遠處我們假設平衡電位為零，而我們考慮平板為不帶電故其邊界條件可寫成(2-36)式，而式(2-37)表示軸對稱邊界條件。 at η=0 and ξ=π or η=0 (2-36) 0 * e = φ 0 * e ₌ ξ ∂ φ ∂ at ξ=0 or ξ=π (2-37) 在外加擾動電位方面，我們假設膠體粒子的表面傳導可以忽略。且一般而言膠體粒子的介電常數通常小於流體，故可以假定膠體粒子為一非導體粒子。所以電場效應不可穿透，故膠體粒表面上之邊界條件就可以寫成(2-38)式。此外我們從無窮遠處沿著 z 軸外加了一擾動電場，故其邊界條件為(2-39)式所示。而式(2-40)表示軸對稱邊界條件。 δφ 0 * = η ∂ δφ ∂ at η=η0 (2-38) 0 * = ξ ∂ δφ ∂ at ξ=0 or ξ=π (2-39) at η=0 and ξ=0 or η=0 (2-40) * * z * ₌₋_E _⋅_r δφ 其中Ez*為外加於z方向之電場，r為方向向量定義為r*⋅ir而r在雙圓柱座標為 2 / 1 * * x cos cosh c r ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ η+ ξ = (2-41)

(28)

而所造成的擾動電位。在流動擾動電位的邊界條件上，在膠體粒子表面需滿足電解質不可穿透粒子表面，故∂n_j ∂n應為零。其中nj表示電解質第j種物種之濃渡， n表示垂直方向。又因平衡電位和外加擾動電位在膠體粒子表面皆是對垂直方向為零[(2-38)式]，故流動擾動電位邊界條件在膠體粒子表面可表示成： 0 g* j ₌ η ∂ ∂ at η=η0 (2-42) 在無窮遠處和在平板表面要滿足擾動為零： at η=0 and ξ=0 or η=0 (2-43) * * j g =−δφ 在對稱邊界條件須滿足軸對稱故寫為 0 g* j ₌ ξ ∂ ∂ at ξ=0 or ξ=π (2-44) 對於流場方程式之邊界條件而言，我們考慮一個膠體粒子垂直於一平板以一速度U 運動。根據 Happel(51)書中指出，在一軸對稱的運動中粒子的表面之邊界條件有二種型式(kinematical 及 dynamical)：

(

v−U

)

⋅n=0 kinematical condition (2-45)

(

v−U

)

⋅s=0 dynamical condition (2-46) 若將上二式轉成流線函數，並無因此化之則可轉變成以下表示式： * 2 * * _U 21 ω − = ψ * 3 2 * * U sinh sin x c _ξ _η = η ∂ ψ ∂ at η=η0 (2-47) 而無遠處溶液以及平板都靜止不動，因可表示成 0 * ₌ ψ 0 * = η ∂ ψ ∂ at η=0 (2-48) 0 r*2 * → ψ at η=0 and ξ=0 (2-49) 軸對稱邊界條件可表示成 0 * ₌ ψ 0 * = ξ ∂ ψ ∂ at ξ=0 or ξ=π (2-50)

(29)

其中ω* =c*sinξ/x而U* ₌U/

(

kT ez₁

)

2_ε/_µa 2. 第三年：液滴對平板之電泳現象 由於液滴必須同時考慮內外兩區流場的變化，因此為了簡化我們的計算，我們討論低表面電位下非傳導液滴對一平板之電泳現象。這裡先說明一點，由於第二年的強電場的處理沒有經過簡化，故在這裡我們將順序作一對調。 (2-15)式考慮的是平衡狀態，無流動的情況下的離子分佈。若產生外加作用力時，流體流動以及外加電場則有可能對電位分佈造成影響。然而在本研究裡吾人假設液滴的表面電位相當低，外加電場的強度遠小於液滴表面離子雲扭曲所造成的影響，也就是說，液滴外部的電解質濃度分佈僅受到液滴表面電位的影響，故這裡假定(2-15)式依舊成立。將(2-7)式取∇2後，可得： e 2 2 Z 2₍_δφ_E ₎₌_∇ _φ₋_∇ _φ ∇ (2-51) 將(2-1)式、(2-11)式、與(2-16)式代入(2-51)，則我們可以得到 0 ) kT e z 1 ( ) kT e z 1 ( en z ) E ( 2 1 j e j e j e 0 j j Z 2 = ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ φ − − φ − ε = δφ ∇

∑

= (2-52) 其對應的邊界條件為 0 = η ∂ δφ ∂ , η=η₀ (2-53) r E_z ⋅ − = δφ , η=0 (2-54) 0 = ξ ∂ δφ ∂ , ξ=0,π (2-55) , 0 = δφ η→∞ (2-56)

(30)

上述之第一個邊界條件代表液滴內部為非傳導體，其介電常數通常遠小於外部流體，故在液滴表面可給定電荷不可穿透之條件；第二個邊界條件則說明外加電場的效果；而第三個邊界條件為滿足系統軸對稱邊界條件；最後(2-26)式則闡述外加擾動電場不會影響內部非傳導流體。　同樣的，接著，我們對(2-3)式取旋度，可消去壓力項得： 0 ) ( ) v (_µ_∇2 ₋_∇_× _ρ_∇_φ ₌ × ∇ v (2-57) 而對於一個2D 的流動問題，我們通常會引入流線函數來減少未知數的數目。在雙球座標系統下，其流線函數之定義為： ξ ∂ ψ ∂ ξ = η sin c x v ₂ 2 (2-58) η ∂ ψ ∂ ξ − = ξ sin c x v ₂ 2 (2-59) 將(2-15)式代入 ρ 的定義裡，並且進一步把(2-28)、(2-29)式代入(2-27)式裡，詳細之計算情形可參考Appendix B，最後可得： ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ ξ ∂ φ ∂ η ∂ ρ ∂ − η ∂ φ ∂ ξ ∂ ρ ∂ µ ξ = ψ c sin x E4 (2-60) 其中運算子E4的定義為： 2 2 4 _E _E E = ； _⎥ ⎦ ⎤ ⎢ ⎣ ⎡ ξ ∂ ∂ ξ η ξ − + η ∂ ∂ η + ξ ∂ ∂ + η ∂ ∂ = sin x cosh cos 1 x sinh c x E ₂ 2 2 2 2 * 2 2 _(2-61) 這樣一來，本來需要解的兩個速度分量以及壓力項變成僅需求解流線函數即可。然而這裡需要注意的是，由於液滴內外均有流體，因此我們必須同時求解液滴內外的流場分佈，而根據式(2-30)，本研究的流場主控方程式為： ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ ξ ∂ φ ∂ η ∂ ρ ∂ − η ∂ φ ∂ ξ ∂ ρ ∂ µ ξ = ψ(e) ₍_e₎ 4 c sin x E (2-62) 0 E4ψ(i) = (2-63)

(31)

其中上標、代表液滴外部以及內部；、代表液滴外部電解質溶液的流線函數及黏度；、代表液滴內部非電解質溶液的流線函數及黏度，而且因為內部是非電解質溶液，故式(2-3)與式(2-60)的電力項不存在，只需要求解緩流的流線分佈。 ) e ( (i) ψ(e) µ(e) ) i ( ψ µ(i) 對於流場方程式之邊界條件，我們考慮一個液體球垂直於一平板以一速度U

運動。根據Wacholder and Weihs[47]指出，在液滴表面上的邊界條件應為：

0 ) (v(e) −U ⋅i_n = (2-64) 0 n ) i ( _{⋅ i} ₌ v (2-65) t (i) _i v i U v(e) − )⋅ _t = ⋅ ( (2-66) ) e ( nt ) i ( nt =π π (2-67) 其中下標符號n 、t 代表「與液滴表面」的垂直方向及切線方向(相當於球座標之 r 、方向)；上標符號、則分別表示液滴內部、外部的特性；代表單位向量；則為液滴表面的切應力。式(2-64)與式(2-65)說明了液滴維持球形不破裂，且內外流體不會自由進出；而因為液滴本身為內外兩種不同流體的界面，故式(2-66)代表交界面內外切線速度連續；式(2-67)則闡述了切應力(tangential stress

component)連續。而根據 Happel 與 Brenner[51]將以上四個邊界條件轉換成流線

函數以及雙球座標的定義： θ (i) (e) i nt π U 2 1 ₂ ) e ( ₌ − _ω ψ ，at η=η₀ (2-68) 0 ) i ( ₌ ψ ，at η=η₀ (2-69)

(32)

U sinh sin x c 2 3 2 ) i ( ) e ( η ξ = η ∂ ψ ∂ − η ∂ ψ ∂ ，at η=η₀ (2-70) ) e ( ) i ( ηξ ηξ =π π ，at η=η₀ (2-71) 其中ω=csinξ/x，代表與對稱軸之水平距離；x =coshη−cosξ；而雙球座標下切應力的定義為：

(

)

⎪ ⎪ ⎪ ⎭ ⎪⎪ ⎪ ⎬ ⎫ ⎪ ⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎨ ⎧ η ∂ ψ ∂ − ∂ ψ ∂ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ + − + + ∂ ψ ∂ ξ − ξ ∂ ψ ∂ ξ ⋅ µ − = π_ηξ ) i ( 3 2 ) i ( 3 2 2 3 2 2 2 ) i ( 2 3 3 2 ) i ( 2 3 3 ) i ( ) i ( ξ sin c η sinh x 3 ξ c x ξ sin c η cosh ξ cos ξ sin 1 x η sin c x sin c x (2-72)

(

)

⎪ ⎪ ⎪ ⎭ ⎪⎪ ⎪ ⎬ ⎫ ⎪ ⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎨ ⎧ η ∂ ψ ∂ − ∂ ψ ∂ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ + − + + ∂ ψ ∂ ξ − ξ ∂ ψ ∂ ξ ⋅ µ − = π_ηξ ) e ( 3 2 ) e ( 3 2 2 3 2 2 2 ) e ( 2 3 3 2 ) e ( 2 3 3 ) e ( ) e ( ξ sin c η sinh x 3 ξ c x ξ sin c η cosh ξ cos ξ sin 1 x η sin c x sin c x (2-73) 另外，在液滴內部的速度是存在的，因此： 0 ) i *( = η ∂ ψ ∂ ，ψ*(i) =0 when η→∞ (2-74) 而在無窮遠處溶液是靜止不動的，因此可表示成 0 r2 ) e ( → ψ at η=0 and ξ=0 (2-75) 軸對稱邊界條件可表示成 ψ(e) =0， (e) =0 ξ ∂ ψ ∂ at ξ=0 and ξ=π (2-76) 最後，針對流體在實體平板上滿足不可滑動(no-slip)邊界條件，亦即： 0 v_η = at η=0 (2-77) 0 v_ξ = at η=0 (2-78)

(33)

而上述兩式可以進一步化簡成： 0 ) e ( ₌ ψ at η=0 (2-79) 0 ) e ( = η ∂ ψ ∂ at η=0 (2-80) 同理，將前述之無因次群代入之前所提到之主控方程式與邊界條件，分別進行無因次化的處理，即可將系統變數無因次化，經過一些計算以及化簡之後，分別整理成以下的表示式；並且為了描述液滴內外不同流體特性，這裡吾人亦將前述之主控方程式整理成液滴內部與液滴外部的形式：根據式(2-16)，在外部電解質溶液的無因次化平衡電位方程式為： (2-81) * e 2 * e 2 * _φ ₌₍_κ_a₎ _φ ∇ 其中: ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ η ∂ ∂ η − ξ ∂ ∂ ξ − η ξ + η ∂ ∂ + ξ ∂ ∂ = ∇ x sinh sin x 1 cosh cos c x 2 2 2 2 2 * 2 2 * _(2-82) * e φ 代表無因次化以後的平衡電位；符號「*」代表無因次化後的變數，對底下的擾動電位以及流場方程式亦然；

κ

a為粒子半徑與電雙層厚度的比值，為一重要

的無因次群。又電雙層厚度為(Debye length),

κ

−1，而其倒數(inversed Debye length)

κ

表示如下 :

( )

1/2 2 1 2 0 / ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ =

∑

= j j j ez kT n ε κ (2-83) 而我們假設液滴本身是一非導體球，換句話說液滴內部的流體不含電解質，且電荷均勻分佈在表面上。因此液滴內部所對應的平衡電位主控方程式為： * 2 *

(34)

其對應之無因次化邊界條件為： 1 * e = φ at η=η₀ (2-85) 0 * e = φ at η=0 (2-86) 0 * e ₌ ξ ∂ φ ∂ at ξ=0,π (2-87) 0 * e = φ at η→∞ (2-88) 由式(2-52)，液滴外部或內部的無因次外加擾動電位都可以寫成： 0 ) ( * 2 * _δφ ₌ ∇ (2-89) 其中運算子∇*2之定義同前面(2-82)式。另外其對應之邊界條件為： * =0 η ∂ δφ ∂ at η=η₀ (2-90) δφ* =−E_z*⋅r* at η=0 (2-91) * =0 ξ ∂ δφ ∂ at ξ=0,π (2-92) δφ* =0 at η→∞ (2-93) 根據式(2-62)與式(2-63)，液滴外部與液滴內部的無因次化流場方程式可以表示成： ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ ξ ∂ δφ ∂ η ∂ φ ∂ − η ∂ δφ ∂ ξ ∂ φ ∂ ξ κ − = ψ *e * *e * * 2 ) e ( * 4 * c sin x ) a ( E (2-94) 0 E*4ψ*(i) = (2-95) 其中 : 2 * 2 * 4 * _E _E E = (2-96)

(35)

⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ ∂ ∂ + ∂ ∂ − + ∂ ∂ + ∂ ∂ = η η ξ ξ η ξ η ξ x x c x E sinh sin cosh cos 1 2 2 2 2 2 * 2 2 * _(2-97) 而式(2-68)到式(2-71)的液滴表面流場邊界條件亦可改成： * 2 * ) e ( * _U 21 ω − = ψ at η=η₀ (2-98) 0 ) i *( ₌ ψ at η=η₀ (2-99) * 2 3 2 * ) i *( ) e *( U sinh sin x c η ξ = η ∂ ψ ∂ − η ∂ ψ ∂ at η=η₀ (2-100) ) e *( ) i *( ηξ ηξ =π π at η=η₀ (2-101) 而內外無因次shear stressπ*(_ηξi)、π*(_ηξe)為：

(

)

⎪ ⎪ ⎪ ⎭ ⎪⎪ ⎪ ⎬ ⎫ ⎪ ⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎨ ⎧ η ∂ ψ ∂ − ∂ ψ ∂ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ + − + + ∂ ψ ∂ ξ − ξ ∂ ψ ∂ ξ ⋅ µ − = π_ηξ ) i ( * 3 * 2 ) i ( * 3 * 2 2 3 * 2 2 2 ) i ( * 2 3 * 3 2 ) i ( * 2 3 * 3 * r ) i ( * ξ sin c η sinh x 3 ξ c x ξ sin c η cosh ξ cos ξ sin 1 x η sin c x sin c x (2-102)

(

)

⎪ ⎪ ⎪ ⎭ ⎪⎪ ⎪ ⎬ ⎫ ⎪ ⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎨ ⎧ η ∂ ψ ∂ − ∂ ψ ∂ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ + − + + ∂ ψ ∂ ξ − ξ ∂ ψ ∂ ξ ⋅ − = π_ηξ ) e ( * 3 * 2 ) e ( * 3 * 2 2 3 * 2 2 2 ) e ( * 2 3 * 3 2 ) e ( * 2 3 * 3 ) e ( * ξ sin c η sinh x 3 ξ c x ξ sin c η cosh ξ cos ξ sin 1 x η sin c x sin c x (2-103) 另外，在液滴內部的邊界條件可從式(2-74)得到： 0 ) i *( = η ∂ ψ ∂ ，ψ*(i) =0 when η→∞ (2-104) 無窮遠處溶液的邊界條件改為： 0 r2 ) e *( → ψ at η=0 and ξ=0 (2-105) 而對稱軸上的邊界條件可化簡成：

(36)

ψ*(e) =0， *(e) =0 ξ ∂ ψ ∂ at ξ=0 and ξ=0 (2-106) 最後，平板上的邊界條件也可化簡成： 0 ) e *( ₌ ψ at η=0 (2-107) 0 ) e *( = η ∂ ψ ∂ at η=0 (2-108) 2. 第二年：強外加電場作用下球對平板之電泳現象 最後，考慮第二年的強電場，由於外加電場不可忽略，換句話說我們必須利用(2-21)式直接求解，亦即不經過任何線性化的處理。因此該年度所討論的平衡電位方程式φe變成： )] e e [( ) 1 ( ) a ( * e r * e r r 2 * e 2 −φ φ ₋ αφ φ φ α + κ − = φ ∇ (2-109) 而擾動電位方程式δφ則可寫成： )] e e ( ) n n [( ) 1 ( ) a ( * e r * e r * 2 * 1 r 2 * e 2 * 2 * 2 ₋ ₋ −φ φ ₋ αφ φ φ α + κ − = φ ∇ − φ ∇ = δφ ∇ (2-110) 而如果把先前的、、、代入離子守恆式中，並且將不同的離子分開，可得底下兩式： φ 2 ∇ ∇2φ₁ ∇2φ₂ n_j ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∇ ⋅ ∇ + δφ ∇ φ + ∇ + δφ ∇ + φ ∇ = ∇ φ ∇ φ − ∇

⋅

1 * 1 * r 1 * e 1 1 e r * 1 2 g ) g ( ) g Pe g g v

(

(2-111) ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∇ ⋅ ∇ + δφ ∇ φ − ∇ + δφ ∇ + φ ∇ = ∇ φ ∇ φ − ∇

⋅

2 * 2 * r 2 * e 2 2 e r * 2 2 g ) g ( ) g Pe g g v

(

(2-112) 而無因次化流線方程式為：

(37)

⎥ ⎥ ⎦ ⎤ ξ ∂ φ ∂ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ η ∂ ∂ α + η ∂ ∂ − ⎢ ⎢ ⎣ ⎡ η ∂ φ ∂ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ ξ ∂ ∂ α + ξ ∂ ∂ ξ α + κ − = ψ * * 2 * e 2 * 1 * e 1 * * 2 * e 2 * 1 * e 1 * 2 * 4 g n g n g n g n c sin x ) 1 ( ) a ( E (2-113) 而對應之邊界條件可參考第一年度的計畫。

(38)

2-5 電泳速度之計算:

根據O’Brien and White[9]之理論，因為第一年與第三年的主控方程式式皆為線性方程式，我們可以將本系統的問題分成兩個子問題(sub-problem)來看。第一個子問題（以下簡稱Problem 1）為膠體粒子不動，在流體裡受到外加電場作用；而第二個子問題（以下簡稱Problem 2）為膠體粒子在沒有外加電場的作用下相對於流體運動。故流場與電場之解，可視作此兩問題個別之解的總和。而值得注意的是，這兩個問題實際上沒有什麼太大的物理意義，只是當我們分別求解上述兩個子問題的時候，其解之合恰為原物理問題之解。無因次化之電泳速度之定義為 * z * * m E U U = ，代表無因次化的終端速度值，為無因次化外加電場值，其物理意義為單位電場下粒子移動之速度。而根據 O’Brien and White[9]，當我們將原系統分解成兩個子問題時，在 Problem 1 中，

液滴所得之總力若稱之為，其值和液滴所受到的外加電場值成正比： * U * z E 1 F E*_z F₁ =f₁'E*_z (2-114) 其中與無關；而在Problem 2 中，膠體粒子所得之總力稱之為，與粒子的運動速度成正比： ' 1 f E*_z F₂ * U * ' 2 2 f U F = (2-115) 其中與無關。當液滴達穩態時，由力平衡可知力總和為零，即，則可得液滴之無因次化電泳速度為 : ' 2 f U* F₁+F₂ =0 ' 2 ' 1 * z * * m f f E U U = = (2-116) 由於與f₁' E*_z無關、f₂' 與U*無關，為了計算方便，我們個別取U* =1以及E*_z =1，

(39)

如此一來、分別代表Problem 1、Problem 2 的合力。換言之，只要我們能計算出兩個子問題的總力，便可經由(2-116)式求出液滴的電泳速度。 1 ' 1 F f = f₂' =F₂ 既然電泳速度之計算為計算施於粒子之合力為零，分析本系統後，我們發現施於粒子上之力可以分成電力FEz和拖曳力FDz，其表示如下 (2-117)

∫

π δ φ −∇ σ π = 0 Ez 2 ( ) s F ∫ δ ∂ φ ∂ ρ ω π − ∫ δ ω ψ ∂ ∂ ω µπ = π π 0 2 0 2 2 3 Dz _n E s _s s F (2-118) 其中σ為液滴表面電荷密度，定義為ε(−∇φ)⋅i_n；而δ 代表液滴表面單位面積，s 表示成δS=2πωδs，其詳細推導過程可以參考Appendix C。 ξ ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ ξ ∂ φ ∂ ξ η − η ∂ φ ∂ ξ η − η ∂ φ ∂ ξ − π

=

_∫

π (1 cosh cos ) sinh sin d

x sin 2 F * * 0 * Ez (2-119)

∫

π π ξ ξ ∂ φ ∂ ξφ π κ + ξ ω ψ η ∂ ∂ ω − π = 0 * * 1 2 2 2 * 2 0 2 * * 2 3 * * Dz sin d x c ) a ( d E F (2-120) 所以，在兩個子問題裡，我們只要將個別子問題裡液滴表面所受合力加總，包括電力及拖曳力的作用，換句話說：

(

F F

)

f ' ,

(

i 1,2 F_zi* = _Ez* + _Dz* _problem_i = _i =

)

(2-121) 而液滴的電泳速度便可將(2-119)到(2-121)式代入(2-116)式加以計算，可以整理成如(2-122)式：

(

)

(

Ez* Dz*

)

problem2 1 problem * Dz * Ez 2 1 * z * * m F F F F ' f ' f E U U + + = − = = (2-112)

(40)

第

三章

數值方法

由於電泳現象之主控方程式為一組聯立之偏微分方程組，在數值方法的選擇上我們採用正交配位法(Orthogonal collocation method)或稱之為假性光譜法 (Pseudospectral method)來處理。正交配位法於近年來特別受到重視，此乃由於其可在較少的計算點數下，獲得較高的準確度；而其缺點為較不適合複雜的幾何形狀，這是因為使用正交配位法的時候，在計算空間為一矩形才有其優勢。而在計畫中我們採用雙極座標，此一座標為一正交系統(Orthogonal coordinate)，經過映射(Mapping)之後可以在計算區域為一矩形，故使用正交配位法相當適合。在本章中，將說明本計畫正交配位法的基本原理以及實際應用[52,53]。

(41)

3-1 正交配位法

配位法是光譜法的其中一種方法，故配位法才又稱為假性光譜法。而本文即是以Chebyshev polynomial為基底的正交配位法所做之研究，正交配位法可以說是加權剩餘法(Method of weighted residuals, MWR)的一種。使用MWR的觀念在於為嘗試函數(Trial function)和試驗函數(Test function)，試驗函數也稱為加權函數 (Weighting function)。以嘗試函數作為一未知待解函數的斷級數展開式的基本函數；而試驗函數則是藉著餘數(Residual)最小化來確定斷數展開式是否接近微分方程式。在配位法中試驗函數為平移delta函數(Shift delta function)δ(x-xi)至某一特別的配置點(Collocation points)xi上，在這逼近法中，必須在配置點上餘數為零。光譜法是將區間內之任意函數以嘗試函數展開，其 N 階近似函數可寫為 ∑ φ = = N 0 k k k N(x) a (x) u (3-1)

其中φk為此區間內的一組完全正交函數集(Complete set of orthogonal function)。假性光譜法中之(3-1)式，在N+1 個配置點上，滿足殘餘數(Resiual)為零。在此所謂的配置點即是一般的計算格點，不同的是它的位置有一定的決定方式，故本文稱之為配置點。所以： ∑ φ = = N 0 k k k j j) a (x ) x ( u (3-2) 其中為所選取的配置點，為配置點上的函數值。在區間上取 N+1 個配置點，由(3-1)至(3-2)式可得 j x u(x_j) x_j ∑ = = N 0 j j j N(x) u(x )q (x) u (3-3) ) x ( q_j 為對應於配置點的內插多項式(Interpolation polynomial)，它與選取的配置點的位置及所採用的嘗試函數直接相關。而關於配置點的選取與內插多項式的推

(42)

導如下：

在區間[-1,1]內之 N 次正交多項式(Orthogonal polynomial)，取其一次微分等

於零的根，加上 1,-1 兩點得到 N+1 個配置點，這種方法所選取的配置點稱為

Lobbato points。本研究以 Chebyshev 多項式作為 N 次正交多項式： )) x ( cos N cos( ) x ( T 1 N − = (3-4)

因此對應於 Chebyshev polynomial 的配置點，亦稱為 Chebyshev-Gauss-Lobatto 法。其表示式如下：- ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ π = N j cos x_j (3-5) 對應於所選取配置點的內插多項式為 ) x x ( N c ) 1 )( x ( T ) x 1 ( ) x ( q j 2 j 1 j ' N 2 j ₋ − − = + (3-6) 其中 ⎩ ⎨ ⎧ − ≤ ≤ = = 1 N j 1 , 1 N , 0 j , 2 c_j (3-7) I.正交配位法在[-1, +1]區間內的一維微分矩陣表示式 若想得到uN(x)的導函數，只需要對內插多項式作微分，(3-3)式對x作微分k 次，則產生： ∑ = = N 0 j k j k j k N k dx ) x ( q d ) x ( u dy ) x ( u d (3-8) 如果是想得到的是在配置點x=x_l上的微分值，則可寫成如下： N ,..., 0 l ; ) x ( u D dx ) x ( q d ) x ( u dx ) x ( u d N 0 j j ) k ( Nlj N 0 j k l j k j k l N k = ∑ = ∑ = = = _(3-9) 其中，D(_Nljk)即是一維的正交配位法微分矩陣，而k=1,2,3,4 分別代表一次、二次、

(43)

三次及四次的微分矩陣。由上述的推導得知，微分運算子_dyd 可以用一維的正交配位法之微分矩陣加以離散化(Discretization)： Φ ≡ (_N1) _lj N) (D ) u ( dy d (3-10) 此時(D(_N1))_lj =q(_j1)(x_l)，Φ=

[

u(x₀),u(x₁),_L_L,u(x_N)

]

T 根據(3-6)式我們可以得到一次微分矩陣(D(_N1))_lj的明顯形式(Explicit form)： ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎧ = = + − = = + ≤ = ≤ − − ≠ − − = + N j i 6 1 N 2 0 j i 6 1 2N 1 -N j i 1 ) x 1 ( 2 x j i x x ) 1 ( c c ) D ( 2 2 2 j j j i j i j i lj ) 1 ( N (3-11) II.正交配位法在[-1,+1]區間內的二維微分矩陣表示式

一個二維邊界值問題(Boundary value problem)(如本文中所運用的座標即為一二維軸對稱的系統)，可利用分離變數法(Separation of variables)。將二維的變數以Chebyshev 多項式的雙向斷級數展開(Double truncated series expansion)：

(3-12) ∑ ∑ = ψ = ψ = = N 1 n M 0 m nm n m NM(x,y) a T (x)T (y) 此時x及y方向可以選擇不同階數(Order)的多項式，各方向仍然選取前述所提的配置點，其分別為(xi,yj)： N , 0, i ) N cos( x i i = LL π = (3-13) M , 0, j ) M cos( y i i = LL π = (3-14)

(44)

應用正交配位法之後，可得： (3-15) ∑ ∑ = ψ = = N 1 i M 0 j nm n i m i i i NM(x ,y ) a T (x )T (y ) 上述的做法相當於將一二維變數以 x 方向以及 y 方向試驗函乘積的組合來近似，如此一來，假光譜法的二維型式可寫成： (3-16)

∑∑

= = = N 1 i M 0 j j i i i NM NM(x,y) f (x ,y )q (x)q (y) f 其中，x方向為N階近似、y方向為M階近似。qi(x)為x方向的內插多項式，qj(y) 為y方向的內插多項式： N , 0, i ) x x ( N c ) 1 )( x ( T ) x 1 ( ) x ( q i 2 i 1 i ' N 2 i = LL − − − = + (3-17) M , 0, j ) y y ( M c ) 1 )( y ( T ) y 1 ( ) x ( q j 2 j 1 j ' M 2 j = LL − − − = + (3-18) 當要做偏微分時，只需對某方向的內插多項式作偏微分。例如，需要對 x 方向作偏微分一次： ∑ ∑ ∂ ∂ = ∂ ∂ = = N 1 i M 0 j j i i i NM NM _q ₍_y₎ x ) x ( q ) y , x ( f x ) y , x ( f (3-19) 以類似一維微分矩陣形成的方法，二維偏微分矩陣： ∑ ∑ = ∑ ∑ ∂ ∂ ∂ ∂ = ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ ∂ ∂ ∂ = = = = = = + N 1 i M 0 j NM i i ) r , p ( klij NM N 1 i M 0 j r l j r p k i p i i NM y y , x x r pNM ) r p ( ) y , x ( f P y ) y ( q x ) x ( q ) y , x ( f y x ) y , x ( f l k (3-20) 在實際數值計算上，二維偏微分矩陣乃是藉著一維微分矩陣各分量的組合來計算，不需另外推導。 ) r , p ( klij NM P

(45)

3-2 空間映射 一般而言我們真正感興趣的物理空間並不是一定介於-1 至 1，而上述推導微分矩陣的區間我們稱之為計算空間卻必須介於-1 至 1 之間，因此吾人須將系統中的物理空間轉換為計算空間。另外基於液滴必須考慮內外兩區流體的影響，因此我們必須採用兩組不同的空間映射方法。通常來說對於一無限大的空間映射是屬於一非線性代數映射，但我們採用的雙球座標可以將無限大的空間座標值以有限值來表示。因此在本文中液滴外部的空間映射只需做一線性的代數映射。對於雙球座標而言外區的座標範圍為 π ≤ ξ ≤ η ≤ η ≤ , 0 0 ₀ ，故其空間轉換關係式為： 2 x 2 0 0 ₊η η = η (3-21) 2 y 2 π + π = ξ (3-22) 而液滴內部的座標範圍為η₀ ≤η<∞, 0≤ξ≤π，在這裡空間映射會遭遇到 ∞ → η 的困擾，但是實際上我們觀察(2-55)式到(2-84)式，我們在液滴內部實際上採用的邊界條件除了表面η=η₀以外，其餘都落在ξ=0、ξ=π、與η→∞三個地方，其中令人困擾的只有η→∞的邊界條件。而從(2-55)、(2-60)、(2-70)、與 (2-71)四條式子可以發現，因為我們假設液滴是一個非傳導體球，所以只有(2-70) 與(2-71)式之流場邊界條件會對內外兩區的結果產生嚴重的影響(因為流場的計算牽涉到兩區聯解，詳見3-3)，這裡我們思考的方向是，如果我們將網格開到十分逼近圖二裡的對稱軸，是否在數值上能達到可行的忽略呢？根據雙球座標與極座標的轉換關係，如下表一與表二所示，我們針對η=11η₀以及η=21η₀做為我們網格η方向的終點，對應到真實極座標空間：

(46)

0 /η η ξ z/c y/c 11 0 1-000033 0 11 0-174533 1-000033 5-8E-06 11 0-349066 1-000031 1-14E-05 11 0-523599 1-000029 1-67E-05 11 0-698132 1-000026 2-15E-05 11 0-785398 1-000024 2-36E-05 11 0-872665 1-000021 2-56E-05 11 0-959931 1-000019 2-74E-05 11 1-047198 1-000017 2-89E-05 11 1-22173 1-000011 3-14E-05 11 1-396263 1-000006 3-29E-05 11 1-570796 1 3-34E-05

Table 1 The values of polar coordinates (y/c, z/c) calculated from bipolar coordinates when

(

η ,ξ

)

η=11η0 0 /η η ξ z/c y/c 21 0 1 0 21 0-174533 1 2-63E-10 21 0-349066 1 5-19E-10 21 0-523599 1 7-58E-10 21 0-698132 1 9-75E-10 21 0-785398 1 1-07E-09

(47)

21 0-872665 1 1-16E-09 21 0-959931 1 1-24E-09 21 1-047198 1 1-31E-09 21 1-22173 1 1-43E-09 21 1-396263 1 1-49E-09 21 1-570796 1 1-52E-09

Table 2 The values of polar coordinates (y/c, z/c) calculated from bipolar coordinates when

(

η ,ξ

)

η=21η0 由以上的測試我們可以發現，因為兩座標系統的轉換幾乎是以指數進行，因此當η>η₀以後，換句話說在液滴內部的網格，(y/c, z/c)就變化得很快。由於我們希望y/c 越接近零越好(at y= ⇒0 z 對稱軸)，所以在本研究裡我們取η=21η0 作為我們內部的網格終點，也就是液滴內部的座標範圍為 π ≤ ξ ≤ η ≤ η ≤ η₀ 21 ₀, 0 ，其空間轉換關係式為： 0 0x 11 10η + η = η (3-23) 2 y 2 π + π = ξ (3-24) 所以，根據以上的關係式，可推導出程式裡液滴外部運算子的轉換關係式： η ∂ ∂ η = ∂ η ∂ η ∂ ∂ = ∂ ∂ 2 x x 0 _(3-25) n n n 0 n n 2 x ∂η ∂ ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ η = ∂ ∂ (3-26) ξ ∂ ∂ π = ∂ ξ ∂ ξ ∂ ∂ = ∂ ∂ 2 y y (3-27) n n n n n 2 y ∂ξ ∂ ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ π = ∂ ∂ (3-28)

(48)

另外，液滴內部運算子的轉換關係式也可推導如下： η ∂ ∂ η = ∂ η ∂ η ∂ ∂ = ∂ ∂ 0 10 x x (3-29)

(

)

n n_n 0 n n 10 x ∂η ∂ η = ∂ ∂ (3-30) ξ ∂ ∂ π = ∂ ξ ∂ ξ ∂ ∂ = ∂ ∂ 2 y y (3-31) n n n n n 2 y ∂ξ ∂ ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ π = ∂ ∂ (3-32) 故只要先處理計算空間中的運算子，之後再乘上轉換係數即可。

(49)

3-3 兩區聯解問題之處理 在液滴的系統中內外兩區的主控方程式不同，我們無法在一個單一的區間中計算這樣的問題。因此，我們必須將這個問題進行領域分割(Domain Decomposition)成兩個區間[52]，也就是我們在液滴系統敘述中所指的內區與外區。以下我們將以一個二階的微分方程式為例說明領域分割的作法。考慮一Helmholtz 方程式： 0 ) ( f 2_φ₋ _φ ₌ ∇ , a≤x≤c (3-33) 現在我們將這個問題分成兩個部分，其中內區的範圍是a≤x≤b，而外區的範圍是b≤x≤c，則原本的(3-33)式可分解為 ⎩ ⎨ ⎧ ≤ ≤ = φ − φ ∇ ≤ ≤ = φ − φ ∇ c x b 0 ) ( f b x a 0 ) ( f 2 2 (3-34) 我們若以光譜法展開上面的式子，則我們可以得到 0 F L_ijφ_j− _ijφ_j= (3-35) 或 0 ) F L ( _ij− _ij φ_j = (3-36) 其中是運算子的微分矩陣而代表函數f 的代數矩陣。而(3-36)式可以寫成矩陣的型態： ij L _∇2 ij F ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ − − − − − + + + + + + + + 1 N , 1 N 1 N , 1 N 1 , 1 N 1 , 1 N 1 N , 1 1 N , 1 12 12 11 11 F L F L F L F L F L L L M O M M O M L ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ φ φ +1 N 1 M M ₌ _(3-37) ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ 0 0 M M 之後我們將內區與外區的矩陣做組合，如圖三所示。虛線表示的區域為兩區域重疊之處，其中疊合的一點使外區的最後一點與內區的第一點相重疊，此舉使

(50)

件，我們將虛線中的微分矩陣，以連續條件的微分矩陣取代之。

Outside region

Inside region

0

(51)

3-4 牛頓-拉福生疊代法 在本研究中，平衡電位必須先單獨求解 Poisson-Boltzmann 方程式。而為了更有效率地求解此條式子，作者採用牛頓- 拉福生疊代法 (Newton-Raphson iteration scheme)。其為求解勘根問題f(x)=0時，將此一函數取其切線而得一收歛速度較為有效並且最著名的數值方法之一。首先簡介一下此一方法[65]：首先對於單一變數 F(x)=0 而言如圖四所示，取函數 F 在處的切線，交 x 軸於，再取函數在的切線交x 軸於。直到不再變化，則即為其解，其數學關係式為: r x 1 r x + xr+1 xr+2 xr xr r r r 1 r dx dF F x x ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ − = + (3-38)

Fig.4 Newton-Raphson method for one dependent variable

Xr Xr+1

(52)

對於n 條代數方程式 n ~ 1 i ; 0 ) x , , x x ( F ) x ( F_i = _i ₁_, ₂ _L _n = = (3-39) 來說，牛頓拉福生法在求解時是以下列的方式進行疊代。 (3-40) ) F(x )δ (xk k =− k ϕ 其中 n ~ 1 j ; n ~ 1 i ; ) x ( _k i ₌ ₌ _(3-41) x F ) x ( F ) (x j k ij k _∂ ∂ = = ϕ k n k 2 k 1 k (3-42) (3-43) 下標 k 代表第 k 次疊代。而只要能給予足夠好的猜值，便能以相當快的速度勘出根值，因為此法的收斂速度乃為二階。為如下所示: T )] x ( F , ), x ( F ), x ( F [ ) F(x = L T 0 n 20 10 0 k 1 k k x x ,x [x ,x , ,x ] δ = ₊ − = _L 若我們將(3-38)式寫為矩陣的型式，則 ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎢ ⎢ ∂ ∂ ∂ ₂ 2 F F F _⎢⎢x2(k+1) x2k ⎥ ⎢F2 ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ n n 2 n 1 n n 2 2 1 2 n 1 1 1 1 x F x F x F x x x x F x F x F L M O M M L L _⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ − − − − − + − + − + + nk ) 1 k ( n k ) 1 n ( ) 1 k ( 1 n k 3 ) 1 k ( 3 k 1 ) 1 k ( 1 x x x x x x x x M = ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ − − n 1 n 3 1 F F F F M (3-44) 我們可將(3-44)式之等號左側拆為兩項，其中一項含有第 k 次所得之值，一項則含有未知變數。然後將已知之項展開併入等號右側之矩陣，此時矩陣圍。所提之方法，我們便可求解 Poisson-Boltzmann 方程式，即(2-51) 。將(2-51)式以(3-38)式的型式展開，乃為： ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ k x 1 k x ₊

形式AX= 可直接以高斯消去法(Gauss elimination method)求解B xk+1，若其最大

相對誤差或最大絕對誤差大於容忍值，則再將x_k₊₁代入x_k重複計算直到誤差小於

容忍範

藉由上面式

膠體粒子對平面之電泳運動(3/3)

行政院國家科學委員會專題研究計畫 成果報告