2324 J. Opt. Soc. Am. A/Vol. 23, No 9/September 2006 Trapping and binding of an arbitrary number of cylindrical particles in an in-plane electromagnetic field Tomasz M. Grzegorczyk, Brandon A Kemp, and Jin Au Kong Research Laboratory of Electronics, Massachusetts Institute of Technology, Cambridge, Massachusetts 02139 Received January 9, 2006; revised March 1, 2006; accepted March 16, 2006: posted April 6, 2006(Doc. ID 66668) The Mie theory and the Foldy-Lax multiple-scattering equations are applied to compute the scattered field of ane incidences. The Max. well stress tensor is then used to compute the force on each cylinder. Trapping and binding forces are studied as a fund of particle size, number, ittivity, and separation. Finally, the formulation is applied to a system of 20 particles, and the results show clear similarities with known experimental reports. The formu lation presented here extends the capabilities of modeling particle interaction and optical matter beyond the imple cases of the Rayleigh regime and two-particle systems. o 2006 Optical Society of America oCIS codes:290.4210,2602 20.7010 1 INTRODUCTION cal considerations have been limited to considering only The manipulation of small dielectric particles, typicall special cases where approximations could be used. Typi dielectric spheres, by electromagnetic waves wa first cally, a widely studied situation involves a single particle demonstrated in Ref. 1 and was immediately followed by in a Rayleigh regime: it has been shown that the force can important experimental observations such as optica be split into the gradient force and scattering force, 3-16 levitation and single-beam trapping. In such experi- respectively referring to the force due to the gradient of ments, the exerting force on the particles is due to the the field intensity and to the force due to the momentum wave to the particles, as was already predicted by James model approximates the scattered field of the particle by a C Maxwell. Levitation and trapping are associated with dipole radiation and therefore has to be used with caution two different physical origins: Levitation is due to the for larger particles or when multiple particle are consI scattering force of the impinging electromagnetic wave, ered in a close proximity. The scattering from large nd trapping is due to the gradient of the electric field in- spherical particles can be directly obtained from the Mie tensity acting on small particles. These two new effect theory; some researchers have used different approaches needed to be displaced without damage. 4 Optical binding study particles in the presence of a ground planel or th found applications in systems where minute particles such as the discrete dipole approximation(DDA) yas later reported,b as a third manifestation of the opti- multiple-multipole method to study cylindrical circular cal forces generated within a system of particles submit- and elliptical objects. However, the computation of the ted to an electromagnetic excitation. Optical binding, force on a system of multiple particles has not received which can become a dominant effect in the near field of much attention yet, despite some work on a particle in a the particles and can dictate their motion, can be under- complex environment such as a rough surface- and on stood as a secondary trapping effect, due not to the sole a system of two particles still in the Rayleigh regime22or incident field but to the scattered field from all the par- using the dda. The advantage of the approach used ticles in the system. Both trapping and binding have been Ref. 23 is the inherent advantage of DDA, namely, the recently experimentally verified on a collection of spheri- flexibility of studying particles of arbitrary shapes and cal polystyrene beads submitted to two different types of permittivities, but the disadvantages of the DDa are that incidence,, and a new trapping regime based on binding large particles require important computer resources forces has been reported in Ref 9 This is particularly true when the fields need to be com- The theoretical understanding of such puted in the close proximity of the particle boundary quires the computation of the force via either the Maxwell here a large number of dipoles have to be used in order stress tensor, or the lorentz force on bound currents to reduce the effects of surface discretization and charges(Ref. 12 and references herein). In both cases In the present paper, we present an exact theoretical e methods require the knowledge of the total electric model based on Mie theory to compute the force on an ar- and magnetic field around or within the particles due to bitrary number of dielectric particles under plane-wave the incident wave and the scattered waves from all the incidences. For the purpose of simplicity, we limit our- particles in the system. With an increasing number of selves to infinite dielectric cylinders and in-plane inci randomly positioned particles, the computation of the to- dences (i.e, with all the wave vectors perpendicular to the tal field quickly becomes nontrivial so that first theoreti- axes of the cylinders)with the electric field parallel to the 7/$15.00 o 2006 Optical Society of America
Trapping and binding of an arbitrary number of cylindrical particles in an in-plane electromagnetic field Tomasz M. Grzegorczyk, Brandon A. Kemp, and Jin Au Kong Research Laboratory of Electronics, Massachusetts Institute of Technology, Cambridge, Massachusetts 02139 Received January 9, 2006; revised March 1, 2006; accepted March 16, 2006; posted April 6, 2006 (Doc. ID 66668) The Mie theory and the Foldy–Lax multiple-scattering equations are applied to compute the scattered field of an arbitrary number of infinite dielectric cylinders of arbitrary size, subject to in-plane incidences. The Maxwell stress tensor is then used to compute the force on each cylinder. Trapping and binding forces are studied as a function of particle size, number, permittivity, and separation. Finally, the formulation is applied to a system of 20 particles, and the results show clear similarities with known experimental reports. The formulation presented here extends the capabilities of modeling particle interaction and optical matter beyond the simple cases of the Rayleigh regime and two-particle systems. © 2006 Optical Society of America OCIS codes: 290.4210, 260.2110, 020.7010. 1. INTRODUCTION The manipulation of small dielectric particles, typically dielectric spheres, by electromagnetic waves was first demonstrated in Ref. 1 and was immediately followed by important experimental observations such as optical levitation2 and single-beam trapping.3 In such experiments, the exerting force on the particles is due to the transfer of momentum from the incident electromagnetic wave to the particles, as was already predicted by James C. Maxwell. Levitation and trapping are associated with two different physical origins: Levitation is due to the scattering force of the impinging electromagnetic wave, and trapping is due to the gradient of the electric field intensity acting on small particles. These two new effects found applications in systems where minute particles needed to be displaced without damage.4 Optical binding was later reported5,6 as a third manifestation of the optical forces generated within a system of particles submitted to an electromagnetic excitation. Optical binding, which can become a dominant effect in the near field of the particles and can dictate their motion, can be understood as a secondary trapping effect, due not to the sole incident field but to the scattered field from all the particles in the system. Both trapping and binding have been recently experimentally verified on a collection of spherical polystyrene beads submitted to two different types of incidence,7,8 and a new trapping regime based on binding forces has been reported in Ref. 9. The theoretical understanding of such experiments requires the computation of the force via either the Maxwell stress tensor10,11 or the Lorentz force on bound currents and charges (Ref. 12 and references herein). In both cases, the methods require the knowledge of the total electric and magnetic field around or within the particles due to the incident wave and the scattered waves from all the particles in the system. With an increasing number of randomly positioned particles, the computation of the total field quickly becomes nontrivial so that first theoretical considerations have been limited to considering only special cases where approximations could be used. Typically, a widely studied situation involves a single particle in a Rayleigh regime: it has been shown that the force can be split into the gradient force and scattering force,13–16 respectively referring to the force due to the gradient of the field intensity and to the force due to the momentum of the electromagnetic wave. However, the Rayleigh model approximates the scattered field of the particle by a dipole radiation and therefore has to be used with caution for larger particles or when multiple particles are considered in a close proximity. The scattering from large spherical particles can be directly obtained from the Mie theory; some researchers have used different approaches such as the discrete dipole approximation17 (DDA) to study particles in the presence of a ground plane18 or the multiple-multipole method to study cylindrical circular and elliptical objects.19 However, the computation of the force on a system of multiple particles has not received much attention yet, despite some work on a particle in a complex environment such as a rough surface20,21 and on a system of two particles still in the Rayleigh regime22 or using the DDA.23 The advantage of the approach used in Ref. 23 is the inherent advantage of DDA, namely, the flexibility of studying particles of arbitrary shapes and permittivities, but the disadvantages of the DDA are that large particles require important computer resources. This is particularly true when the fields need to be computed in the close proximity of the particle boundary, where a large number of dipoles have to be used in order to reduce the effects of surface discretization. In the present paper, we present an exact theoretical model based on Mie theory to compute the force on an arbitrary number of dielectric particles under plane-wave incidences. For the purpose of simplicity, we limit ourselves to infinite dielectric cylinders and in-plane incidences (i.e., with all the wave vectors perpendicular to the axes of the cylinders) with the electric field parallel to the 2324 J. Opt. Soc. Am. A/Vol. 23, No. 9/September 2006 Grzegorczyk et al. 1084-7529/06/092324-7/$15.00 © 2006 Optical Society of America
Grzegorczyk et a VoL23,No9/September2006/J.Opt.Soc.Am.A2325 axis of the cylinders, but the theory is directly generaliz N,(k, p)=ikN(kp)ein, able to oblique incidences and spherical particles as well The interactions between the particles are computed us- where the prime denotes a derivative with respect to the ing the Foldy-Lax exact multiple-scattering equations, argument and where the Hankel functions are replaced yielding the expansion coefficients of the total field in the by the Bessel function for Rg Mn and RgNn. Note that ylindrical coordinate system(P, c, 2). The advantage of Eqs. (2)the coefficients an are unknown and are solved for the method is the ability to model an arbitrary number of upon applying the boundary condition on the total tan tion other than the mode truncation in the cylindrical cylinder. Before this is done, the fields inside the cyln o particles of arbitrary sizes, without additional approxima- gential electric and magnetic fields at the boundary of ti wave expansion. The disadvantages are that the particles also need to be expressed as in Eqs. (2)with another set of are required to be of a canonical shape(cylinders in our unknown coefficients. However, since we shall use in this case)and of identical permittivity. Since the number of paper the approach based on the Maxwell stress tensor articles as well as the number of plane-wave incidences the computation of the force on the cylinder due to the in can be arbitrary, both trapping and binding phenomena cident field requires the knowledge of the field outside the can be predicted, as we shall illustrate hereafter. We cylinder only, hence we do not write the field inside for should note that the limitations of in-plane incidence and brevity. When multiple particles are present, Eqs. (2) identical radius(which we shall use later)are easily need to be generalized so that the scattered field of a lifted. In addition, the case of an electric field perpendicu- given particle includes the effects from the incident field lar to the axis of the cylinders is not discussed here but and from the scattered field of all the other particles, ir could be studied in a very similar manner cluding all interactions. The self-consistent formulation proposed by Foldy2 and Lax26 is used here and can be written as 2. CYLINDRICAL WAVE EXPANSION AND FOLDY-LAX APPROACH TO MULTIPLE E(p)=En(p)+∑ECn(p) SCATTERING The problem we are considering is expressed as follows: a total of L infinite cylinders with their axes parallel to Equation(4)states that the final exciting field of particle are located at positions(x(,yo(e=l, .. L). For the sake g is equal to the incident wave plus the scattered waves to of simplicity, all the cylinders are assumed to be identical, particle q from all the other particles e except q itself. A of radius a, real relative permittivity ee, wavenumber ke, a generalization of Eqs. (2), Escat(p) is written as and nonmagnetic. A plane wave with electric field Einc(p)=iEoeik'p is incident onto the system with k=kk Egat(p)=> asaN,(k, p-Pa xk+yky, where a bold symbol denotes a vector and a at denotes a vector of norm 1. It is already well known where the cylindrical functions Nn have to be evaluated at that such incidence, along with the associated magnetic translated origins for each particle q located at pg field,can be decomposed onto a cylindrical coordinate sys- =(xa, y). The coefficients ana) are the generalized version of an and reduce to an in the case of a single particle. Their expression is obtained from Eine(p)=iEoek'p=> aRgN,(k, p), (1a) where Tn is the known T coefficient of a cylinder for the zxKE given polarization, Hine(p) ×k ke n(kca) Jn(ka) k Hocan hl ohere an=Eoine-indi/k. The scattered field in the regg g tside the cylinder is decomposed in a similar fashion B,=-ckH(ka) ,(ka)-ky(ka)()]- n=一x (7b) and where wn are solved from Hap)=:∑aM(,p), an+ T here n is the impedance of the background medium. Ir Eqs.(1)and(2) with ra being the angle between the x axis and the vec M,(k, p)=p-h(kp)eind-akHn(kp)eine, tor joining the centers or particles (e) and(g). equation (8)has been obtained upon applying the translational
axis of the cylinders, but the theory is directly generalizable to oblique incidences and spherical particles as well. The interactions between the particles are computed using the Foldy–Lax exact multiple-scattering equations, yielding the expansion coefficients of the total field in the cylindrical coordinate system ˆ ,ˆ ,zˆ. The advantage of the method is the ability to model an arbitrary number of particles of arbitrary sizes, without additional approximation other than the mode truncation in the cylindrical wave expansion. The disadvantages are that the particles are required to be of a canonical shape (cylinders in our case) and of identical permittivity. Since the number of particles as well as the number of plane-wave incidences can be arbitrary, both trapping and binding phenomena can be predicted, as we shall illustrate hereafter. We should note that the limitations of in-plane incidence and identical radius (which we shall use later) are easily lifted. In addition, the case of an electric field perpendicular to the axis of the cylinders is not discussed here but could be studied in a very similar manner. 2. CYLINDRICAL WAVE EXPANSION AND FOLDY–LAX APPROACH TO MULTIPLE SCATTERING The problem we are considering is expressed as follows: a total of L infinite cylinders with their axes parallel to zˆ are located at positions x ,y =1,...,L. For the sake of simplicity, all the cylinders are assumed to be identical, of radius a, real relative permittivity c, wavenumber kc, and nonmagnetic. A plane wave with electric field Einc=zˆE0 eik· is incident onto the system with k=kˆ k =xˆ kx+yˆ ky, where a bold symbol denotes a vector and a hat denotes a vector of norm 1. It is already well known that such incidence, along with the associated magnetic field, can be decomposed onto a cylindrical coordinate system as24 Einc = zˆE0 eik· = n=− + anRgNnk,, 1a Hinc = − zˆ kˆ zˆ kˆ E0 0 eik· , 1b where an=E0i n e−ini /k. The scattered field in the region outside the cylinder is decomposed in a similar fashion as Escat = n=− + an s Nnk,, 2a Hscat = 1 i n=− + an s Mnk,, 2b where is the impedance of the background medium. In Eqs. (1) and (2), Mnk, = ˆ in Hn 1 kein − ˆ kHn 1 kein , 3a Nnk, = zˆkNn 1 kein , 3b where the prime denotes a derivative with respect to the argument and where the Hankel functions are replaced by the Bessel function for RgMn and RgNn. Note that in Eqs. (2) the coefficients an s are unknown and are solved for upon applying the boundary condition on the total tangential electric and magnetic fields at the boundary of the cylinder. Before this is done, the fields inside the cylinder also need to be expressed as in Eqs. (2) with another set of unknown coefficients. However, since we shall use in this paper the approach based on the Maxwell stress tensor, the computation of the force on the cylinder due to the incident field requires the knowledge of the field outside the cylinder only; hence we do not write the field inside for brevity. When multiple particles are present, Eqs. (2) need to be generalized so that the scattered field of a given particle includes the effects from the incident field and from the scattered field of all the other particles, including all interactions. The self-consistent formulation proposed by Foldy25 and Lax26 is used here and can be written as24 Eexq = Einc +=1 lq L Escat . 4 Equation (4) states that the final exciting field of particle q is equal to the incident wave plus the scattered waves to particle q from all the other particles except q itself. As a generalization of Eqs. (2), Escat q is written as Escat q = n=− + an sq Nnk, − q, 5 where the cylindrical functions Nn have to be evaluated at translated origins for each particle q located at q =xq ,yq. The coefficients an sq are the generalized version of an s and reduce to an s in the case of a single particle. Their expression is obtained from an sq = Tnwn q , 6 where Tn is the known T coefficient of a cylinder for the given polarization, Tn = kc k Jnkca Hn 1 ka Bn − Jnka Hn 1 ka , 7a Bn = 2ik akc kHn 1 kaJnkca − kcJn kcaHn 1 ka−1, 7b and where wn q are solved from wn q = eik·qan + n=− + =1 lq L Hn−n 1 kl − qe−in−nlqTn wn , 8 with q being the angle between the xˆ axis and the vector joining the centers or particles and q. Equation (8) has been obtained upon applying the translational Grzegorczyk et al. Vol. 23, No. 9/September 2006/J. Opt. Soc. Am. A 2325
2326 J. Opt. Soc. Am. A/Vol. 23, No 9/September 2006 theorem of the Hankel function to each individual par The validation of the previous equations has been p ticle formed by comparing the predicted total field scattered by multiple particles with that obtained by the commercial Hh (elpe-P l)e,'p=> Jn(klPe-paDeinetg package CST MICROWAVE STUDIO. The comparison has been done with a number of particles, ranging from one to five submitted to one, two, and three incident plane waves. In XHn-n, (Alpp-Pale-i(n-n)pg.(9) all cases, the fields(E2, Hr, Hy) have been monitored and verified to agree well, both inside and outside the dielec In practice, the infinite sums are truncated at n tric cylinders. In addition, the force has been compared ax(n)(2+x-2+N), where N has to be properly chosen with the one obtained by the method presented in Ref 12 Ls a function of the number of particles in the system, in the way shown in Refs. 28 and 29, where bound cur- their size, permittivity, and separation, in order to achieve good convergence(values of N are specified in ticle. The agreement between the two calculation meth each case hereafter). The solution of the system of Eq (8) ods has been found to be excellent in all the situations requires the inversion of a square matrix of dimensions studied hereafter. a noticeable advantage of the Maxwell X(2N+1)=21L for N=10. For our purpose stress tensor approach, however, is that it requires fewer Gauss inversion is used, but more efficient algorithms can integration points, since the force is expressed as a line be implemented to improve the computation efficiency of integral instead of a surface integral. On the other hand the calculation. It should be noted that such matrix size is the method of bound currents and charges has the advan significantly smaller than the one obtained from a DDa tage of computing the local force density within the par method applied to particles of similar properties. Hence ticle, which is not accessible from the Maxwell stress ap- losing the ability of modeling arbitrary shapes and per- mittivities, we gain a considerable saving in computer memory and computational time. Once the w a) coefficients are solved for asia) are com- puted from Eq(6), and Eq) (p)is obtained from Eq (5). 3. EVIDENCE OF THE TRAPPING FORCE Finally, the total scattered field is expressed as Trapping of particles has been experimentally verified in Refs. 6 and 7 with a set of three incident gaussian beams E∞a(p)=∑Eq(p (10) The corresponding conclusions are here confirmed theo- retically on a simplified system of cylindrical particles submitted to three plane waves of identical amplitude The interior fields are obtained in a similar manner but Eo=1(V/m)propagating in the (xy) plane with angles are not required for the computation of the force on the f 7/2, 7/6, 11/6)(rad). The three incident plane waves cylinders in have a similar wavelength in free space of \=546 nm, and Knowing the scattered field, and therefore the total the particles are polystyrene cylinders of permittivity e fields (E, H)in outside all cylinders, we c 2.56E0 embedded in water(e=1.69Eo). The radius of the pute the time-averaged Lorentz force from the Maxwell's particles used in Ref 7 is 1 um, which is not a constraint stress tensor as in our case since size is not a limitation. However. the traps created by the three in-plane incident waves are F= del-Re((En)E*)+Re((Hn)H, maller than 1 um, as can be seen from the background pattern in Fig. 1, and therefore 1 um particles would feel the effect of multiple traps. To keep the particle size of the -E·E*n--(H·H*n same order as the traps, we resort to particles of radii =0.15A=819 nm and a=03A=1638mm. In Eq(11), Rel) denotes the real-part operator, n=i Y The trapping force on a single particle is shown in Fig where the background field pattern shows the interfer- he normal vector to the surface of the cylinder, and the ence of the three incident plane waves(only a relevant closed contour C surrounds the particle from the outside portion of the (xy) plane is shown, since the pattern is pe- thus requiring the knowledge of the fields outside the par- riodic in both directions) and the arrows show the force ticle only. In this case, the permittivity e and permeability acting on the particle as it sweeps the(xy) plane. The ar- u are those of the background medium. Equation(11)has rows in Fig. 1 show both the relative magnitude and di lready been extensively used in the literature but often rection of the force, where the starting point of each arrow applied to small particles. In such cases, gradient and corresponds to the location of the center of the cylinder at scattering forces can be separated, and a series of ap- which the computation is performed. It can immediately proximations can be used to simplify Eq (11)into expres- be seen that the trapping ability of the incident pattern sions with direct physical interpretations. However, in the significantly different for the two different sizes of pa case of large particles, or when many closely packed pal ticles: The high field intensity at(, y)=(0, 0)is a point of ticles are considered, no systematic approximation of Eq. stable equilibrium for a=0.15,, while it is a point of un 11)can be used. As we shall see hereafter, the results of stable equilibrium for a=0.3. It should be noted that a the force on large particles or on many particles do not al- similar behavior was found in Ref. 30 for a particle in the ways follow the intuition built from the study of small Rayleigh regime(a<X/20)as a function of its refractive
theorem of the Hankel function to each individual particle: Hn 1 k − peinlp = n=− + Jnk − qeinlq Hn−n 1 kp − qe−in−npq. 9 In practice, the infinite sums are truncated at N =maxn − +→−N +N, where N has to be properly chosen as a function of the number of particles in the system, their size, permittivity, and separation, in order to achieve good convergence (values of N are specified in each case hereafter). The solution of the system of Eq. (8) requires the inversion of a square matrix of dimensions L2N+1=21L for N=10. For our purpose, a direct Gauss inversion is used, but more efficient algorithms can be implemented to improve the computation efficiency of the calculation. It should be noted that such matrix size is significantly smaller than the one obtained from a DDA method applied to particles of similar properties. Hence, losing the ability of modeling arbitrary shapes and permittivities, we gain a considerable saving in computer memory and computational time. Once the wn q coefficients are solved for, an sq are computed from Eq. (6), and Escat q is obtained from Eq. (5). Finally, the total scattered field is expressed as Escat =q=1 L Escat q . 10 The interior fields are obtained in a similar manner but are not required for the computation of the force on the cylinders in our case. Knowing the scattered field, and therefore the total fields E,H in the region outside all cylinders, we compute the time-averaged Lorentz force from the Maxwell’s stress tensor as10,27 F = C d 2 Re E · nˆ E * + 2 Re H · nˆ Hˆ * − 4 E · E * nˆ − 4 H · H * nˆ . 11 In Eq. (11), Re{ } denotes the real-part operator, nˆ =ˆ is the normal vector to the surface of the cylinder, and the closed contour C surrounds the particle from the outside, thus requiring the knowledge of the fields outside the particle only. In this case, the permittivity and permeability are those of the background medium. Equation (11) has already been extensively used in the literature but often applied to small particles. In such cases, gradient and scattering forces can be separated, and a series of approximations can be used to simplify Eq. (11) into expressions with direct physical interpretations. However, in the case of large particles, or when many closely packed particles are considered, no systematic approximation of Eq. (11) can be used. As we shall see hereafter, the results of the force on large particles or on many particles do not always follow the intuition built from the study of small particles. The validation of the previous equations has been performed by comparing the predicted total field scattered by multiple particles with that obtained by the commercial package CST MICROWAVE STUDIO. The comparison has been done with a number of particles, ranging from one to five, submitted to one, two, and three incident plane waves. In all cases, the fields Ez ,Hx ,Hy have been monitored and verified to agree well, both inside and outside the dielectric cylinders. In addition, the force has been compared with the one obtained by the method presented in Ref. 12 in the way shown in Refs. 28 and 29, where bound currents and charges are used within the volume of the particle. The agreement between the two calculation methods has been found to be excellent in all the situations studied hereafter. A noticeable advantage of the Maxwell stress tensor approach, however, is that it requires fewer integration points, since the force is expressed as a line integral instead of a surface integral. On the other hand, the method of bound currents and charges has the advantage of computing the local force density within the particle, which is not accessible from the Maxwell stress approach. 3. EVIDENCE OF THE TRAPPING FORCE Trapping of particles has been experimentally verified in Refs. 6 and 7 with a set of three incident Gaussian beams. The corresponding conclusions are here confirmed theoretically on a simplified system of cylindrical particles submitted to three plane waves of identical amplitude E0=1 V/m propagating in the xy plane with angles /2,7/6,11/6 (rad). The three incident plane waves have a similar wavelength in free space of =546 nm, and the particles are polystyrene cylinders of permittivity c =2.560 embedded in water =1.690. The radius of the particles used in Ref. 7 is 1 m, which is not a constraint in our case since size is not a limitation. However, the traps created by the three in-plane incident waves are smaller than 1 m, as can be seen from the background pattern in Fig. 1, and therefore 1 m particles would feel the effect of multiple traps. To keep the particle size of the same order as the traps, we resort to particles of radii a =0.15 =81.9 nm and a=0.3 =163.8 nm. The trapping force on a single particle is shown in Fig. 1 where the background field pattern shows the interference of the three incident plane waves (only a relevant portion of the xy plane is shown, since the pattern is periodic in both directions) and the arrows show the force acting on the particle as it sweeps the xy plane. The arrows in Fig. 1 show both the relative magnitude and direction of the force, where the starting point of each arrow corresponds to the location of the center of the cylinder at which the computation is performed. It can immediately be seen that the trapping ability of the incident pattern is significantly different for the two different sizes of particles: The high field intensity at x,y=0,0 is a point of stable equilibrium for a=0.15 , while it is a point of unstable equilibrium for a=0.3 . It should be noted that a similar behavior was found in Ref. 30 for a particle in the Rayleigh regime a /20 as a function of its refractive index. 2326 J. Opt. Soc. Am. A/Vol. 23, No. 9/September 2006 Grzegorczyk et al.
Grzegorczyk et a Vol 23, No 9/September 2006/J. Opt. Soc. Am. A 2327 200-150-100 xm501001a0 Fig. 1. Force(represented by the black arrows) linder due to the interference of three plane waves(represented by the ackground pattern)of identical amplitude Eo=1 relength A=546 nm, and incident angles /2, 7/6, 11/6)(rad). Other pa- rameters are e=2.56∈0,∈=169∈0,N=5, and the of the cylinder is as indicated. The background pattern shows the absolute value of the electric field on a scale from 0(V/m)(black)to 3(V/m)(white)(a)a=0.15 (b)a=0.3. In Fig. 1(a), the force pattern is seen to closely obey a gradient law, which can be expected for small particles. In act, calculations have been performed for a series of sizes smaller than a=0.15, and all yielded a similar force pat- tern. The drastic difference with the pattern of Fig. 1(b) indicates the limit of the gradient approximation. With these properties(size and permittivity contrast), the par- ticle is already slightly larger than a trap and feels the nfluence of the periodic repetition of the high field inten sities. Interestingly, calculations with particles of sizes a 0.5A reveal a force pattern very similar to the one of Fig 1(a), while particles of a=0.75A yield a pattern similar to Fig. 1(b). In both ca wever, the particles are very Fig. 2. Two infinite cylinders of radius a and permittivity are arge compared with the trapping sites, and a simple gra- embedded in a background medium of permittivity e=1. 69Eo and dient law is not expected to remain valid. Finally, it subject to an incident field propagating in the y direction as worth mentioning that, for large particles, vortices are formed in the regions of low intensities: The particles are attracted toward these regions and are trapped following rection(due to the transfer or momentum from the elec a spiral motion before coming to rest tromagnetic wave to the particle) would be measured The results of Fig. 1 have been obtained with 11 cylin- Hence, any force in the x direction is due to interaction drical modes(N=5), which we have confirmed to be suffi effects, associated with binding phenomena cient to yield an error of less than 1% compared with our Figure 3 shows the i directed force on the right-hand reference computed with 81 modes (N=40) particle in Fig. 2 as function of relative distance and par- ticle size. Since the system is symmetric, the force on the 4. EVIDENCE OF BINDING FORCE left-hand particle is immediately obtained by symmetry. It can be seen first that the force is oscillatory as a func- Unlike the trapping force, which is immediately felt by a tion of the distance between the two particles. This obser- ingle particle in a properly arranged incident-field con- vation generalizes the result already presented in Ref 23 figuration, binding forces require multiple particles. In for very small particles. Second, it is seen that the principle, binding forces arise as soon as two particles are strength of the force does not depend on the particle size present: the scattered field of the first one modifies the in- in a trivial way: Although the strength is generally larger cident field on the second one, thus modifying the force with the larger particle the strength on the particle with that would have been felt by the second particle if it had a=0.3\ is smaller than on the particle with a=0. 15\. We been alone. In practice, however, the perturbation of the attribute this effect to the size of the particle compared incident field by the scattered field is often weak, and with the trapping site. As has been mentioned already, a pinding forces are easier to measure when a large number particle with a=0. 15A is approximately the same size as of particles scatter together? the trap, while a particle with a=0. 3\ is slightly larger al In theoretical calculations such as we present here ready (i.e, the particle feels a nonnegligible influence however, two particles are sufficient to show the binding from the neighboring traps), and a particle with a =0.5X is effect. The configuration, originally studied in Ref. 23 on influenced by many traps. The same effect that dictates very small particles (a=\/60), is depicted in Fig. 2: two the difference in the trapping force patterns of Fig. 1 die- particles are aligned in the x direction and submitted to tates here the difference between the binding forces. Fur- single incident plane wave propagating in the y direction. thermore, we have confirmed that for particles smaller If the two particles were isolated, only a force in the y than the traps, the binding force diminishes with the size
In Fig. 1(a), the force pattern is seen to closely obey a gradient law, which can be expected for small particles. In fact, calculations have been performed for a series of sizes smaller than a=0.15 , and all yielded a similar force pattern. The drastic difference with the pattern of Fig. 1(b) indicates the limit of the gradient approximation. With these properties (size and permittivity contrast), the particle is already slightly larger than a trap and feels the influence of the periodic repetition of the high field intensities. Interestingly, calculations with particles of sizes a =0.5 reveal a force pattern very similar to the one of Fig. 1(a), while particles of a=0.75 yield a pattern similar to Fig. 1(b). In both cases, however, the particles are very large compared with the trapping sites, and a simple gradient law is not expected to remain valid. Finally, it is worth mentioning that, for large particles, vortices are formed in the regions of low intensities: The particles are attracted toward these regions and are trapped following a spiral motion before coming to rest. The results of Fig. 1 have been obtained with 11 cylindrical modes N=5, which we have confirmed to be suffi- cient to yield an error of less than 1% compared with our reference computed with 81 modes N=40. 4. EVIDENCE OF BINDING FORCE Unlike the trapping force, which is immediately felt by a single particle in a properly arranged incident-field con- figuration, binding forces require multiple particles. In principle, binding forces arise as soon as two particles are present: the scattered field of the first one modifies the incident field on the second one, thus modifying the force that would have been felt by the second particle if it had been alone. In practice, however, the perturbation of the incident field by the scattered field is often weak, and binding forces are easier to measure when a large number of particles scatter together.7 In theoretical calculations such as we present here, however, two particles are sufficient to show the binding effect. The configuration, originally studied in Ref. 23 on very small particles a /60, is depicted in Fig. 2: two particles are aligned in the xˆ direction and submitted to a single incident plane wave propagating in the yˆ direction. If the two particles were isolated, only a force in the yˆ direction (due to the transfer or momentum from the electromagnetic wave to the particle) would be measured. Hence, any force in the xˆ direction is due to interaction effects, associated with binding phenomena. Figure 3 shows the xˆ directed force on the right-hand particle in Fig. 2 as function of relative distance and particle size. Since the system is symmetric, the force on the left-hand particle is immediately obtained by symmetry. It can be seen first that the force is oscillatory as a function of the distance between the two particles. This observation generalizes the result already presented in Ref. 23 for very small particles. Second, it is seen that the strength of the force does not depend on the particle size in a trivial way: Although the strength is generally larger with the larger particle, the strength on the particle with a=0.3 is smaller than on the particle with a=0.15 . We attribute this effect to the size of the particle compared with the trapping site. As has been mentioned already, a particle with a=0.15 is approximately the same size as the trap, while a particle with a=0.3 is slightly larger already (i.e., the particle feels a nonnegligible influence from the neighboring traps), and a particle with a=0.5 is influenced by many traps. The same effect that dictates the difference in the trapping force patterns of Fig. 1 dictates here the difference between the binding forces. Furthermore, we have confirmed that for particles smaller than the traps, the binding force diminishes with the size Fig. 1. Force (represented by the black arrows) on a single cylinder due to the interference of three plane waves (represented by the background pattern) of identical amplitude E0=1 V/m, wavelength =546 nm, and incident angles /2,7/6,11/6 (rad). Other parameters are c=2.560, =1.690, N=5, and the radius a of the cylinder is as indicated. The background pattern shows the absolute value of the electric field on a scale from 0 V/m (black) to 3 V/m (white) (a) a=0.15 , (b) a=0.3 . Fig. 2. Two infinite cylinders of radius a and permittivity c are embedded in a background medium of permittivity =1.690 and subject to an incident field propagating in the yˆ direction as indicated. Grzegorczyk et al. Vol. 23, No. 9/September 2006/J. Opt. Soc. Am. A 2327
2328 J. Opt. Soc. Am. A/Vol. 23, No 9/September 2006 of the particle. Such result could be confirmed by an nalysis in the Rayleigh regime but is outside the scope of this paper. The effect of the permittivity of the cylinders on the binding force is illustrated in Fig. 4. As expected the larger the permittivity contrast between the particle and the background, the stronger the binding. This effect is here confirmed far from the Rayleigh regime and can be directly understood from the Born approximation: the scattered field is directly proportional to (Ec-e) and is therefore reduced in a low-permittivity contrast system Finally, we analyze the amplitude of the force exerted n the right-hand particle as a function of the number of particles on the left. The results for one particle are iden- tical to the ones presented in Fig 4 with the correspond- ing parameters, while the cases of two and four particles x [nm are illustrated in the insets of Fig. 5. The computation is Fig. 5. Force in the i direction for varying number of cylinders performed here with N=20 to ensure good convergence of right). The case of one cylinder corresponds to Fig. 2, the cases of two and four cylinders correspond to the situations shown in the linder case. Two interesting conclusions can be drawn insets. The positions of the fixed cylinders are (x,y) from these results. First, the binding force is stronger the cases of two and four particles, respectively. The parameters two particles. This directly justifies the experimental veri- are A=546 nm, =2.5660.6-1.6960,@=0.3A, and N=20 fication that binding phenomena become measurable when multiple particles are scattering together. Second the binding force is almost identical in the three-particle system and the five-particle system. Hence, this result puts a limitation on the previous conclusion, which is re- lated to the incident field: in the configuration studied here, the scattered fields from two and four particles on the left along the x axis are almost identical. This can be directly understood from a generalization of the Mie theory, which indicates that a large particle has a strong forward scattering and low side scattering. The clusters of two and four particles behave in some sense like a single ∧AA八 large particle, having different forward-scattering fields and similar side-scattering fields. Hence, the similarity in the binding force would obviously not hold for ar x-directed incidence, in which case the particle on the right would be in the forward-scattering region of the clusters and would be affected differently by their scat tered fields Fig 3. Force x direction for the configuration of Fig. 2 as unction of the relative distance between parti eters are A=546 nm, e-2.56Eo, E=1. 69Eo, N=10, and radius a as 5 RANDOM ARRANGEMENT OF 20 ndicated in the labels PARTICLES within the limitation of in-plane incident waves, we re- produce here part of the experimental conclusions of Ref. 7. Toward thi return to the of three plane-wave incidences like those shown in Fig. 1, where all the parameters are as indicated in the caption of the figure. The size of the cylinders is taken to be a=0. 151 since we have seen in Fig. 1(a)that the incident interfer- ence pattern yields good traps. Within a region of space of ApUaNe (5a X 5a), we randomly position 20 particles and compute the force on each due to the incident field and the scat tered field from all the particles in the system. The posi tion of each particle is then adjusted in space by ar amount proportional to the force acting on it. Mathemati cally, if we denote by (x(, yo and vre,fy the position and x [nmI Fig 4. Force in the x direction for the configuration of Fig. 2 as force on particle t at iteration i, respectively, the updated function of the relative distance cles. The positions at iteration(i+1)are obtained from eters are A=546 nm, E=1.69Eo, a= N=10. The relative permittivity of the cylinders is as in the labels (12a)
of the particle. Such result could be confirmed by an analysis in the Rayleigh regime but is outside the scope of this paper. The effect of the permittivity of the cylinders on the binding force is illustrated in Fig. 4. As expected, the larger the permittivity contrast between the particle and the background, the stronger the binding. This effect is here confirmed far from the Rayleigh regime and can be directly understood from the Born approximation: the scattered field is directly proportional to c− and is therefore reduced in a low-permittivity contrast system. Finally, we analyze the amplitude of the force exerted on the right-hand particle as a function of the number of particles on the left. The results for one particle are identical to the ones presented in Fig. 4 with the corresponding parameters, while the cases of two and four particles are illustrated in the insets of Fig. 5. The computation is performed here with N=20 to ensure good convergence of the fields even for closely packed cylinders as in the fourcylinder case. Two interesting conclusions can be drawn from these results. First, the binding force is stronger within the system of three particles than in the system of two particles. This directly justifies the experimental veri- fication that binding phenomena become measurable when multiple particles are scattering together.7 Second, the binding force is almost identical in the three-particle system and the five-particle system. Hence, this result puts a limitation on the previous conclusion, which is related to the incident field: in the configuration studied here, the scattered fields from two and four particles on the left along the xˆ axis are almost identical. This can be directly understood from a generalization of the Mie theory, which indicates that a large particle has a strong forward scattering and low side scattering. The clusters of two and four particles behave in some sense like a single large particle, having different forward-scattering fields and similar side-scattering fields. Hence, the similarity in the binding force would obviously not hold for an xˆ-directed incidence, in which case the particle on the right would be in the forward-scattering region of the clusters and would be affected differently by their scattered fields. 5. RANDOM ARRANGEMENT OF 20 PARTICLES Within the limitation of in-plane incident waves, we reproduce here part of the experimental conclusions of Ref. 7. Toward this purpose, we return to the case of three plane-wave incidences like those shown in Fig. 1, where all the parameters are as indicated in the caption of the figure. The size of the cylinders is taken to be a=0.15 since we have seen in Fig. 1(a) that the incident interference pattern yields good traps. Within a region of space of 5a5a, we randomly position 20 particles and compute the force on each due to the incident field and the scattered field from all the particles in the system. The position of each particle is then adjusted in space by an amount proportional to the force acting on it. Mathematically, if we denote by x i ,y i and f x i ,f y i the position and force on particle at iteration i, respectively, the updated positions at iteration i+1 are obtained from xl i+1 = xl i + f xl i , 12a Fig. 3. Force in the xˆ direction for the configuration of Fig. 2 as a function of the relative distance between particles. The parameters are =546 nm, c=2.560, =1.690, N=10, and radius a as indicated in the labels. Fig. 4. Force in the xˆ direction for the configuration of Fig. 2 as a function of the relative distance between particles. The parameters are =546 nm, =1.690, a=0.3 , and N=10. The relative permittivity of the cylinders is as indicated in the labels. Fig. 5. Force in the xˆ direction for varying number of cylinders on the left (the force and the positions are for the particle on the right). The case of one cylinder corresponds to Fig. 2, the cases of two and four cylinders correspond to the situations shown in the insets. The positions of the fixed cylinders are x,y =0, ±170 nm and x,y= −170 nm, ±170 nm,0, ±170 nm in the cases of two and four particles, respectively. The parameters are =546 nm, c=2.560, =1.690, a=0.3 , and N=20. 2328 J. Opt. Soc. Am. A/Vol. 23, No. 9/September 2006 Grzegorczyk et al.
Grzegorczyk et a ol. 23, No 9/September 2006/J. Opt. Soc. Am. 2329 (12b) where a is arbitrarily chosen such that the particles are not moved too much or too little between the first and sec. ond iterations (we typically tolerate a motion within a fraction of the particle radius). At the new set of positions the forces acting on each particle are computed again, and the process is reiterated until the correction terms on the 鵝翻 positions are negligible(typically less than 1% of the cur- rent position). It should be emphasized that the motion of the particles resulting from this iterative process is not trictly derived from the equations of motions of mechan cs but is merely a visual tool to estimate the evolution of the particles as a function of the forces acting on them when strong damping is present as in the case of a water background < The experiment with 20 particles was repeated mul- ole times, each time with a different set of initial posi- tions. An example of initial positions is given in Table 1 and illustrated in Fig. 6(a). Because of the large number of particles in the relatively small constrained surface the particles are often closely packed over single trapping sites and may compete for a site until an equilibrium is reached. The equilibrium is attained either when the par- ticles rearrange themselves to occupy separate traps, as shown in Fig. 6(b)(with the final positions given in Table 1)or when the total field is such that the effect of neigh- Fig. 6. Positions of 20 dielectric cylinders and field distributions boring traps on a conglomerate of particles produces no s cases:(a) random initial position in a three-plane-wave interference pattern (incident force, as shown close to the origin in Fig. 6(d). In both fi- field shown);(b) organized final position due to trapping and nal positions, it should also be noted that the particles do binding forces (incident field shown); (e)same as case(b)but with not always settle in regions of high field intensity, where the total field shown; (d)organized final position corresponding particle [see Fig. 1(a)]. This is a direct manifestation of a=0.15, set of initial positions different from that in case (a) the traps are predicted to be the strongest for a single 6∈a,∈=169 and N=10. The background patterns show the absolute the binding forces: The scattered field by all the particles value of the electric field (either incident or total field )on a scale from 0(black) to 3 V/m(white) Table 1. Initial Positions(xi, yi) and Final Positions(af, ye) of 20 Particles including all their interactions. modifies the incident-field ThreePlane- Wave Interference pattern pattern and thus modifies the location of the traps. The background pattern shown in Fig. 6(b), representing the yi(nm) incident field only, therefore only serves as a reference to 14891 59.30 00.08 120.94 see where the traps should have been if the scattered field -377.76 -38895 -410.69 -533.34 between the particles had been ignored. Obviously, the fi- -380.12 nal positions do not correspond to the traps due to the in- -335.74 7.45 -486.75 cident field, and therefore the gradient force approxi 106.53 396.91 18.50 tion used in the modeling of small particles cannot be 286.33 used in this case. Figure 6(c)shows the same case as Fig 278.84 301.18 6(b)with the background pattern corresponding to the to- 17949 -427.16 tal field instead of the incident field only. Interestingly, it 276.50 36289 219.12 346.56 is seen that the particles come to rest in regions of high 27925 91.55 204.00 111.49 total field intensity, indicating that for such size the gr dient force due to the total field is dominant Finally, we have confirmed that if the dielectric con 70.09 231 trast between the background medium and the particles is very low(typically, e-e=0.1), the traps due to the in 151.11 349.56 198.72 cident field single particle are still located in the high-field regions as in Fig. 1(a), and their strength is sig- 0637 nificantly increased. The immediate consequence is that 22498 art 416.51 -199.41 168.19 2597 within an area of (5a X5a) the particles evolve to almost 397.59 249.79 633.45 exactly occupy the positions of high field intensities due to he incident waves only, similar to what has been experi The corresponding positions are shown in Figs. 6(a)and 6(b). mentally reported in Refs. 7 and 31. This indicates ne
yl i+1 = yl i + f yl i , 12b where is arbitrarily chosen such that the particles are not moved too much or too little between the first and second iterations (we typically tolerate a motion within a fraction of the particle radius). At the new set of positions, the forces acting on each particle are computed again, and the process is reiterated until the correction terms on the positions are negligible (typically less than 1% of the current position). It should be emphasized that the motion of the particles resulting from this iterative process is not strictly derived from the equations of motions of mechanics but is merely a visual tool to estimate the evolution of the particles as a function of the forces acting on them when strong damping is present as in the case of a water background. The experiment with 20 particles was repeated multiple times, each time with a different set of initial positions. An example of initial positions is given in Table 1 and illustrated in Fig. 6(a). Because of the large number of particles in the relatively small constrained surface, the particles are often closely packed over single trapping sites and may compete for a site until an equilibrium is reached. The equilibrium is attained either when the particles rearrange themselves to occupy separate traps, as shown in Fig. 6(b) (with the final positions given in Table 1) or when the total field is such that the effect of neighboring traps on a conglomerate of particles produces no force, as shown close to the origin in Fig. 6(d). In both fi- nal positions, it should also be noted that the particles do not always settle in regions of high field intensity, where the traps are predicted to be the strongest for a single particle [see Fig. 1(a)]. This is a direct manifestation of the binding forces: The scattered field by all the particles, including all their interactions, modifies the incident-field pattern and thus modifies the location of the traps. The background pattern shown in Fig. 6(b), representing the incident field only, therefore only serves as a reference to see where the traps should have been if the scattered field between the particles had been ignored. Obviously, the fi- nal positions do not correspond to the traps due to the incident field, and therefore the gradient force approximation used in the modeling of small particles cannot be used in this case. Figure 6(c) shows the same case as Fig. 6(b) with the background pattern corresponding to the total field instead of the incident field only. Interestingly, it is seen that the particles come to rest in regions of high total field intensity, indicating that for such size the gradient force due to the total field is dominant. Finally, we have confirmed that if the dielectric contrast between the background medium and the particles is very low (typically, c−0.1), the traps due to the incident field on a single particle are still located in the high-field regions as in Fig. 1(a), and their strength is significantly increased. The immediate consequence is that when a set of 20 such particles are randomly positioned within an area of 5a5a the particles evolve to almost exactly occupy the positions of high field intensities due to the incident waves only, similar to what has been experimentally reported in Refs. 7 and 31. This indicates not Table 1. Initial Positions „xi,yi… and Final Positions „xf,yf… of 20 Particles in a Three-Plane-Wave Interference Patterna xi (nm) yi (nm) xf (nm) yf (nm) −148.91 −59.30 −200.08 −120.94 −377.76 −388.95 −410.69 −533.34 296.79 −380.12 179.12 −356.27 30.74 −335.74 −7.45 −486.75 −106.53 396.91 18.50 454.02 286.33 −85.47 398.75 16.19 278.84 301.18 219.68 374.49 −272.76 179.49 −427.16 245.14 −276.50 362.89 −219.12 346.56 279.25 91.55 204.00 111.49 −287.76 −231.91 −388.75 −230.97 42.83 17.06 −5.79 −2.15 70.09 231.64 −7.01 235.09 −310.92 3.13 −409.60 19.53 −151.11 −349.56 −198.72 −402.18 10.59 −162.72 2.95 −247.57 −101.12 118.56 −206.37 118.08 388.51 −224.98 416.51 −199.41 168.19 −235.86 193.02 −125.97 141.00 397.59 249.79 633.45 a The corresponding positions are shown in Figs. 6a and 6b. Fig. 6. Positions of 20 dielectric cylinders and field distributions (shown in the background) for various cases: (a) random initial position in a three-plane-wave interference pattern (incident field shown); (b) organized final position due to trapping and binding forces (incident field shown); (c) same as case (b) but with the total field shown; (d) organized final position corresponding to another set of initial positions different from that in case (a). In all cases, the parameters are =546 nm, c=2560, =1690, a=0.15 , and N=10. The background patterns show the absolute value of the electric field (either incident or total field) on a scale from 0 (black) to 3 V/m (white). Grzegorczyk et al. Vol. 23, No. 9/September 2006/J. Opt. Soc. Am. A 2329
2330 J. Opt. Soc. Am. A/Vol 23, No 9/September 2006 only the good efficiency of the traps but also that the seat- 10. J. Stratton, Electromagnetic Theory (McGraw-Hill, ring is weak and does not perturb the incident waves 1941) significantly. For this type of dielectric contrast, a born 11. J. A. Kong, Electromagnetic Waue (EMW,2000 approximation can be used to compute the field in the sys- 12. A. R. Zakharian, M. Mansuript "Radiation pressure and the distri tem and confirm these cor inclusions force in dielectric media, "Opt. 321-2336 13. A. Ashkin, Trapping of by resonance radiation 6. CONCLUSIONS sure, Phys. Rev. Lett. 40, 729-732(1978) 14. J. P Gordon. "Radiation forces and momenta in dielectric We have demonstrated in this paper that both trapping media, "Phys. Rev. A 8, 14-21(1973 nd binding forces in a complex arrangement of particles 15. P C Chaumet and M. Nieto-Vesperinas,"Time-averaged an be predicted without approximations on particle size total force on a dipolar sphere in an electromagnetic field, permittivity, or separation. The theoretical model, based Opt.Lett.25,1065-1067(2001) on an extension of Mie theory to cylindrical particles com- 16. R. Arias-Gonzalez and M. Nieto- Vesperinas, " Optical forces on small particles: attractive and repulsive nature bined with the Foldy-Lax multiple-scattering equations d plasmon-resonance conditions, J. Opt. Soc. Am.A nd the Maxwell stress tensor, has been shown to predict 01-1209(2003 some recent experimental results reasonably well with 17. E. M. Purcell and C. R. Pennypacker, "Scattering imited requirements on computer memory. Such model- ing capability represents a step forward in the under- 18. P. C. Chaumet and M. Nieto-Vesperinas, "Coupled dipole standing of optical tweezers, optical matter, and other method determination of the electromagnetic force systems where large particles are manipulated by radia particle over a fat dielectric substrate, " Phys. Rev. B 61 tion pressure 14119-14127(2000) 19. C. Rockstuhl and H. P. Herzig, " Rigorous lied to the analysis of the optical ACKNOWLEDGMENTS J. Opt. A, P 921-931(2004 It is a pleasure to acknowledge many stimulating discus- 20. A. Madrazo and M. Nieto- Vesperinas, Scattering of ns with J -M. Fournier. This work is sponsored by the Department of the U. S. Air Force under Air Force con- tract FA8721-05-C-0002 and by the Chinese Nation 2t onducting plane,"J.Opt. Soc. Am. A 12, 1298-1309 A. Madrazo and M. Nieto-Vesperinas, " Surface structure Foundation under contracts 60371010 and 60531020 Opinions, interpretations, conclusions, and recommenda- lectromagnetic waves from a r in front of a tions are those of the author and are not necessarily en- onducting grating: theory for the reflection photon eling microscope, J. Opt. Soc. Am. A 13, Corresponding author T. M. Grzegorczyk can be 22. F. Depasse and J.-M. Vigoureux, Optical binding force reached by e-mail at tomasz@mit. edu. 23. P C Chaumet and M. Nieto- Vesperinas, "Optical binding of REFERENCES articles with or without the of a fat dielectr urface, Phys. Rev. B 64, 035422(2001) 1. A. Ashkin. "Acceleration and tra of particles by 24. L. Tsang, J. Kong, K. Ding, and C. Ao, Scattering of Electromagnetic Waves: Numerical Simulations 2. A. Ashkin and M. Dziedzic (1971 5. LL Foldy, " The multiple scattering of waves, "Phys. Rev. 3. A. Ashkin and J. M. Dziedzic, "Optical levitation in high 26. M. Lax, "Multiple scattering of waves. Il. The effective field acuum, " Appl. Phys. Lett. 28, 333-335(1976) 4. A. Ashkin, "Applications of laser radiation pressure, 27. M. Lester and M. Nieto-Vesperinas, "Optical forces on Science210,1081-1088(1980) 5. M. M. Burns. J.M. Fournier. and A. Go microparticles in an evanescent laser field, Opt Lett. 24 "Optical binding, "Phys. Rev. Lett. 63, 1233-123 urns. J 28. B A. Kemp, T. M. Grzegorczyk, and J. A. Kong, "Ab initio nd binding in intense study of the radiation media, "Opt tical fields. "Science 249, 749-754(1990) 92809291(2005 7. J-M. Fournier. G 29. B. A. Kemp, T M. Grzegorczyk, and J. A. Kong, Lorentz R. Proe.SPI5514,309-317 Waves Appl.20,827-839(2006) (2004) 30. P. Zemanek, V. Karasek, and A. Sasso, "Optical forces 8. A Casaburi, G. Pesce, P Zemanek, and A Sasso, Two- and 240,401-415(2004) Commun.251,393-404(2005 T. M. Grzegorczyk, B. A. Kemp, and J. A. Kong, "Stable 31. and, Rousaiat e. Assem b inga quet, sco ichanat ices ia optical trapping based on optical binding forces, "Phys. Rev various optical schemes, in Proc. SPIE 5930, 238-247 Let.6,113903(2006) 2005
only the good efficiency of the traps but also that the scattering is weak and does not perturb the incident waves significantly. For this type of dielectric contrast, a Born approximation can be used to compute the field in the system and confirm these conclusions. 6. CONCLUSIONS We have demonstrated in this paper that both trapping and binding forces in a complex arrangement of particles can be predicted without approximations on particle size, permittivity, or separation. The theoretical model, based on an extension of Mie theory to cylindrical particles combined with the Foldy–Lax multiple-scattering equations and the Maxwell stress tensor, has been shown to predict some recent experimental results reasonably well with limited requirements on computer memory. Such modeling capability represents a step forward in the understanding of optical tweezers, optical matter, and other systems where large particles are manipulated by radiation pressure. ACKNOWLEDGMENTS It is a pleasure to acknowledge many stimulating discussions with J.-M. Fournier. This work is sponsored by the Department of the U.S. Air Force under Air Force contract FA8721-05-C-0002 and by the Chinese National Foundation under contracts 60371010 and 60531020. Opinions, interpretations, conclusions, and recommendations are those of the author and are not necessarily endorsed by the U. S. Government. Corresponding author T. M. Grzegorczyk can be reached by e-mail at tomasz@mit.edu. REFERENCES 1. A. Ashkin, “Acceleration and trapping of particles by radiation pressure,” Phys. Rev. Lett. 24, 156–159 (1970). 2. A. Ashkin and J. M. Dziedzic, “Optical levitation by radiation pressure,” Appl. Phys. Lett. 19, 283–285 (1971). 3. A. Ashkin and J. M. Dziedzic, “Optical levitation in high vacuum,” Appl. Phys. Lett. 28, 333–335 (1976). 4. A. Ashkin, “Applications of laser radiation pressure,” Science 210, 1081–1088 (1980). 5. M. M. Burns, J.-M. Fournier, and J. A. Golovchenko, “Optical binding,” Phys. Rev. Lett. 63, 1233–1236 (1989). 6. M. M. Burns, J.-M. Fournier, and J. A. Golovchenko, “Optical matter: crystallization and binding in intense optical fields,” Science 249, 749–754 (1990). 7. J.-M. Fournier, G. Boer, G. Delacrétaz, P. Jacquot, J. Rohner, and R. Salathé, “Building optical matter with binding and trapping forces,” in Proc. SPIE 5514, 309–317 (2004). 8. A. Casaburi, G. Pesce, P. Zemánek, and A. Sasso, “Two- and three-beam interferometric optical tweezers,” Opt. Commun. 251, 393–404 (2005). 9. T. M. Grzegorczyk, B. A. Kemp, and J. A. Kong, “Stable optical trapping based on optical binding forces,” Phys. Rev. Lett. 96, 113903 (2006). 10. J. Stratton, Electromagnetic Theory (McGraw-Hill, 1941). 11. J. A. Kong, Electromagnetic Wave Theory (EMW, 2000). 12. A. R. Zakharian, M. Mansuripur, and J. V. Moloney, “Radiation pressure and the distribution of electromagnetic force in dielectric media,” Opt. Express 13, 2321–2336 (2005). 13. A. Ashkin, “Trapping of atoms by resonance radiation pressure,” Phys. Rev. Lett. 40, 729–732 (1978). 14. J. P. Gordon, “Radiation forces and momenta in dielectric media,” Phys. Rev. A 8, 14–21 (1973). 15. P. C. Chaumet and M. Nieto-Vesperinas, “Time-averaged total force on a dipolar sphere in an electromagnetic field,” Opt. Lett. 25, 1065–1067 (2001). 16. J. R. Arias-González and M. Nieto-Vesperinas, “Optical forces on small particles: attractive and repulsive nature and plasmon-resonance conditions,” J. Opt. Soc. Am. A 20, 1201–1209 (2003). 17. E. M. Purcell and C. R. Pennypacker, “Scattering and absorption of light by nonspherical dielectric grains,” Astrophys. J. 186, 705–714 (1973). 18. P. C. Chaumet and M. Nieto-Vesperinas, “Coupled dipole method determination of the electromagnetic force on a particle over a flat dielectric substrate,” Phys. Rev. B 61, 14119–14127 (2000). 19. C. Rockstuhl and H. P. Herzig, “Rigorous diffraction theory applied to the analysis of the optical force on elliptical nano- and micro-cylinders,” J. Opt. A, Pure Appl. Opt. 6, 921–931 (2004). 20. A. Madrazo and M. Nieto-Vesperinas, “Scattering of electromagnetic waves from a cylinder in front of a conducting plane,” J. Opt. Soc. Am. A 12, 1298–1309 (1995). 21. A. Madrazo and M. Nieto-Vesperinas, “Surface structure and polariton interactions in the scattering of electromagnetic waves from a cylinder in front of a conducting grating: theory for the reflection photon scanning tunneling microscope,” J. Opt. Soc. Am. A 13, 785–795 (1996). 22. F. Depasse and J.-M. Vigoureux, “Optical binding force between two Rayleigh particles,” J. Phys. D 27, 914–919 (1994). 23. P. C. Chaumet and M. Nieto-Vesperinas, “Optical binding of particles with or without the presence of a flat dielectric surface,” Phys. Rev. B 64, 035422 (2001). 24. L. Tsang, J. Kong, K. Ding, and C. Ao, Scattering of Electromagnetic Waves: Numerical Simulations (Wiley, 2000). 25. L. L. Foldy, “The multiple scattering of waves,” Phys. Rev. 67, 107–119 (1945). 26. M. Lax, “Multiple scattering of waves. II. The effective field in dense systems,” Phys. Rev. 85, 261–269 (1952). 27. M. Lester and M. Nieto-Vesperinas, “Optical forces on microparticles in an evanescent laser field,” Opt. Lett. 24, 936–938 (1999). 28. B. A. Kemp, T. M. Grzegorczyk, and J. A. Kong, “Ab initio study of the radiation pressure on dielectric and magnetic media,” Opt. Express 13, 9280–9291 (2005). 29. B. A. Kemp, T. M. Grzegorczyk, and J. A. Kong, “Lorentz force on dielectric and magnetic particles,” J. Electromagn. Waves Appl. 20, 827–839 (2006). 30. P. Zemánek, V. Karásek, and A. Sasso, “Optical forces acting on Rayleigh particle placed into interference field,” Opt. Commun. 240, 401–415 (2004). 31. J.-M. Fournier, J. Rohner, P. Jacquot, R. Johann, S. Mias, and R. Salathé, “Assembling mesoscopic particles by various optical schemes,” in Proc. SPIE 5930, 238–247 (2005). 2330 J. Opt. Soc. Am. A/Vol. 23, No. 9/September 2006 Grzegorczyk et al