毋飛鵬，張 波，劉子龍，唐 禹

(武漢理工大學 理學院 物理系，武漢 430070)

T矩陣方法計算雙層球形粒子的受力

毋飛鵬，張 波，劉子龍*，唐 禹

(武漢理工大學 理學院 物理系，武漢 430070)

為了分析聚焦光束對多層粒子的捕獲效率，結合矢量衍射積分、T矩陣方法以及Maxwell應力張量積分，通過理論推導給出了雙層球形粒子的T矩陣的詳細表達式，并對雙層球形粒子在聚焦光場中的受力進行了數(shù)值計算，詳細分析了內層折射率和內層尺寸對光場捕獲效率的影響。結果表明，只有內層折射率在一定范圍內，聚焦光束對雙層球形粒子才具有捕獲作用，隨著內層折射率增加，最大后向捕獲效率先增加后減小至零，對于空心粒子，內層尺寸越大，聚焦光束對粒子的捕獲作用越弱，且平面波的捕獲作用比高斯光束更強。此雙層球形粒子的受力計算可以拓展到多層的復雜粒子的情形。

激光技術；光學捕獲；T矩陣；雙層球形粒子

引 言

自從ASHKIN在實驗中首次觀察到激光對微粒的輻射力、并利用一束強聚焦激光對微粒實現(xiàn)了穩(wěn)定的3維捕獲以來[1-2]，光鑷技術得到了巨大的發(fā)展，并廣泛運用到原子物理學、生物醫(yī)藥、生命科學等領域[3-5]。光鑷技術已被用于細胞分選[6]、測量細胞彈性[7]、研究染色體轉錄[8]、研究驅動蛋白[9]、測量液體粘滯系數(shù)[10]等。在光鑷技術的相關研究領域中，理論上計算各種微粒的受力具有重要的指導意義。射線模型和瑞利模型可以計算尺寸遠大于和遠小于入射波長的粒子的受力[11-12]，對于尺寸近似波長的粒子，可以利用多種數(shù)值方法來計算散射場，并利用麥克斯韋應力張量積分求得粒子的受力[13-16]。

T矩陣方法于1965年被首次提出[17]，并在1971年得到進一步完善[18]，在1974年被推廣到求多層粒子的散射場[19]，之后T矩陣得到了巨大的發(fā)展，被廣泛運用到計算光鑷對多種微粒的捕獲力上。T矩陣只與粒子的形狀和物理性質有關，適用于任意形狀的粒子，且能夠適用的粒子的尺寸范圍較廣。T矩陣方法可以用來計算介質均勻的橢球粒子在光場中的受力[20]，也可以用于分析桿狀粒子在光場中所受的力與力矩[21]，如今已可以利用工具箱計算常見粒子的受力[22]。光鑷捕獲和操縱生物粒子或鍍有金屬膜層的納米粒子時，這類粒子通常具有多層結構。在計算光鑷的捕獲力時，如果將粒子處理為單層結構，得到的結果是不夠準確的。目前，利用T矩陣方法計算多層粒子在光場中的受力的相關研究比較少，粒子各部分的參量對捕獲力的影響也沒有得到詳細的分析。雙層球形粒子的受力計算可以進一步推廣到任意光源對多層球形粒子的受力和力矩，因此，多層粒子受力的計算是以雙層結構的計算為基礎的。

本文中給出了雙層球形粒子的T矩陣的表達式，利用T矩陣方法結合矢量衍射積分和麥克斯韋應力張量積分計算了雙層球形粒子在聚焦光場中的受力，分析了內層折射率和內層尺寸對粒子受力的影響，討論了空心粒子在光場中的捕獲與運動軌跡。首先介紹了T矩陣方法計算粒子受力的相關理論；然后計算了雙層球形粒子在光場中的受力，并對數(shù)值結果進行了詳細的討論；最后結合理論和數(shù)值結果總結了本文中的工作。

1 理論分析

利用T矩陣方法計算粒子受力的基本原理是將粒子周圍的入射場Einc(r)和散射場Esca(r)都展開如下式所示的矢量球形波函數(shù)的疊加，其中入射場的展開系數(shù)amn,bmn和散射場的展開系數(shù)pmn,qmn可以通過T矩陣聯(lián)系起來，通過amn,bmn和T矩陣可以得到散射場，進而得到粒子周圍的場分布，然后通過麥克斯韋應力張量積分得到粒子的受力。

T矩陣方法中入射場Einc(r)一般是透鏡系統(tǒng)形成的聚焦場。如圖1所示，入射光經(jīng)過一個高數(shù)值孔徑的透鏡聚焦之后入射到一個雙層球形粒子之上，利用Richard-Wolf矢量衍射積分方法以及角譜展開方法將聚焦之后的場展開為一系列平面波的疊加，利用平面波的展開系數(shù)[23]，可得入射場的展開系數(shù)為:

Fig.1 The particles trapped by the focused beam

(2)式給出了粒子表面入射場的展開系數(shù)，為了得到散射場，還需要求得雙層球形粒子的T矩陣。雙層粒子的T矩陣可以表示為:

T12=(-K1,Ⅰ,ⅠK1,Ⅰ,Ⅲ-1+K1,Ⅲ,ⅠK2,Ⅰ,ⅠK2,Ⅰ,Ⅲ-1K1,Ⅰ,Ⅲ-1)×

(1-K1,Ⅲ,ⅢK2,Ⅰ,ⅠK2,Ⅰ,Ⅲ-1K1,Ⅰ,Ⅲ-1)-1=

式中，T1和T2為外層界面S1和內層界面S2構成的單一均勻球形粒子的T矩陣，T12是雙層粒子的T矩陣，K是一個與T矩陣等大小的矩陣[23]，其下標中的第1個阿拉伯數(shù)字表示積分界面的標號，后面兩個羅馬數(shù)字表示求解過程中矢量球形波函數(shù)的形式，矩陣K由4個子塊J1，J2，J3，J4組成:

在雙層球形粒子的情況下，J2≡0，J3≡0，

式中,δmm′nn′為Kronecker函數(shù)(m=m′，且n=n′時取1；否則取0)，相對折射率ny,y-1=ny/ny-1(ny和ny-1分別表示第y層和第y-1層介質的絕對折射率)，定義參量xy=k0ryny-1/n0(ry表示第y層球形界面的半徑)，函數(shù)fu和fv在各種情況下的形式為:

通過雙層粒子的T矩陣與入射場的展開系數(shù)可求得散射場的展開系數(shù)，進一步得到粒子外的場分布，利用麥克斯韋應力張量積分，可以求得粒子的受力F。

式中,〈〉表示求時域平均，D,E,B,H分別為電位移矢量、電場強度、磁感應強度、磁場強度，I為單位張量，中括號的第2項為并矢，dS為積分的面積元矢量，上標*表示取相應矢量的共軛。

2 數(shù)值結果分析

2.1 內層折射率的影響

根據(jù)以上的理論分析，假設周圍環(huán)境的折射率n0=1.33(在水中捕獲微粒)，透鏡的焦距f=1mm，透鏡的數(shù)值孔徑為1.25。通過數(shù)值模擬，分別計算了雙層球形粒子在聚焦平面波和聚焦高斯光束中的受力情況，并繪制了捕獲效率Q與位置的關系曲線，Q=Fc/(n0P)(c為真空中的光速，P為光功率)無量綱。若平面波作為入射光源，真空波長λv=1064nm(對應生物粒子的透明波長)，振幅l0=1(歸一化)，沿x方向偏振。設粒子的外層半徑r1=1μm(第3類粒子)，內層半徑r2=0.5μm，外層折射率n1=1.43，內層折射率分別取空氣、水和聚苯乙烯的折射率，即n2為1,1.33和1.59，分別對應空心粒子、內層折射率小于和大于外層的3種情形。這3種情況下，微粒在光軸上的軸向捕獲效率Qz以及在焦平面x軸上橫向捕獲效率Qx如圖2所示，其中Qb,m為最大后向捕獲效率。

從圖2可以看到，當n2=1(空心粒子，實線)，粒子在軸上有兩個平衡點，分別位于z=-1.02μm和z=1.43μm(后面的討論中將-z軸和+z軸的平衡點分別稱為A和B)，粒子的橫向平衡點處于離軸的位置x=0.78μm。此時，粒子的內層折射率小于外層，甚至小于環(huán)境折射率，由于內層結構產(chǎn)生的梯度力都會將粒子推離焦點[24]，而外層折射率大于環(huán)境，外層產(chǎn)生的梯度力將粒子拉回焦點。由于兩者的共同作用，在光軸上產(chǎn)生了兩個平衡點，同理在橫向上也有兩個橫向平衡點(對稱性導致-x軸上還有一個平衡點，預測y軸上情況類似，而且這里的平衡點只是橫向上受力平衡)。在此處因外層結構產(chǎn)生的梯度力會將粒子拉向光軸，而內層結構產(chǎn)生的力則會將粒子推離光軸，在粒子偏離光軸時，推離的力會隨著光強的減弱而變小，最后與拉向光軸力達到平衡，形成一個偏離光軸的橫向平衡點。當n2=1.33時(點線)，粒子的受力曲線及平衡點產(chǎn)生的原因與n2=1時是相似的。當n2=1.59時(點劃線)，在軸向上只有一個平衡點位于z=0.07μm，而橫向上平衡點位于光軸上，因為此時內層的折射率大于外層(也大于環(huán)境)，光場對內層與外層的梯度力都指向焦點處，所以只有一個軸向平衡點，橫向上由于對稱性導致平衡點位于光軸上。

Fig.2 Trapping efficiency of plane wave with different n2

前面的分析發(fā)現(xiàn)n2=1時在軸上有兩個平衡點，為了判斷兩個平衡點是否都能穩(wěn)定的捕獲粒子，圖3給出了平面波入射時，粒子在y-O-z平面上受力的矢量圖。圖中箭頭的長度代表受力的大小，箭頭的方向代表受力的方向，插圖的箭頭只表示受力的方向而不表示大小。

由圖3可以看出，粒子在焦平面的y軸上的橫向受力和x軸類似，在離軸位置有一個橫向平衡點(只是橫向受力為0，軸向受力不為0)，和前面預期的結果是一致的。而左下角的插圖表明，在平衡點A處，粒子的橫向受力平衡但是不穩(wěn)定，粒子在受到環(huán)境因素擾動的情況下會逐漸遠離A點。右下角的插圖表示粒子在平衡點B處的橫向和軸向都處于穩(wěn)定狀態(tài)。所以，位于B點附近的微粒會直接移動到B點，并被穩(wěn)定捕獲在該位置。而位于A點附近的粒子會先向A點移動，然后由于環(huán)境干擾和布朗運動，沿著圖3實線所示的路徑緩慢移動到B點并被穩(wěn)定捕獲。綜上所述，只有+z軸的B點是穩(wěn)定的平衡點，內層折射率n2為1,1.33,1.59時，聚焦的平面波可以將粒子分別穩(wěn)定地捕獲在z為1.43μm,0.86μm,0.07μm處。

Fig.3 Vector graph of trapping force of particles iny-O-zplane(enlarged: directions of trapping force nearAandB)

如果入射光是高斯光束并且束腰位于出瞳位置，束腰半徑w0=1mm，只將振幅分布改為高斯分布，其它條件不變。如圖4所示，在3種情況下粒子的受力曲線以及各個平衡點的形成原因與平面波相似。需要注意的是：n2=1(實線)時不再形成B平衡點。射線模型可以用來解釋這個現(xiàn)象，光束中的所有光線經(jīng)過透鏡聚焦后，靠外的光線提供了更大的梯度力[11]。高斯光束相比于平面波而言，靠外的光線的功率占比很小，靠內的光線的功率占比較大，導致高斯光束提供的梯度力較小，不能形成平衡點B。當粒子位于-z軸時，內層的空心結構使得粒子被沿-z方向推離焦點，而外層提供的梯度力和散射力都沿+z軸指向焦點，兩者達到平衡，形成了平衡點A。

Fig.4 Trapping efficiency of Gaussian beam with different n2

2.2 內層尺寸的影響

除了內層折射率的影響，內層尺寸也是一個影響捕獲效率的因素。利用平面波和高斯光作為捕獲光源，設外層折射率n1=1.43，內層折射率n2=1，外層半徑r1=1μm，內層半徑r2分別為0.2μm,0.4μm,0.6μm，0.8μm，其它參量不變。微粒軸向和橫向上的捕獲效率如圖5所示。

Fig.5 Trapping efficiency with different r2

a—plane wave, axial b—plane wave, transverse c—Gaussian beam, axial d—Gaussian beam, transverse

如圖5a所示，當r2=0.2μm時(實線)，粒子在+z軸有一個平衡點位于z=0.78μm。參照第2.1節(jié)中對平面波的分析，由于內層結構尺寸更小，粒子更接近于一個均勻的介質球，所以該平衡點是穩(wěn)定的平衡點。而粒子內層結構提供的推離焦點的力過小使得此時沒有-z軸的平衡點A。當r2=0.4μm時(虛線)，粒子的受力與第2.1節(jié)中r2=0.5μm時相似。當r2=0.6μm(點線)和r2=0.8μm(點劃線)時，在光軸上都只有不穩(wěn)定的平衡點A，分別在z=-1.15μm和z=-1.58μm，由于內層結構尺寸過大，外層提供的拉向焦點的力不能平衡內層提供的推離焦點的力，所以不能形成穩(wěn)定的平衡點B。從圖5a中看出,r2=0.6μm時，平衡點B剛好消失，說明當內層尺寸在0μm～0.6μm之間時，聚焦的平面波可以對雙層球形粒子進行穩(wěn)定的3維捕獲。橫向上，內層尺寸為0μm時，粒子相當于一個單層粒子，平衡位置位于光軸上，隨著內層尺寸的增加，焦平面上橫向的平衡位置越來越遠離光軸。

若高斯光束作為捕獲光源，內層尺寸變化對捕獲效率的影響如圖5c和圖5d所示，粒子的受力曲線和平面波的結果非常相似。其中r2=0.4μm時，高斯光束提供的后向捕獲效率已經(jīng)非常小，捕獲作用十分微弱，所以當內層尺寸在0μm～0.4μm之間時，聚焦高斯光束可以對雙層球形粒子進行穩(wěn)定的3維捕獲，內層尺寸繼續(xù)增加時不再穩(wěn)定捕獲。

2.3 對最大后向捕獲效率的影響

前面分析了幾種情況下粒子捕獲效率與位置的關系曲線，現(xiàn)在將進一步分析內層折射率和內層尺寸與最大后向捕獲效率Qb,m的關系。如圖6a所示，只有內層折射率處于一定的范圍內，粒子才能被穩(wěn)定捕獲，參考文獻[25]中在單層粒子的情況下也得到了相似的結論，隨著內層折射率增加，內層界面上功率反射率先減少后增加，相應的散射力也按同樣的規(guī)律變化，導致捕獲效率先增加再減小，如圖6a所示。平面波入射時(實線)，只有內層折射率處于1≤n2≤2.1，粒子(r1=1μm，r2=0.5μm，n1=1.43)才能被聚焦的平面波捕獲，最大后向捕獲效率先增加然后減小至0。而當入射光束為高斯光束時(虛線)，只有內層折射率在1.05～1.9之間時，粒子才能被聚焦的高斯光束穩(wěn)定捕獲，且最大后向捕獲效率變化趨勢與平面波的情況類似，可見，平面波捕獲的粒子的內層折射率變化范圍比高斯光束要大。而對于同一個內層折射率，聚焦平面波的最大后向捕獲效率要比聚焦高斯光束更大，其原因如前面的分析。

Fig.6 The maximum backward trapping efficiency of plane wave and Gaussian wave

a—refractive indexn2of the inner layer b—sizer2of the inner layer

如圖6b所示，若內層折射率為1(即為空心粒子)，Qb,m會隨著內層尺寸增加而減小直至為0。平面波入射時(實線)，內層尺寸在0μm≤r2≤0.6μm范圍內的粒子可以被捕獲；高斯光束入射時，內層尺寸變化范圍較小(0μm～0.45μm)，這和前面的結論是相符的?？梢钥闯?，相同內層尺寸的情況下，平面波的捕獲效果要比高斯光束更好。

3 結 論

通過T矩陣方法給出了雙層球形粒子的T矩陣的表達式，其結果可以推廣至任意多層球形粒子，通過數(shù)值計算給出了聚焦光場中不同位置的粒子所受的捕獲效率，并分析了粒子內層折射率以及內層尺寸對捕獲效率的影響。計算結果表明：粒子內層的折射率和尺寸都對捕獲效率有明顯的影響，平面波和高斯光束都可以捕獲雙層球形粒子，但是粒子的內層折射率需要處在一定范圍內，平面波作為捕獲光源時的范圍更大；隨著內層折射率的增加，最大后向捕獲效率會先增加然后減小到零；粒子的內層折射率一定的情況下，內層尺寸也會明顯影響捕獲效率，對于空心粒子，內層的空心越大，捕獲效果越弱；對同一個粒子，平面波的捕獲效果要比高斯光束更好。

[1] ASHKIN A. Acceleration and trapping of particles by radiation pressure[J]. Physical Review Letters, 1970, 24(4): 156-159.

[2] ASHKIN A, DZIEDZIC J M, BJORKHOLM J E,etal. Observation of a single-beam gradient force optical trap for dielectric particles[J]. Optics Letters, 1986, 11(5): 288-290.

[3] CHU S. Laser manipulation of atoms and particles[J]. Science, 1991, 253(5022): 861-866.

[4] BLOCK S M. Making light work with optical tweezers[J]. Nature, 1992, 360(6403): 493-495.

[5] YIN H, WANG M D, SVOBODA K,etal. Transcription against an applied force[J]. Science, 1995, 270(5242): 1653-1657.

[6] GRIMBERGEN J A, VISSCHER K, GOMES de MESQUITA D S,etal. Isolation of single yeast cells by optical trapping[J]. Yeast, 1993, 9(7): 723-732.

[7] WU J G, LI Y M, LU D,etal. Measurement of the membrane elasticity of red blood cell with osmotic pressure by optical tweezers[J]. Cryoletters, 2009, 30(2): 89-95.

[8] WEN J D, LANCASTER L, HODGES C,etal. Following translation by single ribosomes one codon at a time[J]. Nature, 2008, 452(7187): 598-603.

[9] BLOCK S M, GOLDSTEIN L S,SCHNAPP B J. Bead movement by single kinesin molecules studied with optical tweezers[J]. Nature, 1990, 348(6299): 348-352.

[10] WEI Y. Measurement of viscosity coefficient of liquid micro-area with laser optical tweezers[J]. Laser Technology, 2016, 40(5): 738-741(in Chinese).

[11] ASHKIN A. Forces of a single-beam gradient laser trap on a dielectric sphere in the ray optics regime[J]. Biophysical Journal, 1992, 61(2): 569-582.

[12] HARADA Y, ASAKURA T. Radiation forces on a dielectric sphere in the rayleigh scattering regime[J]. Optics Communications, 1996, 124(5/6): 529-541.

[13] CHEN J, NG J, WANG P,etal. Analytical partial wave expansion of vector bessel beam and its application to optical binding[J]. Optics Letters, 2010, 35(10): 1674-1676.

[14] CHEN J, NG J, WANG P,etal. Analytical partial wave expansion of vector bessel beam and its application to optical binding: erratum [J]. Optics Letters, 2011, 36(7): 1243-1243.

[15] YANG M, WU Y, REN K F,etal. Computation of radiation pressure force exerted on arbitrary shaped homogeneous particles by high-order bessel vortex beams using mlfma[J]. Optics Express, 2016, 24(24): 27979-27992.

[16] MILNE G, DHOLAKIA K, MCGLOIN D,etal. Transverse particle dynamics in a bessel beam[J]. Optics Express, 2007, 15(21): 13972-13987.

[17] WATERMAN P C. Matrix formulation of electromagnetic scattering[J]. Proceedings of the IEEE, 1965, 53(8): 805-812.

[18] WATERMAN P C. Symmetry, unitarity, and geometry in electromagnetic scattering[J]. Physical Review, 1971, D3(4): 825-839.

[19] PETERSON B,STROM S.T-matrix formulation of electromagnetic scattering from multilayered scatterers[J]. Physical Review, 1974, D10(8): 2670-2684.

[20] YAN Sh H,YAO B L. Transverse trapping forces of focused gaussian beam on ellipsoidal particles[J]. Journal of the Optical Society of America, 2007, B24(7): 1596-1602.

[21] BAREIL P B,SHENG Y. Angular and position stability of a nanorod trapped in an optical tweezers[J]. Optics Express, 2010, 18(25): 26388-26398.

[22] NIEMINEN T A, LOKE V L Y, STILGOE A B,etal. Optical tweezers computational toolbox[J]. Journal of Optics, 2007, A9(8): S196-S203.

[23] MISHCHENKO M I, TRAVIS L D,LACIS A A, Scattering, absorption, and emission of light by small particles[M].New York,USA: Cambridge University Press, 2002: 115-190.

[24] ZHAN Q. Radiation forces on a dielectric sphere produced by highly focused cylindrical vector beams[J]. Journal of Optics, 2003, A5(3): 229-232.

[25] KAWAUCHI H, YONEZAWA K, KOZAWA Y,etal. Calculation of optical trapping forces on a dielectric sphere in the ray optics regime produced by a radially polarized laser beam[J]. Optics Letters, 2007, 32(13): 1839-1841.

Calculationoftrappingforceondouble-layersphericalparticlesusingTmatrixmethod

WUFeipeng,ZHANGBo,LIUZilong,TANGYu

(Department of Physics, School of Science, Wuhan University of Technology, Wuhan 430070, China)

In order to analyze the trapping efficiency of multi-layer particles in a focused beam, theTmatrix of the double-layer spherical particle was deduced. The trapping force on double-layer spherical particles was numerically calculated by means of the vector diffraction integral combining theTmatrix method and the Maxwell stress tensor integral. Effect of the refractive index and size of the inner layer on trapping efficiency was discussed in detail. The focused beam can trap the double-layer spherical particles only when the refractive index of the inner layer is appropriate. The maximum backward trapping efficiency will increase at first and decrease to zero finally when the refractive index increases. For hollow particles, the bigger the hollow is, the weaker the trap is. Besides, the trap formed by plane wave is stronger than Gaussian beam. The calculation method of the force on double-layer spherical particles can be expanded to the case of multilayer complex particles.

laser technique; optical trapping;Tmatrix; double-layer spherical particle

1001-3806(2018)01-0072-06

中央高?；究蒲袠I(yè)務費資助項目(WUT:2016-IA-008)

毋飛鵬(1991-)，男，碩士研究生，主要從事光鑷技術的研究。

*通訊聯(lián)系人。E-mail:zlliu_72@whut.edu.cn

2017-03-31；

2017-05-09

O436

A

10.7510/jgjs.issn.1001-3806.2018.01.014