分子動力學(xué)模擬剪切速度對納米間隙中角鯊?fù)榻缑婊频挠绊?/h1>
2016-07-04 03:43:57潘伶高誠輝福州大學(xué)機(jī)械工程及自動化學(xué)院福建福州350108
化工學(xué)報(bào) 2016年4期 關(guān)鍵詞：聚合物
潘伶，高誠輝（福州大學(xué)機(jī)械工程及自動化學(xué)院，福建 福州 350108）
?
分子動力學(xué)模擬剪切速度對納米間隙中角鯊?fù)榻缑婊频挠绊?/p>
潘伶，高誠輝
（福州大學(xué)機(jī)械工程及自動化學(xué)院，福建 福州 350108）
摘要：采用聚合物一致性力場（PCFF），分別在7種剪切速度V和3種油膜厚度h下對納米間隙中潤滑劑角鯊?fù)檫M(jìn)行分子動力學(xué)模擬，分析固液界面的密度、分子和流速的分布，探究納米薄膜潤滑的潤滑機(jī)理和剪切速度對界面滑移的影響。結(jié)果表明，納米間隙中潤滑劑存在分層現(xiàn)象，各層間距相近，并非越遠(yuǎn)離固體壁面層間距越大，層間距約為角鯊?fù)榉肿訂蝹€(gè)C C鍵距離的3～4倍；隨著油膜厚度的減小，納米間隙中潤滑劑層狀分布越明顯，固化層密度越大；當(dāng)油膜厚度為3.44 nm時(shí)，固液界面滑移現(xiàn)象明顯，滑移長度b值隨著V先增大后減小，當(dāng)V 為22.8 m·s-1時(shí)，b達(dá)到最大值4.35 nm；根據(jù)模擬和計(jì)算結(jié)果，給出滑移長度與剪切速度的關(guān)系公式。
關(guān)鍵詞：分子模擬；納米結(jié)構(gòu)；聚合物；界面滑移；薄膜潤滑；角鯊?fù)?/p>
2015-06-01收到初稿，2015-09-10收到修改稿。
聯(lián)系人：高誠輝。第一作者：潘伶（1969—），女，博士，副教授。
Received date: 2015-06-01.
Foundation item: supported by the National Natural Science Foundation of China (51175085)，the Tribology Science Found of State Key Laboratory of Tribology (SKLTKF13A09) and the Natural Science Foundation of Fujian Province (2016J01226).
引 言
無滑移邊界條件是經(jīng)典流體力學(xué)和潤滑理論中的重要假設(shè)，廣泛應(yīng)用于流體的理論分析、實(shí)驗(yàn)和實(shí)際工程。然而，近年來的研究發(fā)現(xiàn)界面滑移在許多情況下有可能發(fā)生，特別是存在于微納米尺度下的流體。界面滑移的相關(guān)研究方法主要有理論數(shù)值分析、實(shí)驗(yàn)和分子動力學(xué)（MD）模擬。理論數(shù)值分析的物理模型主要有線性滑移長度模型[1]、極限剪應(yīng)力模型[2-3]和非線性滑移長度模型[4-5]。實(shí)驗(yàn)方法常用的有原子力顯微鏡（AFM）[6]和表面力儀（SFA）[7]等。AFM方法是通過探針測量擠壓流中流體阻力和膜厚的關(guān)系，再根據(jù)理論公式分析計(jì)算測量結(jié)果，而界面滑移的理論和數(shù)值分析還不完善。
目前，納米薄膜潤滑的研究仍處于初級階段，其潤滑機(jī)理、數(shù)學(xué)模型、微觀分子分布和油膜厚度等的理論分析、模擬計(jì)算和實(shí)驗(yàn)檢測都與工程實(shí)際有較大差距。對薄膜潤滑進(jìn)行深入研究并開發(fā)相應(yīng)的實(shí)際應(yīng)用技術(shù)有著重要的意義。隨著計(jì)算機(jī)計(jì)算能力、并行技術(shù)以及相關(guān)勢函數(shù)的發(fā)展，分子動力學(xué)模擬方法應(yīng)用日益廣泛[8-11]。采用MD模擬納米間隙潤滑劑，可以動態(tài)適時(shí)地顯示全油膜中分子和密度的分布規(guī)律，克服了實(shí)物實(shí)驗(yàn)無法測量高速加載時(shí)潤滑薄膜的密度和原子運(yùn)動的缺點(diǎn)。另外，實(shí)驗(yàn)的方法不能檢測新型的尚未研制出來的潤滑劑，而MD方法則可以靈活準(zhǔn)確地建立各種分子結(jié)構(gòu)模型。
雖然在過去的20年里，納米流體的界面滑移和薄膜潤滑機(jī)理是研究的熱點(diǎn)，但以往的研究，無論是實(shí)驗(yàn)、MD模擬還是數(shù)值分析的結(jié)論多數(shù)是建立在較為簡化的模型和推理上[12-20]。以往的MD研究多數(shù)采用簡單的單粒子硬球模型[5,21-23]，簡單的LJ勢函數(shù)描述分子間的相互作用[5,21-22,24]，或通過改變某一參數(shù)來模擬不同的濕潤性、固液相互作用勢和硬度等[22-23]對界面滑移的影響，而很少涉及液體的具體分子組成、結(jié)構(gòu)和更為準(zhǔn)確的分子力場。
本文采用聚合物一致性力場（PCFF）模擬研究納米間隙內(nèi)角鯊?fù)樵诠腆w壁面不同的剪切速度和油膜厚度下的分子、密度及流速分布，探究納米間隙中固液界面的界面滑移及潤滑機(jī)理。
1 模型的建立
圖1是納米間隙潤滑劑剪切模型，納米級角鯊?fù)闈櫥な芟抻趦善叫械墓腆w壁面。圖中球形表示原子，橙色球表示鋁Al原子，藍(lán)色球表示碳（C）原子，綠色球表示氫（H）原子。固體壁面采用面心立方fcc晶體鋁Al(001)晶面，因?yàn)锳l原子的Lennard-Jones勢函數(shù)中的能量參數(shù)e=3.32320 kcal·mol-1（1 cal=4.1868 J），距離參數(shù)s=2.9964?（1 ?=0.1 nm），其能量參數(shù)比Pt和Fe等原子小，固液之間的相互作用勢較小，使得常壓下界面滑移現(xiàn)象比較明顯，而且Al也是常用的金屬材料。在保證計(jì)算誤差小于3%的情況下，為減少計(jì)算機(jī)集群的計(jì)算量，取Al晶體表面尺寸為6.48 nm×4.05 nm，晶格常數(shù)為0.405 nm。上、下固體壁面施加同樣大小的反向速度V，V分別取7種介于3～30 m·s-1不同的速度。常壓下油膜厚度h分別取3.44、5.13 和7.02 nm。
圖1 不同剪切速度下納米間隙潤滑劑模型Fig.1 Model of confined thin film at different shear velocities
當(dāng)h=3.44 nm時(shí)，模型中總共有13408個(gè)原子，其中有3840個(gè)Al原子、3120個(gè)C原子和6448個(gè)H原子，即潤滑劑中共有104個(gè)角鯊?fù)榉肿?；?dāng)h=5.13 nm時(shí)，模型中總共有18100個(gè)原子，潤滑劑中共有155個(gè)角鯊?fù)榉肿?；?dāng)h=7.02 nm時(shí)，模型中總共有23344個(gè)原子，潤滑劑中共有212個(gè)角鯊?fù)榉肿印?/p>
常用于美容產(chǎn)品和醫(yī)藥保健品的角鯊?fù)?，還可以作為電接觸或機(jī)械接觸潤滑油的基礎(chǔ)油，無毒且生物降解性高，是良好的環(huán)境友好潤滑劑。角鯊?fù)榈姆肿邮綖镃30H62，分子結(jié)構(gòu)和模型如圖2所示。
2 模擬方法
采用大規(guī)模原子/分子并行模擬器（large-scale atomic/molecular massively parallel simulator，LAMMPS[25]）編程模擬納米間隙潤滑劑的剪切流動，分子力場采用聚合物一致性力場（polymerconsistent force field，PCFF），其力場參數(shù)是由大量的實(shí)驗(yàn)數(shù)據(jù)和精確的量子計(jì)算結(jié)果擬合得出，因此計(jì)算量大，模擬也更準(zhǔn)確，適用于計(jì)算聚合物、有機(jī)物和生物分子體系，并含有一些金屬原子的力場參數(shù)。 PCFF勢函數(shù)包含非鍵結(jié)作用項(xiàng)、鍵伸縮項(xiàng)、鍵角彎曲項(xiàng)、二面角扭曲項(xiàng)、離平面振動項(xiàng)和多種交叉作用項(xiàng)，其勢函數(shù)和各參數(shù)的選取詳見文獻(xiàn)[25-28]。
模型中油膜厚度方向即z方向采用自由邊界條件，其他兩個(gè)方向都設(shè)置成周期性邊界條件。由于界面滑移受上、下固體壁面的變形影響較大，故將固體壁面設(shè)置成彈性可變形，壁厚為1.21 nm。非鍵結(jié)作用項(xiàng)中的LJ勢作用項(xiàng)和庫侖作用項(xiàng)的截?cái)喟霃蕉既?2.5 ?，因?yàn)楫?dāng)油膜厚度為5.31 nm時(shí)，其計(jì)算結(jié)果與截?cái)喟霃饺?6 ?時(shí)非常接近，誤差小于0.5%，但在美國SUN公司的Blade60000刀片服務(wù)器上單節(jié)點(diǎn)的計(jì)算時(shí)間可減少近1/3。分子動力學(xué)模擬時(shí)間步長取1 fs。采用共軛梯度法，通過迭代調(diào)整原子的位置以減小模擬體系的能量，確定合理的初始位形，其中能量偏差取10-6kcal·mol-1，力偏差取10-8kcal·(mol·?)-1。先在微正則系綜（NVE）下對模擬體系進(jìn)行0.5 ns的弛豫以達(dá)到293 K恒溫條件下的平衡。然后，在上、下固體壁面分別施加7種x向的方向相反的速度V，使納米間隙中潤滑劑受到1.8 ns的剪切作用。
圖2 角鯊?fù)榉肿覨ig.2 Squalane molecule
圖3 h=7.02 nm、V=15 m·s-1時(shí)吸附層分子Fig.3 Distributions of squalane molecules adjacent to the wall at h=7.02 nm and V=15 m·s-1
根據(jù)模擬結(jié)果，采用簡單實(shí)用的線性滑移長度模型計(jì)算滑移長度。線性滑移長度模型[1]認(rèn)為固液交界處流體與固體的相對速度即滑移速度vs與流體在壁面的剪切速率g呈正比[29]
式中，b為滑移長度，表示固液界面滑移的程度；vx是流體沿x方向的流動速度。
3 結(jié)果與討論
3.1潤滑劑分子和密度的分布
圖3是h=7.02 nm、V=15 m·s-1時(shí)，上固體壁面和角鯊?fù)槲侥し肿釉诓煌羟兴矔r(shí)t的分布。由圖可見，開始剪切時(shí)即t=0時(shí)，有較多的潤滑劑分子接近垂直于固體表面分布，且跨越2～3層以上的潤滑劑聚集層。隨著剪切時(shí)間的增長，這些接近垂直分布的潤滑劑分子出現(xiàn)彎曲和扭轉(zhuǎn)，并與固體表面形成較小的角度或長分子鏈中的一部分與固體表面平行。
然而，以往的研究多數(shù)認(rèn)為流體分子在近壁面處主要是平行于壁面，距離壁面越遠(yuǎn)的分子越垂直于壁面[12,30-31]。以往的研究多采用Legendre多項(xiàng)式表示的方向參數(shù)S( z )
式中，z是粒子質(zhì)心的z坐標(biāo)，θ是鍵方向與z軸的夾角。當(dāng)鍵平行固體壁面時(shí)，q=90°，S( z )=-0.5；當(dāng)鍵垂直于固體壁面時(shí)，q=0°，S( z )= 1。鍵的方向是取兩相鄰成鍵原子之間的鍵的方向[12,30]，或取幾個(gè)連續(xù)成鍵原子中某兩個(gè)原子中心的連線方向[31]，或以分子中功能基團(tuán)定位，通過統(tǒng)計(jì)各區(qū)域的平均S( z )值來判斷分子的方向。方向參數(shù)S( z )只能表明分子中部分鍵的方向，適于研究短鏈分子的方向，而長鏈分子在運(yùn)動過程中普遍出現(xiàn)多次彎曲和扭轉(zhuǎn)等，故S( z )不宜用于判斷長鏈分子的方向。
當(dāng)剪切時(shí)間t=1.8 ns時(shí)，將模擬體系沿z方向的總尺寸L（圖1）分為200等分，以各原子中心所在位置為各原子的質(zhì)心，計(jì)算出各等分中潤滑劑的平均密度，并畫出密度分布曲線，如圖4(a)、(c)、(e)所示。由圖可見，3種油膜厚度下，由于受到固壁的約束，在固液交界處都有一層密度很大的潤滑劑層，即出現(xiàn)固化現(xiàn)象。油膜厚度越小固化層密度越大，固化層對界面滑移的影響也越大。
圖4 沿油膜厚度方向角鯊?fù)橘|(zhì)量密度和部分分子的分布Fig.4 Density and molecular distributions of lubricants in the z direction
當(dāng)h=3.44 nm時(shí)，納米間隙內(nèi)潤滑劑明顯為層狀分布并分為7層，如圖4(a)所示，圖中標(biāo)出了各層之間的相對距離，可以看出各層之間距離沒有明顯的規(guī)律，而以往其他學(xué)者的研究常認(rèn)為離壁面越遠(yuǎn)層間距越大[12-14]。實(shí)際上，在文獻(xiàn)[32-33]中，采用原子力顯微鏡（AFM）測量探針在趨近潤滑劑中的基片時(shí)所受的阻力，仔細(xì)觀察文中的圖就可以發(fā)現(xiàn)力曲線的層間距是不規(guī)則的，相同的現(xiàn)象在采用MD模擬的文獻(xiàn)[34]的圖中也可以觀察到。當(dāng)h=5.13 nm時(shí)，上、下固液交界處附近分別有明顯的3～4層潤滑劑，膜厚中部波動較小，如圖4(c)所示，同樣并非離壁面越遠(yuǎn)層間距越大；隨著油膜厚度的增加，潤滑劑分層現(xiàn)象越來越不明顯。當(dāng)h=7.02 nm時(shí)，只在上、下固液交界處附近有較明顯的3～4層潤滑劑，膜厚中部波動更小且無規(guī)則。另外，3種油膜厚度下，固體壁面附近的潤滑劑層間距基本相等，約為0.507 nm。角鯊?fù)榉肿訂蝹€(gè)鍵的距離約為0.154 nm，各潤滑劑聚集層中心間距離約為角鯊?fù)榉肿訂蝹€(gè)鍵距離的3～4倍。圖4(b)、(d)、(f)是用不同的顏色表示隨機(jī)抽取的少量分子，以清楚地顯示潤滑劑分子在納米間隙中的分布，由圖可見，同一分子常常跨越多層，而且出現(xiàn)彎曲和扭轉(zhuǎn)。
以往的多數(shù)MD研究采用單粒子硬球模型，結(jié)果顯示在垂直于油膜厚度方向（即xy平面），潤滑劑分子在固液交界處呈有序化排列，x向密度分布呈周期性變化[13,19]。然而，研究結(jié)果顯示，雖然納米薄膜潤滑中，z方向即沿油膜厚度方向潤滑劑分層現(xiàn)象明顯，但x方向沒有分層的現(xiàn)象。圖5是固液交界處xy平面上部分角鯊?fù)榈姆肿臃植?。圖6是角鯊?fù)樵趚方向密度分布曲線，可見曲線波動小，密度分布均勻，不存在周期性變化規(guī)律。
圖5 固液交界面角鯊?fù)椴糠址肿拥姆植糉ig.5 Molecular distribution of squalane at solid- liquid interface
圖6 x向角鯊?fù)橘|(zhì)量密度分布Fig.6 Density profile of squalane in x direction
3.2界面滑移
當(dāng)h=3.44 nm時(shí)，將整個(gè)模擬體系沿z方向分成34等分，圖7是當(dāng)V=10 m·s-1時(shí)，整個(gè)模擬體系沿x方向的速度分布。每個(gè)數(shù)據(jù)點(diǎn)的橫坐標(biāo)表示每等分中心的z坐標(biāo)在整個(gè)體系中的相對位置；縱坐標(biāo)表示每等分中原子x方向的平均速度。圖7中上、下兩條水平數(shù)據(jù)點(diǎn)分別表示上、下兩固體壁面的速度；中間的數(shù)據(jù)點(diǎn)是角鯊?fù)榱黧w的速度分布。由圖可見，固液交界處有明顯的相對速度即出現(xiàn)界面滑移，潤滑劑流速vx呈非線性分布。固體壁面的速度V與流體在固液交界處的速度vx之差即為滑移速度vs，圖中流體速度分布曲線在固液交界處的變化率除以L即為剪切速率g。
圖7 V=10 m·s-1時(shí)x向速度分布Fig.7 Velocity in x direction at V=10 m·s-1
圖8是h=3.44 nm、V=10 m·s-1時(shí)流體在3個(gè)方向的分段平均流速vx、vy和vz分布。與vx相比，vy和vz變化很小，特別是vz變化更小，說明在宏觀上流體屬于層流。圖8中的插圖是放大的vz分布圖，可以看出，在靠近固液交界處vz幾乎為0，而油膜中心附近vz稍大。
圖8 V=10 m·s-1時(shí)角鯊?fù)?個(gè)方向的速度分布Fig.8 Velocity in three directions at V=10 m·s-1
在固體壁面上分別施加7種不同的速度V，并取流體上半部分（z/L≥0.5）進(jìn)行分析，如圖9所示。由圖可見，當(dāng)V=3～30 m·s-1時(shí)，固液交界處的流體速度vx和速度變化率隨著V的增大而增大；當(dāng)V=25 m·s-1后，固液交界處的速度變化率明顯增大；同時(shí)，隨著V的增大，流體速度vx呈非線性分布越明顯。
圖9 不同V時(shí)納米間隙上半部分流體的速度分布Fig.9 Velocities of upper fluids at different shear velocities
圖10(a)、(b)中的離散數(shù)據(jù)點(diǎn)分別是根據(jù)圖9的數(shù)據(jù)計(jì)算出的固液界面的相對速度vs和剪切速率g，圖中曲線分別是對離散點(diǎn)進(jìn)行線性和3次多項(xiàng)式擬合后的曲線。圖10(c)中的離散數(shù)據(jù)點(diǎn)是根據(jù)圖10(a)、(b)中的離散數(shù)據(jù)點(diǎn)，由線性滑移長度模型的式（1）計(jì)算得出的滑移長度b，圖中曲線是二次多項(xiàng)式擬合曲線，其表達(dá)式為
式（3）方便地表達(dá)了滑移長度b與剪切速度V的關(guān)系。由圖10可以看出，隨著V的增大，vs值基本呈線性增大，剪切速率g呈非線性增大且在V＞20 m·s-1后增大更明顯，從而使得滑移長度b值隨著V先增大后減小。當(dāng)V=22.8 m·s-1時(shí)，b達(dá)最大值4.35 nm。
4 結(jié) 論
在PCFF力場下，采用MD方法模擬不同剪切速度下納米間隙角鯊?fù)榈募羟辛鲃?。分析了流體分子、密度和流速的分布，計(jì)算出固液界面的滑移長度，得出以下結(jié)論。
（1）固液交界處液體分子形態(tài)各異，長鏈分子中的部分平行于固體壁面，部分垂直于固壁，還有部分與固壁呈一定角度。方向參數(shù)S(z)不宜作為長鏈分子方向的判定依據(jù)。
（2）各種油膜厚度下，在固體壁面附近，都有一層密度很大的潤滑劑固化層。隨著油膜厚度的減小，納米間隙中潤滑劑層狀分布越明顯，固化層密度越大，固化層對界面滑移的影響也越大。層間距基本相等，并非離壁面越遠(yuǎn)層間距越大。層間距約為角鯊?fù)榉肿訂蝹€(gè)C C鍵距離的3～4倍。
圖10 各參數(shù)隨V的變化Fig.10 Factors as functions of shear velocity
（3）由于Al和squalane原子之間相互作用勢較小時(shí)，在常壓下就出現(xiàn)界面滑移現(xiàn)象。當(dāng)油膜厚度為3.44 nm時(shí)，隨著剪切速度V的增大，流速vx明顯呈非線性分布，固液交界處的相對速度vs基本呈線性增大，剪切速率g呈非線性增大且在V＞20 m·s-1后增大更明顯，從而使得滑移長度b值隨著V先增大后減小。根據(jù)模擬和計(jì)算結(jié)果，給出了滑移長度與剪切速度的關(guān)系公式。
本文的研究克服了實(shí)物實(shí)驗(yàn)無法測量納米間隙中流體的密度及分子運(yùn)動規(guī)律等的缺點(diǎn)，為采用分子動力學(xué)方法定量地、動態(tài)適時(shí)地預(yù)測納米間隙流體的密度和分子分布，進(jìn)一步探究納米薄膜潤滑機(jī)理、研究應(yīng)用界面滑移和開發(fā)新型潤滑劑提供了可靠的依據(jù)。
符號說明
b ——滑移長度，nm
h ——油膜厚度，nm
L ——模擬體系沿z向的總尺寸，nm
t ——剪切時(shí)間，ns
V ——剪切速度，m·s-1
vs——滑移速度，m·s-1vx，vy，vz——分別為流體在x、y和z方向的分速度，
m·s-1
g——剪切速率，ns-1
ε ——能量參數(shù)，kcal·mol-1
ρ ——密度，kg·m-3
σ ——距離參數(shù)，?
References
[1] OLDSTEIN S. Modern Developments in Fluid Dynamics [M]. Vol Ⅱ. Oxford: Clarendon Press，1957: 676-680.
[2] BEAVERS G S，JOSEPH D D. Boundary conditions at a naturally permeable wall [J]. Journal of Fluid Mechanics，1967，30 (1): 197-207. DOI: 10.1017/S0022112067001375
[3] MA G J，WU C W，ZHOU P. Wall slip and hydrodynamics of two-dimensional journal bearing [J]. Tribology International，2007，40: 1056-1066. DOI: 10.1016/j.triboint.2006.10.003.
[4] THOMPSON P A，TROIAN S M. A general boundary condition for liquid flow at solid surfaces [J]. Nature，1997，389: 360-362. DOI: 10.1038/38686.
[5]RIEZJEV N V，TROIAN S M. Molecular origin and dynamic behavior of slip in sheared polymer films [J]. Physical Review Letters，2004，92 (1): 018302. DOI: 10.1103/PhysRevLett.92.018302.
[6] JING D，BHUSHAN B. Boundary slip of superoleophilic，oleophobic，and superoleophobic surfaces immersed in deionized water，hexadecane，and ethylene glycol [J]. Langmuir，2013，29 (47): 14691-14700. DOI: 10.1021/la4030876.
[7] ESPINOSA-MARZAL R M，ARCIFA A，ROSSI A，et al. Microslips to “avalanches” in confined，molecular layers of ionic liquids [J]. The Journal of Physical Chemistry Letters，2013，5 (1): 179-184. DOI: 10.1021/jz402451v.
[8] 陳其樂，孔憲，盧滇楠，等. 外壁電荷性質(zhì)對雙壁碳納米管中水分子運(yùn)動行為的影響 [J]. 化工學(xué)報(bào)，2014，65 (1): 319-327. DOI: 10.3969/j.issn.0438-1157.2014.01.042. CHEN Q L，KONG X，LU D N，et al. Molecular simulation of outer surface charge on water transport through double-wall carbon nanotube [J]. CIESC Journal，2014，65 (1): 319-327. DOI: 10.3969/j.issn.0438-1157.2014.01.042.
[9] 張程賓，趙沐雯，陳永平，等. 流體密度對納通道內(nèi)流動滑移的影響 [J]. 化工學(xué)報(bào)，2012，63 (S1): 12-16. DOI: 10.3969/j.issn.0438-1157.2012.zl.003. ZHANG C B，ZHAO M W，CHEN Y P，et al. Effects of fluid density on velocity slip in nanochannels [J]. CIESC Journal，2012，63 (S1): 12-16. DOI: 10.3969/j.issn.0438-1157.2012.zl.003.
[10] 潘伶，高誠輝. 納米間隙潤滑劑季戊四醇四酯的壓縮性能分子動力學(xué)模擬 [J]. 機(jī)械工程學(xué)報(bào)，2015，51 (5): 76-82. DOI: 10.3901/JME.2015.05.076. PAN L，GAO C H. Molecular dynamics simulation on the compressibility of pentaerythritol tetra in nanogap [J]. Journal of Mechanical Engineering，2015，51 (5): 76-82. DOI: 10.3901/JME.2015.05.076.
[11] PAN L，GAO C H. Confined fluid density of a pentaerythritol tetraheptanoate lubricant investigated using molecular dynamics simulation [J]. Fluid Phase Equilibria，2015，385: 212-218. DOI: 10.1016/j.fluid.2014.11.014.
[12] TSIGE M，PATNAIK S S. An all-atom simulation study of the ordering of liquid squalane near a solid surface [J]. Chemical Physics Letters，2008，457 (4/5/6): 357-361. DOI: 10.1016/j.cplett.2008.04.026.
[13] KARNIADAKIS G，BESKOK A，ALURU N. Microflows and Nanoflows: Fundamentals and Simulation [M]. Berlin Heidelberg，New York: Springer，2004: 1-23.
[14] PERTSIN A J，GRUNZE M. Long-ranged solvation forces in a fluid with short-ranged interactions [J]. Journal of Chemical Physics，2003，118 (17): 8004-8009. DOI: 10.1063/1.1564051.
[15] LOI S，SUN G，F(xiàn)RANZ V，et al. Rupture of molecular thin films observed in atomic force microscopy (II): Experiment [J]. Physical Review E，2002，66 (3): 031602. DOI: 10.1103/PhysRevE.66.031602.
[16] FRANZ V，BUTT H-J. Confined liquids:? solvation forces in liquid alcohols between solid surfaces [J]. The Journal of Physical Chemistry B，2002，106 (7): 1703-1708. DOI: 10.1021/jp012541w.
[17] 陳天星. 利用原子力顯微鏡對近壁面受限液體性質(zhì)的研究 [D].北京: 清華大學(xué)，2011. CHEN T X. Study on properties of the confined liquids at solid surface with AFM[D]. Beijing: Tsinghua University，2011.
[18] SIVEBAEK I M，SAMOILOV V N，PERSSON B N J. Squeezing molecularly thin alkane lubrication films: layering transitions and wear [J]. Tribology Letters，2004，16 (3): 195-200. DOI: 10.1023/B:TRIL.0000009730.31175.82.
[19] 胡元中，王慧，郭炎. 超薄油膜潤滑的分子動力學(xué)模擬 [J]. 摩擦學(xué)學(xué)報(bào)，1995，15 (2): 138-144. DOI: 10.16078/j.tribology.1995.02.007. HU Y Z，WANG H，GUO Y. Molecular dynamics simulation of ultra thin film lubrication (I): Rigid molecule model [J]. Tribology，1995，15 (2): 138-144. DOI: 10.16078/j.tribology.1995.02.007.
[20] 王慧，胡元中，郭炎. 超薄潤滑膜界面滑移現(xiàn)象的分子動力學(xué)研究 [J]. 清華大學(xué)學(xué)報(bào) (自然科學(xué)版)，2000，40 (4): 107-110. WANG H，HU Y Z，GUO Y. Molecular dynamics study of interfacial slip behavior of ultrathin lubricating films [J]. J. Tsinghua Univ. (Sci. & Tech.)，2000，40 (4): 107-110.
[21] NAGAYAMA G，CHENG P. Effects of interface wettability on microscale flow by molecular dynamic simulation [J]. International Journal of Heat and Mass Transfer，2004，47: 501-513. DOI: 10.1016/j.ijheatmasstransfer.2003.07.013.
[22] ASPROULIS N，DRIKAKIS D. Boundary slip dependency on surface stiffness [J]. Physical Review E，2010，81 (6): 061503. DOI: 10.1103/PhysRevE.81.061503.
[23] ZHANG H，ZHANG Z，YE H. Molecular dynamics-based prediction of boundary slip of fluids in nanochannels [J]. Microfluidics and Nanofluidics，2012，12 (1-4): 107-115. DOI: 10.1007/s10404-011-0853-y.
[24] NOORIAN H，TOGHRAIE D，AZIMIAN A R. Molecular dynamics simulation of poiseuille flow in a rough nano channel with checker surface roughnesses geometry [J]. Heat and Mass Transfer，2014，50 (1): 105-113. DOI: 10.1007/s00231-013-1232-x.
[25] PLIMPTON S. LAMMPS molecular dynamics simulator[EB/OL]. [2015-5-15]. http://lammps.sandia.gov.
[26] PLIMPTON S. Fast parallel algorithms for short-range molecular dynamics [J]. Journal of Computational Physics，1995，117 (1): 1-19. DOI: 10.1006/jcph.1995.1039.
[27] SUN H. Ab initio calculations and force field development for computer simulation of polysilanes [J]. Macromolecules，1995，28 (3): 701-712. DOI: 10.1021/ma00107a006.
[28] SUN H. COMPASS:? an ab initio force-field optimized for condensed-phase applications overview with details on alkane and benzene compounds [J]. Journal of Physical Chemistry B，1998，102 (38): 7338-7364. DOI: 10.1021/jp980939v.
[29] NAVIER C L M H. Memoire sur les du movement des fluids [J]. Mem l'Acad Roy Sci l’Inst France，1823，6: 389-440.
[30] VADAKKEPATT A，DONG Y，LICHTER S，et al. Effect of molecular structure on liquid slip [J]. Physical Review E，2011，84 (6): 066311. DOI: 10.1103/PhysRevE.84.066311.
[31] MENDONCA A C F，F(xiàn)OMIN Y D，MALFREYT P，et al. Novel ionic lubricants for amorphous carbon surfaces: molecular modeling of the structure and friction [J]. Soft Matter，2013，9 (44): 10606-10616. DOI: 10.1039/C3SM51689J.
[32] LIM R，O’SHEA S J. Solvation forces in branched molecular liquids [J]. Physical Review Letters，2002，88 (24): 246101. DOI: 10.1103/PhysRevLett.88.246101.
[33] LIM R，LI S F Y，O’SHEA S J. Solvation forces using sample-modulation atomic force microscopy [J]. Langmuir，2002，18 (16): 6116-6124. DOI: 10.1021/la011789+.
[34] ZHENG X，ZHU H，KOSASIH B，et al. A molecular dynamics simulation of boundary lubrication: the effect of n-alkanes chain length and normal load [J]. Wear，2013，301 (1/2): 62-69. DOI: 10.1016/j.wear.2013.01.052.
Molecular dynamics simulation of boundary slip in nanogap: effect of shear velocity
PAN Ling，GAO Chenghui
(School of Mechanical Engineering and Automation，F(xiàn)uzhou University，F(xiàn)uzhou 350108，F(xiàn)ujian，China)
Abstract:Molecular dynamics (MD) simulations using the polymer consistent force field (PCFF) were adopted to investigate the density，molecular and velocity distributions of lubricant squalane in nanogap at 293 K，three different film thicknesses and a wide range of shear velocities. The lubrication mechanism and boundary slip were analyzed. The results showed that the lubricant atoms tended to form layers parallel to the confining wall. The distances between the layers of lubricant atoms were irregular rather than broadening far away from the walls and were about three to four times the length of C C bond in the squalane. The tendency of lubricant atoms to form layers and the density of solid-like layer increased with decreasing film thickness. It was clearly to find the boundary slip at the solid-liquid interface from the velocity profile. The slip lengths increased with increasing velocity of substrates at the beginning，and then decreased. When the film thickness was 3.44 nm，the maximum slip length was 4.35 nm at the substrate velocity of 22.8 m·s-1. According to the simulations，the relationship between the slip length and the shear velocity was given.
Key words:molecular simulation; nanostructure; polymers; boundary slip; thin film lubrication; squalane
DOI：10.11949/j.issn.0438-1157.20150755
中圖分類號：TH 117.2
文獻(xiàn)標(biāo)志碼：A
文章編號：0438—1157（2016）04—1440—08
基金項(xiàng)目：國家自然科學(xué)基金項(xiàng)目（51175085）；清華大學(xué)摩擦學(xué)國家重點(diǎn)實(shí)驗(yàn)室開放基金項(xiàng)目（SKLTKF13A09）；福建省自然科學(xué)基金項(xiàng)目（2016J01226）。
Corresponding author:Prof. GAO Chenghui，gch@fzu.edu.cn
猜你喜歡
注聚污水處理技術(shù)調(diào)研報(bào)告中國科技博覽(2017年5期)2017-04-25 23:36:47 聚合物改善水泥石滲透性評價(jià)方法研究價(jià)值工程(2017年11期)2017-04-18 17:58:11 高溫高鹽油藏用化學(xué)驅(qū)油劑的研究山東工業(yè)技術(shù)(2016年19期)2016-10-29 16:12:54 電化學(xué)法血糖儀電極的改性專利技術(shù)綜述中國科技博覽(2016年18期)2016-10-19 07:34:50 環(huán)保型抗高溫抗鹽聚合物降粘劑的合成及評價(jià)科技傳播(2016年10期)2016-07-15 23:11:56 采油生產(chǎn)過程中聚合物分散技術(shù)的應(yīng)用分析中國科技博覽(2016年13期)2016-07-13 01:53:14 淺談聚合物驅(qū)油區(qū)塊轉(zhuǎn)油站加藥措施中國科技博覽(2016年10期)2016-04-29 04:18:00 聚合物驅(qū)站場施工安全管理措施分析中國科技博覽(2016年1期)2016-04-25 12:04:06 復(fù)合型高分子導(dǎo)電材料淺析科技與創(chuàng)新(2016年5期)2016-03-17 02:33:51 磷酸鎂水泥的聚合物改性研究湖南大學(xué)學(xué)報(bào)·自然科學(xué)版(2014年7期)2014-11-28 17:33:37
潘伶，高誠輝（福州大學(xué)機(jī)械工程及自動化學(xué)院，福建 福州 350108）
?
分子動力學(xué)模擬剪切速度對納米間隙中角鯊?fù)榻缑婊频挠绊?/p>
潘伶，高誠輝
（福州大學(xué)機(jī)械工程及自動化學(xué)院，福建 福州 350108）
摘要：采用聚合物一致性力場（PCFF），分別在7種剪切速度V和3種油膜厚度h下對納米間隙中潤滑劑角鯊?fù)檫M(jìn)行分子動力學(xué)模擬，分析固液界面的密度、分子和流速的分布，探究納米薄膜潤滑的潤滑機(jī)理和剪切速度對界面滑移的影響。結(jié)果表明，納米間隙中潤滑劑存在分層現(xiàn)象，各層間距相近，并非越遠(yuǎn)離固體壁面層間距越大，層間距約為角鯊?fù)榉肿訂蝹€(gè)C C鍵距離的3～4倍；隨著油膜厚度的減小，納米間隙中潤滑劑層狀分布越明顯，固化層密度越大；當(dāng)油膜厚度為3.44 nm時(shí)，固液界面滑移現(xiàn)象明顯，滑移長度b值隨著V先增大后減小，當(dāng)V 為22.8 m·s-1時(shí)，b達(dá)到最大值4.35 nm；根據(jù)模擬和計(jì)算結(jié)果，給出滑移長度與剪切速度的關(guān)系公式。
關(guān)鍵詞：分子模擬；納米結(jié)構(gòu)；聚合物；界面滑移；薄膜潤滑；角鯊?fù)?/p>
2015-06-01收到初稿，2015-09-10收到修改稿。
聯(lián)系人：高誠輝。第一作者：潘伶（1969—），女，博士，副教授。
Received date: 2015-06-01.
Foundation item: supported by the National Natural Science Foundation of China (51175085)，the Tribology Science Found of State Key Laboratory of Tribology (SKLTKF13A09) and the Natural Science Foundation of Fujian Province (2016J01226).
引 言
無滑移邊界條件是經(jīng)典流體力學(xué)和潤滑理論中的重要假設(shè)，廣泛應(yīng)用于流體的理論分析、實(shí)驗(yàn)和實(shí)際工程。然而，近年來的研究發(fā)現(xiàn)界面滑移在許多情況下有可能發(fā)生，特別是存在于微納米尺度下的流體。界面滑移的相關(guān)研究方法主要有理論數(shù)值分析、實(shí)驗(yàn)和分子動力學(xué)（MD）模擬。理論數(shù)值分析的物理模型主要有線性滑移長度模型[1]、極限剪應(yīng)力模型[2-3]和非線性滑移長度模型[4-5]。實(shí)驗(yàn)方法常用的有原子力顯微鏡（AFM）[6]和表面力儀（SFA）[7]等。AFM方法是通過探針測量擠壓流中流體阻力和膜厚的關(guān)系，再根據(jù)理論公式分析計(jì)算測量結(jié)果，而界面滑移的理論和數(shù)值分析還不完善。
目前，納米薄膜潤滑的研究仍處于初級階段，其潤滑機(jī)理、數(shù)學(xué)模型、微觀分子分布和油膜厚度等的理論分析、模擬計(jì)算和實(shí)驗(yàn)檢測都與工程實(shí)際有較大差距。對薄膜潤滑進(jìn)行深入研究并開發(fā)相應(yīng)的實(shí)際應(yīng)用技術(shù)有著重要的意義。隨著計(jì)算機(jī)計(jì)算能力、并行技術(shù)以及相關(guān)勢函數(shù)的發(fā)展，分子動力學(xué)模擬方法應(yīng)用日益廣泛[8-11]。采用MD模擬納米間隙潤滑劑，可以動態(tài)適時(shí)地顯示全油膜中分子和密度的分布規(guī)律，克服了實(shí)物實(shí)驗(yàn)無法測量高速加載時(shí)潤滑薄膜的密度和原子運(yùn)動的缺點(diǎn)。另外，實(shí)驗(yàn)的方法不能檢測新型的尚未研制出來的潤滑劑，而MD方法則可以靈活準(zhǔn)確地建立各種分子結(jié)構(gòu)模型。
雖然在過去的20年里，納米流體的界面滑移和薄膜潤滑機(jī)理是研究的熱點(diǎn)，但以往的研究，無論是實(shí)驗(yàn)、MD模擬還是數(shù)值分析的結(jié)論多數(shù)是建立在較為簡化的模型和推理上[12-20]。以往的MD研究多數(shù)采用簡單的單粒子硬球模型[5,21-23]，簡單的LJ勢函數(shù)描述分子間的相互作用[5,21-22,24]，或通過改變某一參數(shù)來模擬不同的濕潤性、固液相互作用勢和硬度等[22-23]對界面滑移的影響，而很少涉及液體的具體分子組成、結(jié)構(gòu)和更為準(zhǔn)確的分子力場。
本文采用聚合物一致性力場（PCFF）模擬研究納米間隙內(nèi)角鯊?fù)樵诠腆w壁面不同的剪切速度和油膜厚度下的分子、密度及流速分布，探究納米間隙中固液界面的界面滑移及潤滑機(jī)理。
1 模型的建立
圖1是納米間隙潤滑劑剪切模型，納米級角鯊?fù)闈櫥な芟抻趦善叫械墓腆w壁面。圖中球形表示原子，橙色球表示鋁Al原子，藍(lán)色球表示碳（C）原子，綠色球表示氫（H）原子。固體壁面采用面心立方fcc晶體鋁Al(001)晶面，因?yàn)锳l原子的Lennard-Jones勢函數(shù)中的能量參數(shù)e=3.32320 kcal·mol-1（1 cal=4.1868 J），距離參數(shù)s=2.9964?（1 ?=0.1 nm），其能量參數(shù)比Pt和Fe等原子小，固液之間的相互作用勢較小，使得常壓下界面滑移現(xiàn)象比較明顯，而且Al也是常用的金屬材料。在保證計(jì)算誤差小于3%的情況下，為減少計(jì)算機(jī)集群的計(jì)算量，取Al晶體表面尺寸為6.48 nm×4.05 nm，晶格常數(shù)為0.405 nm。上、下固體壁面施加同樣大小的反向速度V，V分別取7種介于3～30 m·s-1不同的速度。常壓下油膜厚度h分別取3.44、5.13 和7.02 nm。
圖1 不同剪切速度下納米間隙潤滑劑模型Fig.1 Model of confined thin film at different shear velocities
當(dāng)h=3.44 nm時(shí)，模型中總共有13408個(gè)原子，其中有3840個(gè)Al原子、3120個(gè)C原子和6448個(gè)H原子，即潤滑劑中共有104個(gè)角鯊?fù)榉肿?；?dāng)h=5.13 nm時(shí)，模型中總共有18100個(gè)原子，潤滑劑中共有155個(gè)角鯊?fù)榉肿?；?dāng)h=7.02 nm時(shí)，模型中總共有23344個(gè)原子，潤滑劑中共有212個(gè)角鯊?fù)榉肿印?/p>
常用于美容產(chǎn)品和醫(yī)藥保健品的角鯊?fù)?，還可以作為電接觸或機(jī)械接觸潤滑油的基礎(chǔ)油，無毒且生物降解性高，是良好的環(huán)境友好潤滑劑。角鯊?fù)榈姆肿邮綖镃30H62，分子結(jié)構(gòu)和模型如圖2所示。
2 模擬方法
采用大規(guī)模原子/分子并行模擬器（large-scale atomic/molecular massively parallel simulator，LAMMPS[25]）編程模擬納米間隙潤滑劑的剪切流動，分子力場采用聚合物一致性力場（polymerconsistent force field，PCFF），其力場參數(shù)是由大量的實(shí)驗(yàn)數(shù)據(jù)和精確的量子計(jì)算結(jié)果擬合得出，因此計(jì)算量大，模擬也更準(zhǔn)確，適用于計(jì)算聚合物、有機(jī)物和生物分子體系，并含有一些金屬原子的力場參數(shù)。 PCFF勢函數(shù)包含非鍵結(jié)作用項(xiàng)、鍵伸縮項(xiàng)、鍵角彎曲項(xiàng)、二面角扭曲項(xiàng)、離平面振動項(xiàng)和多種交叉作用項(xiàng)，其勢函數(shù)和各參數(shù)的選取詳見文獻(xiàn)[25-28]。
模型中油膜厚度方向即z方向采用自由邊界條件，其他兩個(gè)方向都設(shè)置成周期性邊界條件。由于界面滑移受上、下固體壁面的變形影響較大，故將固體壁面設(shè)置成彈性可變形，壁厚為1.21 nm。非鍵結(jié)作用項(xiàng)中的LJ勢作用項(xiàng)和庫侖作用項(xiàng)的截?cái)喟霃蕉既?2.5 ?，因?yàn)楫?dāng)油膜厚度為5.31 nm時(shí)，其計(jì)算結(jié)果與截?cái)喟霃饺?6 ?時(shí)非常接近，誤差小于0.5%，但在美國SUN公司的Blade60000刀片服務(wù)器上單節(jié)點(diǎn)的計(jì)算時(shí)間可減少近1/3。分子動力學(xué)模擬時(shí)間步長取1 fs。采用共軛梯度法，通過迭代調(diào)整原子的位置以減小模擬體系的能量，確定合理的初始位形，其中能量偏差取10-6kcal·mol-1，力偏差取10-8kcal·(mol·?)-1。先在微正則系綜（NVE）下對模擬體系進(jìn)行0.5 ns的弛豫以達(dá)到293 K恒溫條件下的平衡。然后，在上、下固體壁面分別施加7種x向的方向相反的速度V，使納米間隙中潤滑劑受到1.8 ns的剪切作用。
圖2 角鯊?fù)榉肿覨ig.2 Squalane molecule
圖3 h=7.02 nm、V=15 m·s-1時(shí)吸附層分子Fig.3 Distributions of squalane molecules adjacent to the wall at h=7.02 nm and V=15 m·s-1
根據(jù)模擬結(jié)果，采用簡單實(shí)用的線性滑移長度模型計(jì)算滑移長度。線性滑移長度模型[1]認(rèn)為固液交界處流體與固體的相對速度即滑移速度vs與流體在壁面的剪切速率g呈正比[29]
式中，b為滑移長度，表示固液界面滑移的程度；vx是流體沿x方向的流動速度。
3 結(jié)果與討論
3.1潤滑劑分子和密度的分布
圖3是h=7.02 nm、V=15 m·s-1時(shí)，上固體壁面和角鯊?fù)槲侥し肿釉诓煌羟兴矔r(shí)t的分布。由圖可見，開始剪切時(shí)即t=0時(shí)，有較多的潤滑劑分子接近垂直于固體表面分布，且跨越2～3層以上的潤滑劑聚集層。隨著剪切時(shí)間的增長，這些接近垂直分布的潤滑劑分子出現(xiàn)彎曲和扭轉(zhuǎn)，并與固體表面形成較小的角度或長分子鏈中的一部分與固體表面平行。
然而，以往的研究多數(shù)認(rèn)為流體分子在近壁面處主要是平行于壁面，距離壁面越遠(yuǎn)的分子越垂直于壁面[12,30-31]。以往的研究多采用Legendre多項(xiàng)式表示的方向參數(shù)S( z )
式中，z是粒子質(zhì)心的z坐標(biāo)，θ是鍵方向與z軸的夾角。當(dāng)鍵平行固體壁面時(shí)，q=90°，S( z )=-0.5；當(dāng)鍵垂直于固體壁面時(shí)，q=0°，S( z )= 1。鍵的方向是取兩相鄰成鍵原子之間的鍵的方向[12,30]，或取幾個(gè)連續(xù)成鍵原子中某兩個(gè)原子中心的連線方向[31]，或以分子中功能基團(tuán)定位，通過統(tǒng)計(jì)各區(qū)域的平均S( z )值來判斷分子的方向。方向參數(shù)S( z )只能表明分子中部分鍵的方向，適于研究短鏈分子的方向，而長鏈分子在運(yùn)動過程中普遍出現(xiàn)多次彎曲和扭轉(zhuǎn)等，故S( z )不宜用于判斷長鏈分子的方向。
當(dāng)剪切時(shí)間t=1.8 ns時(shí)，將模擬體系沿z方向的總尺寸L（圖1）分為200等分，以各原子中心所在位置為各原子的質(zhì)心，計(jì)算出各等分中潤滑劑的平均密度，并畫出密度分布曲線，如圖4(a)、(c)、(e)所示。由圖可見，3種油膜厚度下，由于受到固壁的約束，在固液交界處都有一層密度很大的潤滑劑層，即出現(xiàn)固化現(xiàn)象。油膜厚度越小固化層密度越大，固化層對界面滑移的影響也越大。
圖4 沿油膜厚度方向角鯊?fù)橘|(zhì)量密度和部分分子的分布Fig.4 Density and molecular distributions of lubricants in the z direction
當(dāng)h=3.44 nm時(shí)，納米間隙內(nèi)潤滑劑明顯為層狀分布并分為7層，如圖4(a)所示，圖中標(biāo)出了各層之間的相對距離，可以看出各層之間距離沒有明顯的規(guī)律，而以往其他學(xué)者的研究常認(rèn)為離壁面越遠(yuǎn)層間距越大[12-14]。實(shí)際上，在文獻(xiàn)[32-33]中，采用原子力顯微鏡（AFM）測量探針在趨近潤滑劑中的基片時(shí)所受的阻力，仔細(xì)觀察文中的圖就可以發(fā)現(xiàn)力曲線的層間距是不規(guī)則的，相同的現(xiàn)象在采用MD模擬的文獻(xiàn)[34]的圖中也可以觀察到。當(dāng)h=5.13 nm時(shí)，上、下固液交界處附近分別有明顯的3～4層潤滑劑，膜厚中部波動較小，如圖4(c)所示，同樣并非離壁面越遠(yuǎn)層間距越大；隨著油膜厚度的增加，潤滑劑分層現(xiàn)象越來越不明顯。當(dāng)h=7.02 nm時(shí)，只在上、下固液交界處附近有較明顯的3～4層潤滑劑，膜厚中部波動更小且無規(guī)則。另外，3種油膜厚度下，固體壁面附近的潤滑劑層間距基本相等，約為0.507 nm。角鯊?fù)榉肿訂蝹€(gè)鍵的距離約為0.154 nm，各潤滑劑聚集層中心間距離約為角鯊?fù)榉肿訂蝹€(gè)鍵距離的3～4倍。圖4(b)、(d)、(f)是用不同的顏色表示隨機(jī)抽取的少量分子，以清楚地顯示潤滑劑分子在納米間隙中的分布，由圖可見，同一分子常常跨越多層，而且出現(xiàn)彎曲和扭轉(zhuǎn)。
以往的多數(shù)MD研究采用單粒子硬球模型，結(jié)果顯示在垂直于油膜厚度方向（即xy平面），潤滑劑分子在固液交界處呈有序化排列，x向密度分布呈周期性變化[13,19]。然而，研究結(jié)果顯示，雖然納米薄膜潤滑中，z方向即沿油膜厚度方向潤滑劑分層現(xiàn)象明顯，但x方向沒有分層的現(xiàn)象。圖5是固液交界處xy平面上部分角鯊?fù)榈姆肿臃植?。圖6是角鯊?fù)樵趚方向密度分布曲線，可見曲線波動小，密度分布均勻，不存在周期性變化規(guī)律。
圖5 固液交界面角鯊?fù)椴糠址肿拥姆植糉ig.5 Molecular distribution of squalane at solid- liquid interface
圖6 x向角鯊?fù)橘|(zhì)量密度分布Fig.6 Density profile of squalane in x direction
3.2界面滑移
當(dāng)h=3.44 nm時(shí)，將整個(gè)模擬體系沿z方向分成34等分，圖7是當(dāng)V=10 m·s-1時(shí)，整個(gè)模擬體系沿x方向的速度分布。每個(gè)數(shù)據(jù)點(diǎn)的橫坐標(biāo)表示每等分中心的z坐標(biāo)在整個(gè)體系中的相對位置；縱坐標(biāo)表示每等分中原子x方向的平均速度。圖7中上、下兩條水平數(shù)據(jù)點(diǎn)分別表示上、下兩固體壁面的速度；中間的數(shù)據(jù)點(diǎn)是角鯊?fù)榱黧w的速度分布。由圖可見，固液交界處有明顯的相對速度即出現(xiàn)界面滑移，潤滑劑流速vx呈非線性分布。固體壁面的速度V與流體在固液交界處的速度vx之差即為滑移速度vs，圖中流體速度分布曲線在固液交界處的變化率除以L即為剪切速率g。
圖7 V=10 m·s-1時(shí)x向速度分布Fig.7 Velocity in x direction at V=10 m·s-1
圖8是h=3.44 nm、V=10 m·s-1時(shí)流體在3個(gè)方向的分段平均流速vx、vy和vz分布。與vx相比，vy和vz變化很小，特別是vz變化更小，說明在宏觀上流體屬于層流。圖8中的插圖是放大的vz分布圖，可以看出，在靠近固液交界處vz幾乎為0，而油膜中心附近vz稍大。
圖8 V=10 m·s-1時(shí)角鯊?fù)?個(gè)方向的速度分布Fig.8 Velocity in three directions at V=10 m·s-1
在固體壁面上分別施加7種不同的速度V，并取流體上半部分（z/L≥0.5）進(jìn)行分析，如圖9所示。由圖可見，當(dāng)V=3～30 m·s-1時(shí)，固液交界處的流體速度vx和速度變化率隨著V的增大而增大；當(dāng)V=25 m·s-1后，固液交界處的速度變化率明顯增大；同時(shí)，隨著V的增大，流體速度vx呈非線性分布越明顯。
圖9 不同V時(shí)納米間隙上半部分流體的速度分布Fig.9 Velocities of upper fluids at different shear velocities
圖10(a)、(b)中的離散數(shù)據(jù)點(diǎn)分別是根據(jù)圖9的數(shù)據(jù)計(jì)算出的固液界面的相對速度vs和剪切速率g，圖中曲線分別是對離散點(diǎn)進(jìn)行線性和3次多項(xiàng)式擬合后的曲線。圖10(c)中的離散數(shù)據(jù)點(diǎn)是根據(jù)圖10(a)、(b)中的離散數(shù)據(jù)點(diǎn)，由線性滑移長度模型的式（1）計(jì)算得出的滑移長度b，圖中曲線是二次多項(xiàng)式擬合曲線，其表達(dá)式為
式（3）方便地表達(dá)了滑移長度b與剪切速度V的關(guān)系。由圖10可以看出，隨著V的增大，vs值基本呈線性增大，剪切速率g呈非線性增大且在V＞20 m·s-1后增大更明顯，從而使得滑移長度b值隨著V先增大后減小。當(dāng)V=22.8 m·s-1時(shí)，b達(dá)最大值4.35 nm。
4 結(jié) 論
在PCFF力場下，采用MD方法模擬不同剪切速度下納米間隙角鯊?fù)榈募羟辛鲃?。分析了流體分子、密度和流速的分布，計(jì)算出固液界面的滑移長度，得出以下結(jié)論。
（1）固液交界處液體分子形態(tài)各異，長鏈分子中的部分平行于固體壁面，部分垂直于固壁，還有部分與固壁呈一定角度。方向參數(shù)S(z)不宜作為長鏈分子方向的判定依據(jù)。
（2）各種油膜厚度下，在固體壁面附近，都有一層密度很大的潤滑劑固化層。隨著油膜厚度的減小，納米間隙中潤滑劑層狀分布越明顯，固化層密度越大，固化層對界面滑移的影響也越大。層間距基本相等，并非離壁面越遠(yuǎn)層間距越大。層間距約為角鯊?fù)榉肿訂蝹€(gè)C C鍵距離的3～4倍。
圖10 各參數(shù)隨V的變化Fig.10 Factors as functions of shear velocity
（3）由于Al和squalane原子之間相互作用勢較小時(shí)，在常壓下就出現(xiàn)界面滑移現(xiàn)象。當(dāng)油膜厚度為3.44 nm時(shí)，隨著剪切速度V的增大，流速vx明顯呈非線性分布，固液交界處的相對速度vs基本呈線性增大，剪切速率g呈非線性增大且在V＞20 m·s-1后增大更明顯，從而使得滑移長度b值隨著V先增大后減小。根據(jù)模擬和計(jì)算結(jié)果，給出了滑移長度與剪切速度的關(guān)系公式。
本文的研究克服了實(shí)物實(shí)驗(yàn)無法測量納米間隙中流體的密度及分子運(yùn)動規(guī)律等的缺點(diǎn)，為采用分子動力學(xué)方法定量地、動態(tài)適時(shí)地預(yù)測納米間隙流體的密度和分子分布，進(jìn)一步探究納米薄膜潤滑機(jī)理、研究應(yīng)用界面滑移和開發(fā)新型潤滑劑提供了可靠的依據(jù)。
符號說明
b ——滑移長度，nm
h ——油膜厚度，nm
L ——模擬體系沿z向的總尺寸，nm
t ——剪切時(shí)間，ns
V ——剪切速度，m·s-1
vs——滑移速度，m·s-1vx，vy，vz——分別為流體在x、y和z方向的分速度，
m·s-1
g——剪切速率，ns-1
ε ——能量參數(shù)，kcal·mol-1
ρ ——密度，kg·m-3
σ ——距離參數(shù)，?
References
[1] OLDSTEIN S. Modern Developments in Fluid Dynamics [M]. Vol Ⅱ. Oxford: Clarendon Press，1957: 676-680.
[2] BEAVERS G S，JOSEPH D D. Boundary conditions at a naturally permeable wall [J]. Journal of Fluid Mechanics，1967，30 (1): 197-207. DOI: 10.1017/S0022112067001375
[3] MA G J，WU C W，ZHOU P. Wall slip and hydrodynamics of two-dimensional journal bearing [J]. Tribology International，2007，40: 1056-1066. DOI: 10.1016/j.triboint.2006.10.003.
[4] THOMPSON P A，TROIAN S M. A general boundary condition for liquid flow at solid surfaces [J]. Nature，1997，389: 360-362. DOI: 10.1038/38686.
[5]RIEZJEV N V，TROIAN S M. Molecular origin and dynamic behavior of slip in sheared polymer films [J]. Physical Review Letters，2004，92 (1): 018302. DOI: 10.1103/PhysRevLett.92.018302.
[6] JING D，BHUSHAN B. Boundary slip of superoleophilic，oleophobic，and superoleophobic surfaces immersed in deionized water，hexadecane，and ethylene glycol [J]. Langmuir，2013，29 (47): 14691-14700. DOI: 10.1021/la4030876.
[7] ESPINOSA-MARZAL R M，ARCIFA A，ROSSI A，et al. Microslips to “avalanches” in confined，molecular layers of ionic liquids [J]. The Journal of Physical Chemistry Letters，2013，5 (1): 179-184. DOI: 10.1021/jz402451v.
[8] 陳其樂，孔憲，盧滇楠，等. 外壁電荷性質(zhì)對雙壁碳納米管中水分子運(yùn)動行為的影響 [J]. 化工學(xué)報(bào)，2014，65 (1): 319-327. DOI: 10.3969/j.issn.0438-1157.2014.01.042. CHEN Q L，KONG X，LU D N，et al. Molecular simulation of outer surface charge on water transport through double-wall carbon nanotube [J]. CIESC Journal，2014，65 (1): 319-327. DOI: 10.3969/j.issn.0438-1157.2014.01.042.
[9] 張程賓，趙沐雯，陳永平，等. 流體密度對納通道內(nèi)流動滑移的影響 [J]. 化工學(xué)報(bào)，2012，63 (S1): 12-16. DOI: 10.3969/j.issn.0438-1157.2012.zl.003. ZHANG C B，ZHAO M W，CHEN Y P，et al. Effects of fluid density on velocity slip in nanochannels [J]. CIESC Journal，2012，63 (S1): 12-16. DOI: 10.3969/j.issn.0438-1157.2012.zl.003.
[10] 潘伶，高誠輝. 納米間隙潤滑劑季戊四醇四酯的壓縮性能分子動力學(xué)模擬 [J]. 機(jī)械工程學(xué)報(bào)，2015，51 (5): 76-82. DOI: 10.3901/JME.2015.05.076. PAN L，GAO C H. Molecular dynamics simulation on the compressibility of pentaerythritol tetra in nanogap [J]. Journal of Mechanical Engineering，2015，51 (5): 76-82. DOI: 10.3901/JME.2015.05.076.
[11] PAN L，GAO C H. Confined fluid density of a pentaerythritol tetraheptanoate lubricant investigated using molecular dynamics simulation [J]. Fluid Phase Equilibria，2015，385: 212-218. DOI: 10.1016/j.fluid.2014.11.014.
[12] TSIGE M，PATNAIK S S. An all-atom simulation study of the ordering of liquid squalane near a solid surface [J]. Chemical Physics Letters，2008，457 (4/5/6): 357-361. DOI: 10.1016/j.cplett.2008.04.026.
[13] KARNIADAKIS G，BESKOK A，ALURU N. Microflows and Nanoflows: Fundamentals and Simulation [M]. Berlin Heidelberg，New York: Springer，2004: 1-23.
[14] PERTSIN A J，GRUNZE M. Long-ranged solvation forces in a fluid with short-ranged interactions [J]. Journal of Chemical Physics，2003，118 (17): 8004-8009. DOI: 10.1063/1.1564051.
[15] LOI S，SUN G，F(xiàn)RANZ V，et al. Rupture of molecular thin films observed in atomic force microscopy (II): Experiment [J]. Physical Review E，2002，66 (3): 031602. DOI: 10.1103/PhysRevE.66.031602.
[16] FRANZ V，BUTT H-J. Confined liquids:? solvation forces in liquid alcohols between solid surfaces [J]. The Journal of Physical Chemistry B，2002，106 (7): 1703-1708. DOI: 10.1021/jp012541w.
[17] 陳天星. 利用原子力顯微鏡對近壁面受限液體性質(zhì)的研究 [D].北京: 清華大學(xué)，2011. CHEN T X. Study on properties of the confined liquids at solid surface with AFM[D]. Beijing: Tsinghua University，2011.
[18] SIVEBAEK I M，SAMOILOV V N，PERSSON B N J. Squeezing molecularly thin alkane lubrication films: layering transitions and wear [J]. Tribology Letters，2004，16 (3): 195-200. DOI: 10.1023/B:TRIL.0000009730.31175.82.
[19] 胡元中，王慧，郭炎. 超薄油膜潤滑的分子動力學(xué)模擬 [J]. 摩擦學(xué)學(xué)報(bào)，1995，15 (2): 138-144. DOI: 10.16078/j.tribology.1995.02.007. HU Y Z，WANG H，GUO Y. Molecular dynamics simulation of ultra thin film lubrication (I): Rigid molecule model [J]. Tribology，1995，15 (2): 138-144. DOI: 10.16078/j.tribology.1995.02.007.
[20] 王慧，胡元中，郭炎. 超薄潤滑膜界面滑移現(xiàn)象的分子動力學(xué)研究 [J]. 清華大學(xué)學(xué)報(bào) (自然科學(xué)版)，2000，40 (4): 107-110. WANG H，HU Y Z，GUO Y. Molecular dynamics study of interfacial slip behavior of ultrathin lubricating films [J]. J. Tsinghua Univ. (Sci. & Tech.)，2000，40 (4): 107-110.
[21] NAGAYAMA G，CHENG P. Effects of interface wettability on microscale flow by molecular dynamic simulation [J]. International Journal of Heat and Mass Transfer，2004，47: 501-513. DOI: 10.1016/j.ijheatmasstransfer.2003.07.013.
[22] ASPROULIS N，DRIKAKIS D. Boundary slip dependency on surface stiffness [J]. Physical Review E，2010，81 (6): 061503. DOI: 10.1103/PhysRevE.81.061503.
[23] ZHANG H，ZHANG Z，YE H. Molecular dynamics-based prediction of boundary slip of fluids in nanochannels [J]. Microfluidics and Nanofluidics，2012，12 (1-4): 107-115. DOI: 10.1007/s10404-011-0853-y.
[24] NOORIAN H，TOGHRAIE D，AZIMIAN A R. Molecular dynamics simulation of poiseuille flow in a rough nano channel with checker surface roughnesses geometry [J]. Heat and Mass Transfer，2014，50 (1): 105-113. DOI: 10.1007/s00231-013-1232-x.
[25] PLIMPTON S. LAMMPS molecular dynamics simulator[EB/OL]. [2015-5-15]. http://lammps.sandia.gov.
[26] PLIMPTON S. Fast parallel algorithms for short-range molecular dynamics [J]. Journal of Computational Physics，1995，117 (1): 1-19. DOI: 10.1006/jcph.1995.1039.
[27] SUN H. Ab initio calculations and force field development for computer simulation of polysilanes [J]. Macromolecules，1995，28 (3): 701-712. DOI: 10.1021/ma00107a006.
[28] SUN H. COMPASS:? an ab initio force-field optimized for condensed-phase applications overview with details on alkane and benzene compounds [J]. Journal of Physical Chemistry B，1998，102 (38): 7338-7364. DOI: 10.1021/jp980939v.
[29] NAVIER C L M H. Memoire sur les du movement des fluids [J]. Mem l'Acad Roy Sci l’Inst France，1823，6: 389-440.
[30] VADAKKEPATT A，DONG Y，LICHTER S，et al. Effect of molecular structure on liquid slip [J]. Physical Review E，2011，84 (6): 066311. DOI: 10.1103/PhysRevE.84.066311.
[31] MENDONCA A C F，F(xiàn)OMIN Y D，MALFREYT P，et al. Novel ionic lubricants for amorphous carbon surfaces: molecular modeling of the structure and friction [J]. Soft Matter，2013，9 (44): 10606-10616. DOI: 10.1039/C3SM51689J.
[32] LIM R，O’SHEA S J. Solvation forces in branched molecular liquids [J]. Physical Review Letters，2002，88 (24): 246101. DOI: 10.1103/PhysRevLett.88.246101.
[33] LIM R，LI S F Y，O’SHEA S J. Solvation forces using sample-modulation atomic force microscopy [J]. Langmuir，2002，18 (16): 6116-6124. DOI: 10.1021/la011789+.
[34] ZHENG X，ZHU H，KOSASIH B，et al. A molecular dynamics simulation of boundary lubrication: the effect of n-alkanes chain length and normal load [J]. Wear，2013，301 (1/2): 62-69. DOI: 10.1016/j.wear.2013.01.052.
Molecular dynamics simulation of boundary slip in nanogap: effect of shear velocity
PAN Ling，GAO Chenghui
(School of Mechanical Engineering and Automation，F(xiàn)uzhou University，F(xiàn)uzhou 350108，F(xiàn)ujian，China)
Abstract:Molecular dynamics (MD) simulations using the polymer consistent force field (PCFF) were adopted to investigate the density，molecular and velocity distributions of lubricant squalane in nanogap at 293 K，three different film thicknesses and a wide range of shear velocities. The lubrication mechanism and boundary slip were analyzed. The results showed that the lubricant atoms tended to form layers parallel to the confining wall. The distances between the layers of lubricant atoms were irregular rather than broadening far away from the walls and were about three to four times the length of C C bond in the squalane. The tendency of lubricant atoms to form layers and the density of solid-like layer increased with decreasing film thickness. It was clearly to find the boundary slip at the solid-liquid interface from the velocity profile. The slip lengths increased with increasing velocity of substrates at the beginning，and then decreased. When the film thickness was 3.44 nm，the maximum slip length was 4.35 nm at the substrate velocity of 22.8 m·s-1. According to the simulations，the relationship between the slip length and the shear velocity was given.
Key words:molecular simulation; nanostructure; polymers; boundary slip; thin film lubrication; squalane
DOI：10.11949/j.issn.0438-1157.20150755
中圖分類號：TH 117.2
文獻(xiàn)標(biāo)志碼：A
文章編號：0438—1157（2016）04—1440—08
基金項(xiàng)目：國家自然科學(xué)基金項(xiàng)目（51175085）；清華大學(xué)摩擦學(xué)國家重點(diǎn)實(shí)驗(yàn)室開放基金項(xiàng)目（SKLTKF13A09）；福建省自然科學(xué)基金項(xiàng)目（2016J01226）。
Corresponding author:Prof. GAO Chenghui，gch@fzu.edu.cn
猜你喜歡
中國科技博覽(2017年5期)2017-04-25 23:36:47
價(jià)值工程(2017年11期)2017-04-18 17:58:11
山東工業(yè)技術(shù)(2016年19期)2016-10-29 16:12:54
中國科技博覽(2016年18期)2016-10-19 07:34:50
科技傳播(2016年10期)2016-07-15 23:11:56
中國科技博覽(2016年13期)2016-07-13 01:53:14
中國科技博覽(2016年10期)2016-04-29 04:18:00
中國科技博覽(2016年1期)2016-04-25 12:04:06
科技與創(chuàng)新(2016年5期)2016-03-17 02:33:51
湖南大學(xué)學(xué)報(bào)·自然科學(xué)版(2014年7期)2014-11-28 17:33:37
