光聚Rayleigh 粒子射流

牛振宇, 周哲瑋, 黃 凱, 張金松, 張建華, 王志亮

(1.上海大學上海市應用數學和力學研究所,上海200072；2.上海大學上海市力學在能源工程中的應用重點實驗室,上海200072；3.上海大學機電工程與自動化學院,上海200444；4.上海大學新型顯示技術及應用集成教育部重點實驗室,上海200444)

由顆粒物質組成的射流廣泛存在于多種自然現象及工程過程之中,如噴墨打印、噴丸工藝、粉煤氣化等.目前,在研究顆粒射流中逐步發展出多種產生顆粒射流的方法.Thoroddsen等[1]使用由直徑為1.34 cm的高速大球撞擊由直徑為80 μm的玻璃小球堆積而成的基底,證明顆粒物質仍然可以在微尺度下通過基底的變形-反彈產生高密度的射流.通過高速氣流載運顆粒也可以形成顆粒射流,在撞擊壁面后得到類似液膜的流動結構.Shi等[2]利用這種方法發現了顆粒流液膜的厚度與來流的速度、顆粒大小無關,液膜的形態演化方式與來流直徑和粒徑之比有關.

近年來,對微電子機械系統(microelectronic mechanical system,MEMS)的研究不斷增多.美國哈佛大學的研究人員就可以使用直徑為100 nm的粉末顆粒,通過3D打印技術成功制造出局部精度可達1 μm的鋰離子微電池[3]；麻省理工學院微系統技術實驗室基于電噴印技術發展出一套低成本加工方法,該方法利用包含納米顆粒的溶液制造出面積僅有0.03 mm2的微傳感器,其內部電路的金屬線僅10 μm寬,最小間隔50 μm[4],微納尺度下的增材制造(additive manuafcturing)技術也因此得到高速發展；Fuller等[5]利用壓電噴墨打印技術載運納米顆粒成功構造出多種微電子元件,展示出噴墨打印技術在MEMS領域的巨大應用潛力；作為噴墨打印技術的一種替代選擇,Huang等[6]利用氣溶膠噴印技術成功地在室溫下的A4紙質基底上完成了導電銀質薄膜的低成本快速打印.

對于微納米顆粒在射流中的輸運,顆粒射流研究著眼于直徑不小于1 μm的顆粒,流體環境影響往往可以忽略.而對于納米顆粒,除噴墨打印和氣溶膠噴印技術外,還存在多種應用方法.Xiao等[7]利用靜電紡絲技術將納米金屬顆粒分散嵌入纖維材料中,減少了顆粒團聚的出現并提高了材料性能；花銀群等[8]通過直流磁控濺射技術在單晶硅上制備鋁納米薄膜,并研究工藝參數對材料性能的影響,發現利用該技術制備的薄膜對特定的晶體結構顆粒有擇優取向；余潔意等[9]利用熱等離子法則成功制備出高純度、高催化效率的SiC納米顆粒.

由于納米顆粒的射流在輸運過程中往往需要各種液體、氣體或其他物質作為載體,流體作用不可忽略,使得在實際應用中受到限制.本工作期望通過借鑒光鑷技術的光場約束原理,構造出載體影響小或無需載體的微尺度顆粒射流.

光鑷技術研究源于20世紀70年代,Ashkin[10]在實驗室里偶然觀察到激光對溶液中膠體顆粒存在擾動現象.之后,Ashkin與其合作者進行了利用激光進行原子/分子、顆粒捕捉的理論實驗研究[11-12].直到1986年,Ashkin等[13]在實驗室中第一次成功利用單光束梯度力阱實現了對從25 nm到10 μm范圍膠體顆粒的穩定捕捉,才標志著光鑷技術的正式誕生.光鑷技術目前已被廣泛應用于染色體分離、微粒間相互作用的測量等微小尺度研究工作[14-15].

現有的光鑷技術文獻大多集中于對單個微粒操控的研究,對大量顆粒的群體行為的影響研究罕見于文獻專利之中.在研究有機發光二極管(organic light-emitting diode,OLED)真空蒸發鍍膜技術時,受光鑷光阱約束顆粒在光軸附近現象的啟發,本工作嘗試探討利用光鑷技術對顆粒流進行收束控制的可行性[16].對于球形微粒在電磁場中的受力計算方法,在Lorenz,Mie和Debye的工作中就已給出,即Lorenz-Mie理論,但該方法計算繁瑣.在光鑷技術誕生后,為了滿足光鑷設計的需求,先后發展出多種方法來簡化計算過程.通常根據激光波長對顆粒半徑的比值λ/a來區分應用于顆粒受力計算的方法.當λ/a＞20時,顆粒被稱為Rayleigh粒子,可以使用電偶極子模型得到十分精確的結果[17]；而當λ/a＜1/5時,顆粒被稱為Mie粒子,可以通過幾何光學進行近似計算[10].更為通用的GLMT(generalized Lorenz-Mie theory)方法則可用于各種尺度的球形顆粒受力計算[18].

在先前的工作中通過理論分析與數值模擬的手段驗證了Mie粒子的計算及收束效果[19-20],其結果適用于微米及納米尺度的顆粒.本工作的重點在于納米尺度Rayleigh粒子在激光作用下的行為.

1 物理模型與模擬方法

1.1 物理模型

本工作使用Lammps軟件進行了數值模擬.Lammps(large-scale atomic/molecular massively parallel simulator)是美國Sandia國家實驗室開發的分子動力學軟件,可以對諸如原子、生物分子、顆粒流或其他粗粒化系統在不同力場與邊界條件下的運動進行建模模擬.同時該軟件架構內集成了消息傳遞接口(message passing interface,MPI)接口,簡化了并行計算的實現過程[21].可以從其官方網站(http：//lammps.sandia.gov)獲得該軟件源代碼,并根據需求進行修改.

最初的物理模型如圖1(a)所示,該模型來自于真空蒸鍍中常見的克努森盒-沉積襯板組合.令激光穿過顆粒流出孔,指向襯板或收集裝置,在真空環境中誘導分子團簇或顆粒的運動.

在初步嘗試從蒸發腔內一直模擬到外部的沉積襯板位置后,發現孔口前的流動模擬需要占用大量的計算資源.另外,雖然在模擬中孔口處克努森數(平均分子自由程與孔徑之比)最小為0.003,已經處于粘滯流狀態,但顆粒出流流量和速度分布形式仍然受孔長度、平均自由程、孔徑等參數影響,其流動形態類似分子流,導致影響模擬結果的變量過多[22].由于最終關注的是激光對顆粒的收束效果,故綜合以上考慮,問題可以進一步簡化為孔口顆粒源與外部的激光作用域的組合,簡化后模型如圖1(b)所示.圖中r為計數平面上一點距光軸的距離,ω0為激光的束腰半徑,h為計數平面在z向的坐標,z0為激光束腰的z向坐標.

本工作將選擇點源與面源來分別討論激光對納米顆粒流動的收束效果.

1.2 Rayleigh粒子所受激光作用力

根據Harada的工作,Rayleigh粒子受激光作用產生的散射力Fscat和梯度力Fgrad可以通過式(1)計算得到[17]：

式中：b z為z方向單位矢量,即激光傳播方向；波數k=2π/λ,λ為激光波長；a為顆粒半徑；m為顆粒對環境的相對折射率；n為顆粒所處環境介質的折射率；c為光速；I(x,y,z)為平均坡印廷矢量.

本工作選擇了常見的線偏振基模高斯光束作為激光源,其平均坡印廷矢量在圖1(b)中的坐標系中為

式中,P為激光器功率,ω0為激光的束腰半徑.

1.3 模型的數值模擬

本工作中設置模擬區域大小為100 μm×100 μm×81 μm.不同的顆粒源放置在如圖1(b)中的坐標原點處,向z軸正向發射顆粒.

對于模擬中的納米顆粒,本工作均等效為理想球體,未引入其具體形貌的影響.納米顆粒在近似真空的環境中運動,流體環境作用較小,運動過程中作用較大的主要是納米顆粒相互之間與納米顆粒-激光之間的作用.

本工作研究的是一個介觀尺度問題,而這一尺度常用的粒子化研究方法有多粒子碰撞動力學(multiparticle collision dynamics,MPC)、耗散粒子動力學(dissipative particle dynamics,DPD)、直接蒙托卡羅(direct simulation Monte Carlo,DSMC)方法等[23].DPD方法經常被用來研究大分子在流體中的運動,通過對流體和分子進行粗粒化建模并應用統計方法得到的作用勢,可以減少計算量,加速模擬過程.粗粒化的過程使得DPD方法的作用勢中自然包含有反映流體環境的粘性作用及熱運動的部分[24].MPC和DSMC方法與之類似,存在由流體環境引入的力項.這些方法通常針對的問題與本工作要模擬的目標差別較大.

參考文獻[25]可知,納米尺度均質球形顆粒之間的吸引力會在遠小于其直徑的表面間距內快速衰減,并且在極近距離與碰撞過程中其相互吸引力會由于顆粒的變形而難以計算,而排斥力整體表現為顆粒碰撞時的彈塑性變形.目前還沒有找到適合本問題的顆粒吸引力的描述方法.根據Hamaker[26]的工作估算本工作中顆粒間的相互吸引力,在顆粒表面間距小于直徑時吸引力才能接近激光作用力.對于納米顆粒間的排斥作用,考慮到在顆粒流研究中,由于著眼于尺度大于1 μm的顆粒,顆粒間由范德華力產生的相互吸引作用可以忽略,故顆粒碰撞過程按照擠壓變形進行計算,其中塑性變形、質量變化等影響通過碰撞耗散來體現[27].因此,本工作在模擬中使用了由Silbert和Brilliantov發展的顆粒流碰撞作用勢來模擬納米顆粒間的排斥作用[28-29],但未模擬納米顆粒之間較弱的相互吸引力.本工作對于納米顆粒間的相互作用的模擬并不完善,還需要進行進一步的研究討論.

此時,按照顆粒流碰撞作用勢(即式(3))計算球形顆粒i,j之間的相互作用：

式中,kn為顆粒間的法向彈性系數(球心連線方向),δnij為顆粒中心距離的變化量,meff=mimj/(mi+mj)為顆粒碰撞時的有效質量,γn為法向阻尼系數,vn為2顆粒相對速度在法向上的分量,kt為切向彈性系數,Δst為2球的切向位移向量,γs切向阻尼系數,vs為相對速度的切向分量.這里,彈性系數和阻尼系數需要人工設置.由于本工作尚不考慮碰撞耗散及顆粒轉動的影響,故在模擬中設置kt=0,kn=0.005 nN/μm,并關閉了阻尼系數的計算.

模擬時按照式(4)和(5),即Velocity-Verlet算法進行時間積分,得到顆粒運動軌跡[30]：

式中,r為顆粒位置矢量,v為顆粒速度矢量,F為顆粒受力,t為當前時刻,Δt為時間步長.

本工作選取了甘油這一常見的試驗物質來作為數值模擬時顆粒的參數來源,同時沿用了真空蒸鍍的環境假定.顆粒密度ρ=1.261×103kg/m3,相對折射率m=1.474 6,顆粒半徑a取0.03 μm.環境折射率n=1(真空).

2 不同顆粒源模擬與結果

2.1 點源模擬

本工作首先使用點源模型開始研究,初步驗證程序及激光收束顆粒的可行性.顆粒的初始速度需要滿足一定的速度分布,為此本工作參考了真空蒸鍍技術中對于各種蒸發源的敘述.通常在很多真空蒸鍍過程中,除在蒸發源出口附近顆粒在飛行過程中幾乎不發生碰撞[31].本工作選擇了容易得到速度分布的點源作為研究起點.點顆粒源出發時的顆粒速度概率密度分布具有各向同性.基于這個特性和統計力學的基本原理,給出在顆粒平均初動能為常數的情況下,其速度矢量v的概率密度函數f(v)在極坐標速度空間上所滿足的方程組(6),進一步推導得到由點源出發顆粒初始速率的概率密度分布(式(7)),由各向同性可以得到各速度分量滿足的概率密度分布[32]：

當參考Oring的關于真空蒸鍍沉積厚度的內容后,可以得到在如圖1(b)所示的計數平面上,無激光收束時距光軸距離r處穿過平面的顆粒數密度滿足的無量綱分布[31]：

式中,N為距光軸r處顆粒的穿越數密度,N0光軸處顆粒的穿越數密度.

30萬個顆粒從坐標原點位置在0.42 s時間范圍內隨機注入計算區域,顆粒速度在進入模擬區域時按照速度分布式(8)賦予.雖然點顆粒源發射顆粒是全向的,但本工作僅模擬了如圖1(b)所示的向上半球運動的顆粒.此外,模擬時需要設置發射距離間隔下限,以保證新注入顆粒不會在注入時刻就與已出發的顆粒發生碰撞,避免由于過近的距離導致碰撞計算時能量不守恒.為了減少模擬步數,允許時間步長Δt在模擬過程中,根據顆粒的最大速度與限定的單步最大移動距離在0.000 5 μs到0.050 0 μs區間上自動調整變化.

在模擬點顆粒源時激光器波長λ=632.8 nm,束腰位置z0=5 μm,束腰半徑ω0=5 μm,模擬點顆粒源時的平均初動能取=0.038 8 eV.不同功率激光作用下在h=10 μm處得到的穿越密度分布如圖2所示.

圖2 不同功率激光作用下在h=10 μm處理論與模擬結果密度分布比較Fig.2 Comparisons of theoretical and simulation result density distributions under different laser powers at h=10 μm

由圖2可以發現,隨著激光功率增大,顆粒的分布向激光光軸(r=0)附近集中；而在無激光作用(P=0)時,圖2(b)顯示出模擬結果回歸到如式(9)計算得出的無收束情況下應遵循的無量綱分布形式.因此,對于點源,激光確實可以在一定參數下有效地對顆粒在空間的運動軌跡進行收束,改善其在指定平面上的分布,同時模擬過程中少量碰撞后的顆粒運動未出現異常,初步驗證了該程序可以正常模擬顆粒的碰撞與運動.

2.2 光聚納米顆粒射流的理論分析和模擬

Lum[33]在1984年給出了的顆粒流的流體力學形式控制方程：

式中,ρ為局部的平均質量密度,T為顆粒溫度,Q為顆粒流的熱流矢量,P為顆粒流的壓力張量,Γ為單位體積的能量耗散速率,F為單位質量所受體力.

顆粒流的溫度T并不是指物質的熱力學溫度,而是指脈動速度的平均度量：

式中,f(v,r,t)為顆粒的速度分布函數,n(r,t)為顆粒的數密度,V(r,t)為速度場分布.

Lum構造的Q,P,Γ均為體積分數函數[33].(射流質量密度與顆粒密度之比)與顆粒溫度的

在之前的工作中,通過體力項F,即式(12)引入了激光的作用：

式中,T為電磁場的動量流密度,g為電磁場的動量密度[19].E,B分別為電場與磁場向量.y是一個單位張量.

然后,構造了如下的顆粒流來進行數值分析.在真空環境下,顆粒射流由半徑R0=3 mm的噴口以速率U=100 m/s流出.射流的出口體積分數為ν=0.1,最大體積分數νm=0.55.顆粒的直徑為100 nm,密度為1.5 g/cm3,折射率為1.3.顆粒間相互作用接近于彈性碰撞,彈性恢復系數η=0.95.射流的顆粒溫度為10 K.由于關注的是激光對顆粒射流的影響,引入無粘假設(Re,Reλ→∞)來簡化后續分析[19].

參考文獻[34],使用正則模態法將擾動變量ρ,(bVz),(bVr),bT寫為如式(15)的形式[34],其中Ω為時間增長因子,k為軸向波數.當Ω＜0時流動穩定,Ω=0時流動中型穩定,Ω＞0時流動隨時間增加會失穩：

代入方程組(13)中并利用無粘假定進行簡化,得到一個2階常微分方程：

利用譜配置法對新方程(16)再次進行變換,成為如式(17)的線性方程組,以便數值求解.方程組為一齊次行列式,其存在非平凡解的條件是矩陣G(Ω,k,···)的行列式為0[19]：

利用Mu¨ller法對方程組(17)進行了求解[35].圖3反映了計算得到的沿徑向收束力和顆粒溫度的變化趨勢.圖中徑向力始終指向光軸,并隨半徑增大而增大,使得顆粒向光軸聚攏,射流存在收束的趨勢.由黃凱[19]推導得到的的解析表達式還表明,如果激光功率增強,約束力就增大.而顆粒溫度的分布則顯示出顆粒在邊界附近的活躍程度快速降低,內部的變化相對較小.圖4是在不同電歐拉數Eu下,軸向波數k與增長因子Ωr的關系.電歐拉數Eu沿圖中箭頭方向減小,增長因子隨之增大.由圖中的情況變化可知,隨電歐拉數Eu減小(如Eu≤0.01)后,增長因子Ωr＞0,按照正則模態法的分析系統將失穩,不能維持穩定射流；而當電歐拉數Eu≥0.5時,系統保持穩定.而增大電歐拉數Eu,最直接的方法就是提高激光功率,這與本工作利用激光可以穩定射流的構想基本一致.而圖4還反映出激光對射流的長波與短波失穩的同時進行了抑制,這一點不同于常見的穩定性理論中典型表面力(即表面張力)對擾動波的選擇性,通常表面張力對射流表面長波有失穩作用,對短波則有抑制作用[19].

圖3 歸一化顯示徑向力與顆粒溫度Fig.3 Reduced axial forceand granular temperature

圖4 電歐拉數對流動穩定性的影響Fig.4 Electrical Euler number effects on the flow stability

2.3 納米顆粒射流的顆粒動力學模擬

在實際科研生產過程中,顆粒源可能為噴管、坩堝等組件,有一定出口面積.對于宏觀尺度下氣體射流或射流中顆粒在噴管出口處的速度分布,雖然Anderson等[36]和Nathanson等[37]做過一定的理論與實驗測量工作,但這些工作大多集中在噴管軸向上結合一維理論結果進行討論,缺少更完整的理論描述.而在微納尺度下對噴管進行模擬的過程中發現,幾何參數的調整會引起出口流動速度分布發生復雜的變化.綜合以上考慮,本工作通過在顆粒流進入模擬區域時疊加一個擾動速度分布的方法,來近似模擬發散的面源.在以下模擬結果中,擾動速度使用上述點顆粒源相同的概率密度分布函數.擾動速度分布參數設=0.004 eV,激光束腰位置z0調整為20 μm,其他參數與點顆粒源模擬相同.

為了減少模擬時的計算量,本工作參考了Moseler等[38]的顆粒注入方法.首先模擬了一群在直徑為1 μm,長度為20 μm的圓柱管道內隨機運動的30萬個顆粒作為液柱,同時固定顆粒間相對位置.將液柱以固定速度vz=0.1 m/s注入,對進入模擬區域的顆粒解除相對位置的固定并在vz基礎上疊加擾動速度.

圖5為不同功率(P=0,0.5,1.0 W)激光收束下,顆粒流某一時刻空間分布的x方向投影.圖中虛線為激光束的邊界.從圖中可以看到,顆粒的空間分布密度有明顯變化,隨著功率的增大向外散逸的顆粒密度逐漸降低,可見提高激光功率(電歐拉數)可以使系統形成穩定的射流,同時射流的直徑在模擬區域中沒有發生明顯變化,即沒有發生類似傳統射流理論中的失穩現象.這與黃凱得到的激光會同時抑制長波與短波失穩的結論相一致.

圖6為P=1.0 W時,50個相鄰編號顆粒的飛行軌跡.本工作在模擬結束后通過計算顆粒的2個相鄰時間步內飛行軌跡的角度改變量來判別是否發生碰撞.在如圖5(c)所示的P=1.0 W的模擬中,以5°為判別依據時顆粒在模擬域中的平均碰撞次數為12.2次,10°時為5.66次,平均碰撞頻率分別約為15.3×103Hz和7.08×103Hz(空氣中分子碰撞的頻率為109～1010Hz量級,即102量級速度與10-7量級分子自由程).從光束范圍外飛出的顆粒占總顆粒數的2.36%.

圖5 不同激光功率模擬下實時分布側視圖Fig.5 Run-time distribution side views from simulations with different laser powers

圖6 P=1.0 W時50個顆粒的運動軌跡Fig.6 50 particale traces when P=1.0 W

在h=80 μm的平面上,統計不同功率下的累計分布密度(見圖7),發現在本工作模擬的功率下光束范圍內的顆粒分布區別不大.使用如式(18)定義的均方半徑來表征整個平面上粒子位置的分散程度,以觀察激光的收束效果(見圖8).通過圖中擬合的曲線可以看出,隨著電歐拉數的線性增長,分散程度呈類似負指數下降的趨勢,

圖7 不同激光功率下h=80 μm處無量綱分布Fig.7 Reduced distributions with different laser powers at h=80 μm

圖8 不同Eu下在h=80 μm處的集中程度比較Fig.8 Concentration comparison with different Euat h=80 μm

3 結束語

本工作基于Lammps軟件實現了由點源與射流產生的納米顆粒流動疊加激光光場后的動力學模擬,同時應用理論分析的方法對激光作用下的納米顆粒射流的穩定性進行了分析.點源模擬結果與經典真空蒸發沉積理論分布的比對后發現,激光功率P的增強可以提高顆粒分布的集中程度,初步證明激光是可以有效影響顆粒的運動軌跡,從而改善納米尺度顆粒的沉積分布.

對于光聚納米顆粒射流穩定性的理論分析表明,提高電歐拉數Eu可以提高顆粒射流的穩定性形成穩定射流,同時發現激光對于顆粒射流的長波與短波穩定性同時進行了抑制.而顆粒射流的模擬結果表明,通過增大激光功率P來提高電歐拉數Eu確實可以更加有效地約束顆粒形成射流,并且由均方半徑隨激光功率的增大而減小表明,顆粒的集中程度也隨電歐拉數的提高而增強,表現出與理論預測一致的穩定性變化趨勢.