不均勻環(huán)境中二聚體微納米馬達的動力學(xué)性質(zhì)

金英俊，董康杰，溫正城

(杭州電子科技大學(xué)理學(xué)院，浙江 杭州 310018)

0 引 言

自驅(qū)動微納米馬達是一種尺寸在微納米量級的微小器件，可以將其他形式的能量轉(zhuǎn)化為動能，實現(xiàn)自主運動[1]。目前，自驅(qū)動微納米馬達的探索尚處于早期階段，譬如Paxton等[2]研究以H2O2為燃料的Au/Pt雙元金屬棒馬達在H2O2中的運動；Gao等[3]先后制備了以H2O2為驅(qū)動力的Al/Ga和Mg微米馬達。馬達動力學(xué)性質(zhì)的研究往往在較簡單的均勻環(huán)境中進行，主要關(guān)注的是馬達的運動速度、驅(qū)動機理、功能應(yīng)用等方面，燃料時空動力學(xué)特性引起馬達運動性質(zhì)變化的相關(guān)研究甚少。

眾所周知，生物系統(tǒng)中的環(huán)境絕大部分是不均勻耗散體系，展示復(fù)雜的時空動力學(xué)現(xiàn)象。比如，在不均勻的環(huán)境中，細菌往往趨向有營養(yǎng)物質(zhì)區(qū)域，遠離有毒性物質(zhì)區(qū)域，展示趨利避害的重要性質(zhì)。隨著相關(guān)研究的深入，不均勻環(huán)境對馬達動力學(xué)的影響引起研究者關(guān)注[4-5]。本文通過構(gòu)造燃料分布不均勻環(huán)境，運用混合的分子動力學(xué)和多粒子碰撞動力學(xué)(Molecular Dynamics-Multi-Particle Collision Dynamics，MD-MPC)方法，探究自驅(qū)動二聚體微納米馬達的動力學(xué)性質(zhì)，以及馬達間相互作用對其形成團簇的影響。

1 系統(tǒng)模型與仿真

本文主要針對自驅(qū)動二聚體微納米馬達在不均勻環(huán)境中的動力學(xué)性質(zhì)和群集效應(yīng)展開研究。建立尺寸為V=Lx×Ly×Lz=100×100×20的三維系統(tǒng)進行仿真模擬。系統(tǒng)中，隨機放置大量的燃料粒子F以及不具有化學(xué)活性的惰性溶液粒子S，粒子的總數(shù)目NT=NF+NS，其中NF為燃料粒子F的數(shù)目，NS為惰性粒子S的數(shù)目。同時，系統(tǒng)中放入ND個二聚體微納米馬達，每個二聚體微納米馬達由催化球C和非催化球N組成。當(dāng)隨機運動的溶劑粒子碰到馬達的催化球C時，C表面發(fā)生不可逆的催化反應(yīng)，F(xiàn)+C→P+C，將燃料粒子F轉(zhuǎn)化為產(chǎn)物粒子P；而在非催化球N處，沒有催化反應(yīng)發(fā)生。二聚體微納米馬達的構(gòu)型如圖1所示，主要包括催化球C和非催化球N。

圖1 二聚體馬達構(gòu)型

將溶液中的F粒子、S粒子以及反應(yīng)生成的P粒子作為質(zhì)量均為m的質(zhì)點。為了限制馬達在Z方向上的運動，采用9-3蘭納-瓊斯(Lennard-Jones，LJ)排斥勢進行描述，具體表達式如下：

(1)

式中，r為2個粒子之間的距離，εW表示馬達與上下墻壁之間相互作用的能量參數(shù)，σW表示馬達與上下墻壁之間相互作用的距離參數(shù)。令σW=LZ/2，其中LZ表示系統(tǒng)在Z方向上的尺寸，由于這個勢的存在，馬達只能在XY平面上做準二維運動[6]。因為壁面的存在，系統(tǒng)中的溶液點粒子不斷通過回彈碰撞與坐標Z=0和Z=LZ處的壁相互作用，與壁面碰撞后，小球的速度不變，并沿碰撞前的鏡面對稱方向彈回。對于含有多達1023量級的微觀粒子的宏觀系統(tǒng)，通過對整個微觀系統(tǒng)的演化來獲得宏觀物理量并不容易。為了獲得較為準確的物理量，在系統(tǒng)的X方向和Y方向上，采用周期性邊界條件來正確描述二聚體馬達的動力學(xué)性質(zhì)。

對于溶液中的粒子與二聚體微馬達之間的相互作用，采用12-6形式的LJ排斥勢進行描述，具體表達式如下：

(2)

式中，λ為C或N，表示馬達的C球或N球，μ為F，S或P，表示溶液中的粒子，rc=21/6σ為馬達的截斷距離，θ(x)為Heaviside函數(shù)。當(dāng)r≤rc時，可以通過式(2)計算得出粒子與馬達之間的作用勢；當(dāng)r>rc時，該作用勢被截斷[7]。除了產(chǎn)物粒子P對N球的能量參數(shù)不同外，其他所有溶液粒子與C球和N球的相互作用能量參數(shù)相同，即εCF=εCS=εNF=εNS=εF≠εP，其中εP=εNP。

對于含有多個二聚體微馬達的系統(tǒng)，不同二聚體的單體之間也通過排斥LJ勢相互作用，其中能量參數(shù)為εD，距離參數(shù)為σD=σλJ+σλK，J=1,…,ND和K=1,…,ND表示系統(tǒng)中的二聚體微馬達。

采用MD-MPC方法對整個系統(tǒng)進行時間演化，MPC-MD方法由分子動力學(xué)(Molecular Dynamics，MD)和多粒子碰撞動力學(xué)(Multi-Particle Collision Dynamics，MPC)構(gòu)成[8-9]。采用MD模擬溶液粒子運動的過程，通過牛頓運動方程演化溶液粒子與馬達的相互作用、馬達與馬達之間的相互作用，其步長為τMD。采用MPC模擬溶液粒子之間的相互作用[10]。經(jīng)過n步MD模擬后，采用MPC進行1次溶液粒子隨機碰撞模擬，其步長τMPC=nτMD。用MPC模擬溶液粒子隨機碰撞過程中，本文把系統(tǒng)分為尺寸a=1的眾多立方體小格子，在每個單獨的方格內(nèi)，溶液粒子發(fā)生碰撞。標記為ξ號方格內(nèi)的i號粒子的速度為：

(3)

數(shù)值模擬時，系統(tǒng)中所有的物理量均為無量綱LJ單位，分別基于能量εα，質(zhì)量m0以及距離σd進行設(shè)置：r/σ→r，t(ε/(mσ2)1/2)→t和KBT/ε→T。具體參數(shù)如下：粒子F，S和P的質(zhì)量均為m0=1。馬達的C球、N球的直徑分別為dC=2σC=2×2=4，dN=2σN=2×4=8，為了保證馬達懸浮在溶液中，其質(zhì)量隨各自直徑變化而變化。本文模擬中，mC=335，mN=2 688，每個馬達的C球、N球之間的核間距固定為R=σC+σN+ξ=2+4+1=7。系統(tǒng)溫度設(shè)定為KBT=1/6。采用MD模擬溶液粒子運動過程中，時間步長ΔτMD=0.01。采用MPC模擬溶液粒子之間相互作用過程中，納米馬達的旋轉(zhuǎn)角為π/2，ΔτMPC=50×0.01=0.5，即n=50步進行1次碰撞，同時伴隨1次能量轉(zhuǎn)換。系統(tǒng)中的粒子總數(shù)NT=2.0×106，平均密度為10，即每個小方格子內(nèi)含有10個溶液粒子。能量參數(shù)為εS=εF=5，εP=0.1，2個馬達之間的εD=5。

初始化系統(tǒng)時，確保二聚體微馬達位于系統(tǒng)中的一側(cè)，放置規(guī)則如下：系統(tǒng)含有單個馬達情形下，N球球心在X和Y方向上需滿足10≤XN≤90,10≤YN≤90，C球相對于N球的角度范圍為[0,2π]，同時滿足5≤XC≤95，5≤YC≤95。當(dāng)系統(tǒng)含有多個馬達時，還需滿足如下規(guī)則：任意2個馬達α、β的N球之間的距離dNαNβ≥2×21/6σN，C球之間的距離dCαCβ≥2×21/6σC，C球和N球之間的距離dCαNβ≥2×21/6(σC+σN)。溶液粒子需滿足條件如下：與任意馬達C球的距離dSiCα≥21/6σC，與任意馬達N球之間的距離dSiNα≥21/6σN。

2 模擬結(jié)果與分析

2.1 單個馬達的自驅(qū)動行為

為了研究馬達在不均勻環(huán)境下的群集效應(yīng)與動力學(xué)性質(zhì)，分別針對環(huán)境中只有惰性溶液粒子S和燃料粒子F這2種情形進行模擬實驗，且只存在單個二聚體微馬達。

當(dāng)環(huán)境中只含惰性溶液粒子S時，溶液粒子S與馬達的C球、N球之間沒有化學(xué)反應(yīng)，馬達周圍沒有產(chǎn)生相應(yīng)的梯度差，又粒子S和C球、N球之間的能量參數(shù)設(shè)置為εSC=εSN=5，所以，此時沒有驅(qū)動力促使馬達做定向運動，馬達在系統(tǒng)中只做隨機的布朗運動。當(dāng)環(huán)境中只含燃料粒子F時，馬達的C球表面發(fā)生催化反應(yīng)生成產(chǎn)物粒子P的梯度場，雖然C球周圍的P粒子場是均勻的，但在N球周圍的P粒子梯度場是不均勻的，P粒子濃度沿著C-N的軸向降低，加上能量參數(shù)εF=5和εP=0.1(εF>εP)的不同，產(chǎn)物粒子P對二聚體N球產(chǎn)生的作用力小于其他溶液粒子對N球的推動，從而在N球上產(chǎn)生N-C指向的凈力，使馬達在系統(tǒng)中產(chǎn)自驅(qū)動。本文研究的二聚體微馬達基于自擴散泳的機理，文獻[12]通過模擬驗證了自擴散泳機理，并且制備出這類二聚體微馬達。

圖2 馬達的均方位移ΔL2隨時間t變化的曲線

2.2 馬達在單個圓形燃料區(qū)中的自驅(qū)動與群集

為了探究馬達在不均勻環(huán)境下的行為，對單個馬達(ND=1)在單個圓形燃料場中的行為進行模擬實驗。以坐標(25,25)為圓心，r=25為半徑建立1個圓形燃料場，區(qū)域內(nèi)的溶液粒子均為具有化學(xué)活性的粒子F。初始時刻，系統(tǒng)中均勻分布大量惰性溶液粒子S，由于粒子S做無規(guī)則布朗運動，為了確保圓形燃料場內(nèi)均為F粒子，將運動到區(qū)域邊界的S粒子受到外界催化轉(zhuǎn)變?yōu)榱Ｗ覨。單個馬達按照第1節(jié)所述規(guī)則進行放置。

圖3 二聚體微馬達的單個燃料場中運動及軸向速度

從圖3(b)可以看出，t在0～810和16 634～17 475這2個時間段中，馬達的軸向速度峰值遠大于其他時間段的軸向速度，馬達受到額外驅(qū)動力(FFN>FPN)的推動，在燃料區(qū)產(chǎn)生自驅(qū)動。

為了進一步探究馬達在受限環(huán)境中的行為，在上述單個燃料場系統(tǒng)中放入15個二聚體微馬達，即ND=15。除了按多個馬達放置規(guī)則來放置馬達外，還要確保惰性S溶液區(qū)和燃料區(qū)均有馬達。為了觀察更多二聚體微馬達的構(gòu)型，將馬達的運行步數(shù)增加到3.0×107。圖4為15個馬達在單個燃料場系統(tǒng)中的瞬態(tài)構(gòu)型。

圖4 單個燃料場內(nèi)，馬達在不同時刻的瞬態(tài)構(gòu)型

從圖4可以看出，初始時刻時，15個微納米馬達不均勻地分布在系統(tǒng)中，大部分馬達分布在惰性溶液粒子區(qū)，少數(shù)馬達位于燃料區(qū)內(nèi)。惰性溶液粒子區(qū)內(nèi)的馬達只做布朗運動，即緩慢進行自擴散，當(dāng)擴散到S溶液粒子區(qū)與燃料區(qū)的交界處時，由催化反應(yīng)產(chǎn)生的N球指向C球的力推動馬達進入燃料區(qū)。當(dāng)然，圓形燃料區(qū)內(nèi)的馬達也會跑出燃料區(qū)進入非燃料區(qū)。經(jīng)過一段時間演化后，多個二聚體微馬達在燃料區(qū)呈群聚狀態(tài)。當(dāng)燃料區(qū)中跑入多個馬達后，由于化學(xué)反應(yīng)，馬達附近產(chǎn)生P粒子的梯度場，考慮到εF>εP，其他馬達探測到該梯度場后，其N球向該馬達的C球靠近，從而形成團簇。C球和N球的聚集排列展現(xiàn)出穩(wěn)定的構(gòu)型，跑出燃料區(qū)的馬達以及原本在溶液區(qū)的馬達則不會產(chǎn)生群集，這是一個十分有趣的現(xiàn)象。

在探究二聚體微馬達形成團簇的過程中，對15個馬達進行標記，在隨機數(shù)為5時，對標號為4的二聚體微馬達進行追蹤，得到其軸向速度如圖5所示。

從圖5可以看出，馬達4的初始時刻軸向速度較大，說明位于系統(tǒng)的燃料區(qū)內(nèi)；隨著模擬實驗的推進，馬達4跑到惰性粒子溶液區(qū)，此時的軸向速度減?。辉?6 256～16 666時間段內(nèi)，其軸向速度處于高值，并于t=16 666時的軸向速度達到最大，說明此時馬達4已重新進入燃料區(qū)。以馬達4為中心，通過計算其他馬達與該馬達之間的距離發(fā)現(xiàn)，這些馬達相互靠近并形成團簇。最終，系統(tǒng)中多個二聚體微馬達以團簇運動的速度繼續(xù)運動。圖6為馬達4的軸向速度分布直方圖，通過計算得到的數(shù)學(xué)期望μ=0.002 5,σ=0.009 6。

圖5 隨機數(shù)為5時，4號二聚體微馬達的軸向速度

圖6 馬達4的軸向速度分布直方圖

2.3 馬達在多個圓形燃料區(qū)環(huán)境中的自驅(qū)動與群集

通過改變?nèi)剂蠄龅臄?shù)目，繼續(xù)研究馬達在燃料區(qū)中的動力學(xué)性質(zhì)與群集行為。在大小為Lx×Ly=100×100的準二維平面系統(tǒng)中進行模擬實驗，設(shè)置5個大小相同的圓形燃料區(qū)，圓心分別為(20，20)，(20，80)，(80，20)，(80，80)，(50，50)，半徑r=15，15個納米二聚體微馬達不均勻分布在整個系統(tǒng)中，15個馬達在系統(tǒng)中的構(gòu)型如圖7所示。

圖7 5個燃料區(qū)內(nèi)，馬達在不同時刻的瞬態(tài)構(gòu)型

從圖7可以看出，初始時刻時，燃料區(qū)和非燃料區(qū)均有馬達。隨著時間的推移，多個燃料區(qū)發(fā)生催化反應(yīng)，t=30 000時，系統(tǒng)呈現(xiàn)相對穩(wěn)定的馬達團簇構(gòu)型。由于單個燃料區(qū)域無法容納全部馬達，所以馬達在不同的燃料區(qū)內(nèi)形成團簇。團簇的形成區(qū)域是隨機的，與燃料區(qū)的排布和大小有關(guān)。

對于富含燃料區(qū)和非燃料區(qū)的不均勻環(huán)境，馬達在燃料豐富區(qū)運動時，呈現(xiàn)自驅(qū)動行為，有較高的活性。需要指出的是，由于本文研究的馬達是化學(xué)場自驅(qū)動，產(chǎn)生的非均勻梯度場驅(qū)動馬達，因此，馬達之間可通過各自的梯度場發(fā)生相互耦合，這種耦合效應(yīng)展示了與其他不同類型微納米馬達的不同性質(zhì)。由于耦合梯度場的存在，馬達之間相互吸引，在富集燃料區(qū)形成團簇。與其他類型馬達不同，馬達在快速運動區(qū)域反而出現(xiàn)群集現(xiàn)象。當(dāng)燃料區(qū)較大時，可以聚集所有的馬達；當(dāng)燃料區(qū)較小時，可以在不同燃料區(qū)分別聚集密集的團簇；在缺乏燃料區(qū)時，馬達四處散開隨機運動，形成有趣的相分離現(xiàn)象。

3 結(jié)束語

本文運用混合的分子動力學(xué)和多粒子碰撞動力學(xué)MPC-MD方法研究由催化球和非催化球構(gòu)成的二聚體微馬達在不均勻燃料環(huán)境中的單體和多體動力學(xué)行為。研究發(fā)現(xiàn)，在燃料富集區(qū)域，微納米馬達展現(xiàn)自驅(qū)動行為并形成團簇；在非燃料區(qū)，馬達呈布朗運動，具有分散-聚集的內(nèi)在物理機制。為面向生物醫(yī)療功能的納米機器的設(shè)計提供借鑒。但是，本文未對馬達形成聚集的條件進行分析，后續(xù)將對此展開進一步研究，通過改變?nèi)剂蠄鱿嚓P(guān)參數(shù)來調(diào)控二聚體微納米馬達的群集區(qū)域。