運用四元數分析橢球微粒所受的光阱力?

張書赫 梁振 周金華

(安徽醫科大學生物醫學工程系,合肥 230032)

運用四元數分析橢球微粒所受的光阱力?

張書赫 梁振?周金華?

(安徽醫科大學生物醫學工程系,合肥 230032)

(2016年10月3日收到;2016年10月17日收到修改稿)

在光鑷的射線模型中,追跡光線在界面的折射反射光是較基本的也是較復雜的問題.傳統的幾何光學法計算光線的方位,對于一些不規則的物體來說,存在一定的困難度.本文提出使用四元數簡化空間光線追跡,從而可計算非球形顆粒的光阱力.以入射光線和界面法線的外積確定入射面的法線為旋轉軸;根據折射定律確定由入射光線旋轉到反射光線和折射光線的角度.將入射光線以四元數表示,根據四元數旋轉,可獲得反射光線和折射光線的空間矢量.根據光子的動量變化,求得光阱對微粒的作用力.本文以球差影響下的不同形變的橢球為計算示例,模擬了橢球在光阱下的動力學行為,結果表明橢球的橫向和軸向俘獲效率受到各方向變形系數的影響;球差增大降低了軸向的最大俘獲效率,穩定俘獲位置也隨之朝負軸偏離中心更遠;在固定球差作用下,最大軸向俘獲效率與軸向變形系數相關,在特定的變形下軸向俘獲效率變得較大.由此驗證了四元數方法的正確性、實用性與普遍性.

四元數,光線追跡,橢球微粒,光阱力

1 引 言

目前光鑷廣泛應用于單分子檢測馬達蛋白[1,2],DNA和RNA力學特性分析[3]、測量細胞膜彈性[4]和探測溶液的黏彈性[5]等方面.在光鑷應用中,光阱對不同顆粒的作用,以及如何低功率下獲得高光阱力一直吸引許多研究者不斷探索[6].在前人的光阱模擬和實驗中,光束的類型[7]、束腰半徑[8]、顆粒大小[9]和聚焦光束的數值孔徑[10]均會影響光阱力的大小.光阱力的模擬方法一直備受關注,可幫助我們更深入理解光阱俘獲小球的動力學行為,還可指導我們設計光鑷提升俘獲效率.

光阱力的模擬方法大致有兩類,一類是幾何光學模型(或者射線模型,ray optics,RO);另一類方法為電磁模型(electromagnetic model,EM).RO模型的算法原理簡單,適于模擬顆粒尺寸遠大于波長的光阱力,即使在顆粒尺寸與波長相近時,也能取得與實驗相近的結果,在光鑷發展的30年內不斷得到應用和推廣.

在RO模型中,Ashkin[11]提出的射線追蹤的模型可以很好地分析球形顆粒在光阱中的受力,并可進一步模擬玻璃-水界面引起球差對光阱力的影響[12,13],以及雙光纖夾持顆粒的受力[14].然而采用歐氏幾何求解向量空間角的方法,在計算橢球形顆粒的受力存在一定困難.我們使用空間向量表示光線,通過空間解析幾何的方法進行光線追跡,用向量的夾角表示入射光和反射光平面的角度關系,根據光子在界面反射和折射的動量改變,進而在三維方向上計算光子對橢球的作用力[15,16].盡管如此,采用解析幾何方法,使用矩陣變換的方式計算折射光、反射光的空間方位角度涉及到多次坐標變換,在后續計算中依然較為復雜[16].

四元數可以很方便地表示空間旋轉[17],它在狹義相對論、電磁學、量子力學得到運用[18],目前還廣泛用于計算機圖形學[19]、飛行器姿態描述、慣性導航領域[20].將四元數旋轉引入折射光、反射光空間方位角的計算中,將會簡化空間光線追跡的復雜程度.本文介紹如何在RO模型中使用四元數法計算強聚焦光束對橢球微粒的作用力,并通過該模擬方法展示不同形變的橢球形顆粒在光阱中的動力學行為.

2 四元數與旋轉

愛爾蘭數學家哈密頓(Hamilton)于1843年提出四元數的數學概念,提供了三維空間中運用這種超復數實現坐標變換的方法[17].

2.1 四元數的表示方法

四元數q可用下列超復數形式表達為q=a+bi+cj+dk,其中i,j,k為虛數單位,且滿足i2=j2=k2=ijk=-1,ij=-ji=k,jk=-kj=i,ki=-ik=j.若將bi+cj+dk視為向量V=(b,c,d),純量S=a,則四元數q可表示為q=S+V,還可以寫成q=[S,V]這樣的形式.

2.2 四元數的運算規則

設有兩個四元數分別為q1=[S1,V1],q2=[S2,V2],則兩個四元數加法規則滿足q1+q2=[S1+S2,V1+V2].若這兩個四元數相乘,則

可見四元數相乘的結果還是四元數,它不滿足乘法交換律,且不同于四元數的點乘與叉乘的規則.對于四元數q=[S,V]的共軛,滿足

四元數q=[S,V]的模為

若一個四元素的模等于1,則該四元數為單位四元數.對于四元數q的逆記為q-1,滿足qq-1=1.四元數的逆可由(2)和(3)式求出：

對于單位四元數,它的逆等于其共軛四元數.

2.3 四元數與旋轉

三維空間平面Σ上已知點A為(xA,yA,zA),點O′為(xO′,yO′,zO′),圖1表示將向量O′A繞點O′旋轉角θ到向量O′B.過點O′的Σ面法線方向單位向量u=(xu,yu,zu).用四元數表示向量旋轉時,旋轉角度滿足右手螺旋為正,反之為負.向量O′A=(xA-xO′,yA-yO′,zA-zO′),被旋轉向量O′A擴展為只有向量部分的四元數：

以單位向量u和旋轉角度θ構建的旋轉軸的四元數為

它為單位四元數.向量O′B對應的四元數qO′B滿足

根據(1)—(7)式,計算所得的qO′B的向量部分即為向量O′B.

圖1 空間平面Σ上O′A以法線u為軸旋轉θ角到O′BFig.1.RotatingO′AtoO′Bwith angleθaround the normal axisuin planeΣ.

3 四元數在RO模型中的應用

3.1 運用四元數表示折射和反射光線

在光鑷的RO模型中,運用四元數法能追蹤入射光線在界面的各種偏折.以光線由光疏介質入射到光密介質為例,如圖2所示,界面的法線N取光密介質指向光疏介質為正方向,旋轉軸向量v滿足

在平面Ω0內以點O為極點,法線N為極軸,以繞向量v按右手螺旋為極角正方向建立平面極坐標系.則L1的極角記為φ1=π+α,L3的極角為φ3=π+β,L2的極角為φ2=2π-α.由L1旋轉到L2為φ2-φ1,L1旋轉到L3為φ3-φ1.

圖2 光線在界面的折射和反射示意圖Fig.2.Sketch of refractive and reflected lights in the interface.

根據(5)式將L1擴展成四元數qL1=[0,L1],根據(6)式用于旋轉到反射光線L2的旋轉軸的四元數為

用于旋轉到折射光線L3的旋轉軸的四元數為

當光線從光密介質進入光疏介質時,同樣取光密指向光疏介質的界面法線,旋轉軸用(8)式表示.根據右手螺旋確定角度正方向,則旋轉軸的四元數依然能用(9)和(10)式表示,則可確定反射光線的四元數qL2與折射光的四元數qL3,進而求得L2和L3的方向向量.根據四元數旋轉向量的方法,很方便地計算出光線通過任意界面產生的反射光與折射光向量.

3.2 單根光線對界面折射和反射產生的作用力

當光子經過折射率不同的界面時,光子由于反射和折射導致動量發生改變,從而對界面產生力的作用.光子從光疏介質進入光密介質的界面發生折射和反射,設入射角為α,折射角為β,根據動量守恒,界面受到光子的作用力為[15]

其中R和T為能流Fresnel的反射和透射系數,nm為光線從光源發出傳播所在介質的折射率,nr=ns/nm,ns為微粒的折射率,dP為光線的功率元.當光子從光密介質進入光疏介質時,光子對界面的作用力為

由(11)和(12)式可知微粒所受到的力始終沿著法線,由微粒內部指向外部.我們采用俘獲效率Q=Fc/(nmP)表示微粒在光阱中的受力情況.

4 運用四元數法模擬橢球微粒所受的光阱力

由于高斯光束被高數值孔徑物鏡強聚焦后,束腰在焦平面上且其半徑非常小,從物鏡后瞳入射面到焦點的光束可近似作為錐形光束傳播.物鏡后瞳處高斯光束強度分布為I(r)=I0exp(-2r2/ω20),ω0為高斯光束的束腰半徑,I0為束腰中心光強,r為光線從后瞳入射點到光軸的距離.在以下計算中,設定高斯光束恰好充滿物鏡后瞳.

4.1 高斯光束的單根光線追跡

以物鏡焦點為坐標原點O,物鏡光軸為z軸建立空間直角坐標,橢球微粒中心E位于(xE,yE,zE),如圖3所示.ng為玻片折射率,NA為物鏡數值孔徑.物鏡匹配油的折射率與玻璃折射率相近,則物鏡焦距fobj=ngRobj/NA.物鏡后瞳任意一點H到物鏡中心的距離為r,從H點入射的光線經過物鏡偏折為IL0指向焦點,且IL0與z軸的夾角為γ1.根據正弦定理sinγ1=r/fobj.當玻片位于z=Zcg處,在玻璃-水界面折射率失配引起的球差作用下,光線在玻璃-水界面上進一步偏折偏離焦點形成IL,IL與z軸交于點Zl,它們的夾角為γ2.σ為O′H與x’軸的夾角,IL=(cosσsinγ2,sinσsinγ2,cosγ2),Zl=Zcg+|Zcg|·tanγ1/tanγ2. 當光線通過玻片之后,光線可看作由Zl點沿IL出射的光線.

圖3 單根光線入射至橢球的光線追跡Fig.3.Tracing a single ray striking an ellipsoidal bead.

設橢球微粒E的方程為

光線IL在橢球微粒E表面反射光線為K0,折射光線為L0,與橢球交于點M0,如圖4.在微粒內部,L0在顆粒-溶液界面依次反射Ln,折射為Kn,并依次與顆粒表面相交于Mn,n為自然數.根據(1)—(10)式,由點Mn與Ln-1可計算得出Kn,Ln,其中Mn可由Ln-1直線方程與(13)式聯立求解確定[16].

圖4 光線在橢球微粒表面的折射和反射Fig.4.Refractive and reflected rays on the interface of an ellipsoidal bead.

根據(11)式和(12)式,光線Ln偏折對微粒的作用力為FIL,n.則該光線對微粒產生的總的作用力為最后再將所有光線對微粒產生的作用力求和,即可得到光束對小球的合力.

4.2 橢球微粒受到的光阱力

微球所在介質折射率nm=1.33,微球折射率ns=1.59,NA=1.25,Robj=3 mm,ng=1.51,玻片位于Zcg=-10μm處.當aE=bE=cE時,橢球方程退化為球方程,本文中取小球半徑rbead=3μm.引入變形系數來描述橢球相對于球體的變形,δx=aE/rbead,δy=bE/rbead,δz=cE/rbead,它們分別表示橢球沿著x,y和z方向的變形.當δx=δy=δz=1時,微粒為半徑3μm的球.在球差的影響下,我們采用四元數計算球形微粒的受力曲線Qz和Qx與我們之前空間坐標旋轉方法獲得的結果一致[16].

圖5(a)表示顆粒的橫向俘獲率.δx=1時微球作為參考曲線,當δx=0.33時小球在z=0處為非穩定狀態,一旦小球發生橫向偏離便會受到橫向的推力加速離開中心位置.對于δx=0.5的微粒,當它橫向偏離光阱中心時,微粒所受到的回復力快速增大,并且最大回復力|Qmax|比其他情況更大.從圖5(b)中可知,非球形顆粒的軸向穩定性與δx有關.在球差的影響下,δx=0.33的微粒的俘獲率曲線與Qz=0沒有交點,此時微粒被推離光阱.δx=0.5時微粒受到的最大軸向回復力較δx=1稍小,比δx=2和δx=3的情況稍大.

圖5(c)表示軸向變形顆粒的橫向俘獲率.當顆粒軸向變形拉長時,最大回復力增大.當變形達到一定程度后再拉長顆粒,最大回復力反而變小,如δz=3的微粒比δz=2的微粒受到的光阱力要小.微粒軸向偏離光阱中心時,它所受到的光阱力的曲線如圖5(d)所示.微粒軸向變形越大,它在平衡位置時受力曲線斜率為負,且其絕對值越大,表明微粒的軸向穩定性越好.

光阱中顆粒離玻片距離|Zcg|不同時,玻璃-水界面球差對橢球顆粒的軸向俘獲效率存在影響.我們計算了不同δz的橢球沿z軸負方向最大俘獲效率Qz0,如表1所列.當δz相同時,|Qz0|隨著|Zcg|增大而減小.當玻片位置相同時,橢球微粒|Qz0|隨δz變化并非單調增大或減小,例如當Zcg=-10,-15或-20μm時,δz=2的橢球微粒的|Qz0|最大,然而當Zcg=-30μm時δz=4的橢球微粒的|Qz0|最大.

圖5 (網刊彩色)不同形變橢球微粒的光阱俘獲效率 (a)小球沿x方向變形后沿x方向偏離光阱中心時的光阱俘獲效率;(b)小球沿x方向變形后沿z方向(光軸方向)偏離光阱中心時的光阱俘獲效率;(c)小球沿z方向變形并橫向偏離光阱中心時的光阱俘獲效率;(d)小球沿z方向變形并軸向偏離平衡位置時的光阱俘獲效率Fig.5.(color online)Trapping efficiencyQof an ellipsoidal bead with different deformations：(a)Transverse trapping efficiencyQxfor the bead with deformation ofδx;(b)axial trapping efficiencyQzfor the bead with deformation ofδx;(c)Qxfor the bead with deformation ofδz;(d)Qzfor the bead with deformation ofδz.

表1 玻片位置和軸向變形系數對橢球軸沿z軸負方向最大俘獲率的影響Table 1.The maximum negative trapping efficiency affected by different positions of coverslip and axial deformations.

表2 玻片位置和軸向變形系數對橢球軸向平衡位置的影響Table 2.The bead’s equilibrium position affected by different positions of coverslip and axial deformations.

表2展示了玻片不同位置對橢球的軸向平衡位置Zp的影響.當橢球軸向變形系數相同時,橢球微粒軸向平衡位置|Zcg|增大而向負方向偏移,這是因為球差增大,光線通過玻片匯聚的位置朝著負方向移動,從而導致平衡位置隨著光線匯聚位置的變化而變化.

4.3 討 論

四元數能方便地表述空間的旋轉.在光線追跡中,對于任意界面光線折射和反射,使用傳統歐氏幾何的方法求解空間光線的方位角尤為復雜,尤其不適于非球形的顆粒.采用矢量光線追跡,我們可以通過坐標旋轉將空間光線的折射和反射簡化到入射平面上處理,再通過旋轉坐標系的方式回到原坐標系中.目前我們采用四元數法的向量旋轉的方式,直接在空間中旋轉向量,避免了多次矩陣旋轉,簡化了光線追跡的困難.

四元數在表述光從光密介質到光疏介質或是光從光疏介質到光密介質時,其公式在表達形式上完全統一,在實際編程計算上更加簡潔方便.此外,四元數在光線追跡中不單可以用來求解折射光線與反射光線,還能對目標剛體進行旋轉.這意味著可以使用四元數的方式模擬任意擺放的非球形顆粒.

5 結 論

在光鑷的RO模型中,追跡光線向量在界面上的反射和折射是計算光阱力的最基本問題.在前人的計算過程中,更多地依賴球形微粒的對稱性求解光阱力,本文中提出四元數法旋轉入射光向量,可以計算任意界面入射光線的反射光線和折射光線,不依賴微粒的幾何特征,比以往的方法適用范圍更寬,使得光鑷的RO模型的應用面更廣泛.

[1]Schnitzer M J,Block S M 1997Nature388 386

[2]Qian H,Chen H,Yan J 2016Acta Phys.Sin.65 188706(in Chinese)[錢輝,陳虎,嚴潔2016物理學報 65 188706]

[3]Fazal F M,Block S M 2011Nature Photon.5 318

[4]Xia P,Zhou J H,Song X Y,Wu B,Liu X,Li D,Zhang S Y,Wang Z K,Yu H J,Ward T,Zhang J C,Li Y M,Wang X N,Chen Y,Guo Z,Yao X B 2014J.Mol.Cell Biol.6 240

[5]Ouyang H D,Wei M T 2010Annu.Rev.Phys.Chem.61 421

[6]Zhou J H,Ren H L,Cai J,Li Y M 2008Appl.Opt.47 6307

[7]Xu S H,Li Y M,Lou L R 2006Chin.Phys.15 1391

[8]Wright W H,Sonek GJ,Berns M W 1994Appl.Opt.33 1735

[9]Stilgoe A B,Nieminen T A,Knoner G,Heckenberg N R,Rubinsztein-Dunlop H 2008Opt.Express16 15039

[10]Gu Y Q,Gong Z,Lou L R,Li Y M 2007Appl.Laser27 98(in Chinese)[谷勇強,龔鏨,樓立人,李銀妹 2007應用激光27 98]

[11]Ashkin A 1992Biophys.J.61 569

[12]F?llman E,Axner O 2003Appl.Opt.42 3915

[13]Reihani S N S,Oddershede L B 2007Opt.Lett.32 1998

[14]Sidick E,Collins S D,Knoesen A 1997Appl.Opt.36 6423

[15]Bareil P B,Sheng Y,Chiou A 2006Opt.Express14 12503

[16]Zhou J H,Zhong M C,Wang Z Q,Li Y M 2012Opt.Express20 14928

[17]Kuipers J B 1999Geometry,Intergrability and Quantization(New York：Coral Press)pp127-143

[18]Xu F G 2012Physics with Quaternions(Beijing：Peking University Press)pp16-186(in Chinese)[許方官2012四元數物理學(北京：北京大學出版社)第16—186頁]

[19]Pletinckx D 1989Visual Comput.5 2

[20]Zhang R H,Jia H G,Chen T,Zhang Y 2008Opt.Precis.Eng.16 1965(in Chinese)[張榮輝,賈宏光,陳濤,張躍2008光學精密工程16 1965]

PACS:87.80.Cc,42.15.—i DOI:10.7498/aps.66.048701

Using quaternions to analyze the trapping force of an ellipsoidal bead?

Zhang Shu-He Liang Zhen?Zhou Jin-Hua?

(Department of Biomedical Engineering,Anhui Medical University,Hefei 230032,China)

3 October 2016;revised manuscript

17 October 2016)

In the ray-optics(RO)model of optical tweezers,tracing refractive and reflected rays with vectors play important roles in calculating the trapping forces.Traditional ray-tracing method with solid geometry,to some extent,is complicated in determining the orientations of those refractive and reflected rays according to spatial incident rays.It is difficult to calculate the trapping forces for irregular particles.In this paper,quaternion is proposed to rotate ray vectors for simplifying the traces of all kinds of spatial rays.Then,it is appropriate to calculate the trapping force of an ellipsoid bead.Based on the algorithm of quaternion and the convention between the interface normal and angular directions,the direction of normal always points from optically denser medium to thinner medium.The rotation axis is the cross product of the incident ray and the interface normal.And the positive angular direction can be determined by right-hand rule based on the orientation of the rotation axis.According to Snell’law,the rotation angle between the incident ray and refractive/reflected ray can be determined.The quaternion for rotation consists of rotation axis and angle.So the refractive and reflected rays are both determined by quaternions of incident ray and rotation based on rotation rules.Furthermore,the force on interface can also be calculated according to momentum changes of the photon before and after the interface refraction and reflection.The quaternion method is used to analyze the effects of coverslip position and deformation ratio on the trapping efficiency of ellipsoid particles.Our simulative results show that the lateral and axial trapping efficiencies are obviously affected by the deformation of the ellipsoid itself.No matter whether the bead deforms transversely or axially,the transverse and axial trapping efficiencies both become larger at a specific deformation.Meantime,the increase of the spherical aberration reduces the maximum axial trapping efficiency,and the equilibrium position of the bead becomes farther away from the center.Using quaternion method,the calculation of refractive lightvector can be simplified in comparison with by using the method of Euclidean geometry or transformation matrix.Theoretically,this quaternion can be used to trace rays on any irregular geometric surfaces.In conclusion,the method of quaternion can make ray tracing easier and extend the applications of RO model.

quaternion,ray tracing,ellipsoidal bead,trapping force

:87.80.Cc,42.15.—i

10.7498/aps.66.048701

?國家自然科學基金青年基金項目(批準號：31400943)、安徽高校自然科學研究重點項目(批準號：KJ2016A361)和安徽醫科大學博士科研資助基金(批準號：XJ201518)資助的課題.

?通信作者.E-mail：liangzhen@foxmail.com

?通信作者.E-mail：zhoujinhua@ahmu.edu.cn

*Project supported by the Young Scientists Fund of the National Natural Science Foundation of China(Grant No.31400943),the Key Project of Natural Science Foundation of the Anhui Higher Education Institutions,China(Grant No.KJ2016A361),and the Grants for Scientific Research of BSKY from Anhui Medical University,China(Grant No.XJ201518).

?Corresponding author.E-mail：liangzhen@foxmail.com

? Corresponding author.E-mail：zhoujinhua@ahmu.edu.cn