Computation of radiation pressure force exerted on arbitrary shaped homogeneous particles by high-order Bessel vortex beams using MLFMA

Minglin Yang; Yueqian Wu; Kuan Fang Ren; Xinqing Sheng

doi:10.1364/OE.24.027979

1. Introduction

A particle illuminated by a beam of light experiences the radiation pressure force (RPF), and under the action of the RPF, the particle could be trapped and moved to a prescribed location [1, 2]. The computation of RPF is of growing interest due to its importance in real applications, such as biological cell trapping [3], microfluidic devices [4], and laser based measurement techniques [5]. Optical force calculation plays an important role in the analysis and understanding of optical micromanipulation.

Traditional optical tweezers use focused Gaussian beams as trapping light sources. Theories for computing RPF have already been developed by using different approaches. When the particle is much smaller than the wavelength, the Rayleigh regime is concerned and RPF can be calculated by using the Rayleigh theory [6]. For particles of size not very small compared to the wavelength, the rigorous electromagnetic theory such as the generalized Lorenz-Mie theory (GLMT) has been developed [7] for studies of RPF or torque exerted on a particle of simple shape such as homogeneous sphere [8, 9], multilayered sphere [10], spheroid [11] and a infinite long cylinder of circular cross section [12]. For complex shaped particles, numerical techniques are possible ways, and researches based on the discrete dipole approximation (DDA) [13–15], the T-matrix method [16, 17] and the finite difference time domain method (FDTD) [18–20] for RPF on non-spherical particles are reported. However, the T-matrix method usually relies on central expansions of the electromagnetic field in VSWF [21]; DDA and FDTD are not efficient for large particles, especially for high refractive index particles [22]. The multilevel fast multipole algorithm (MLFMA) enhanced surface integral equation (SIE) method has shown its great power for solving electromagnetic scattering problems [23–25]. Recently, we have presented MLFMA for computation of RPF exerted on large arbitrary shaped homogeneous particles by a Gaussian beam [26, 27]. This algorithm has been further improved recently [28].

On the other hand, due to special characteristics of nondiffraction and self reconstruction, the Bessel beams have attracted wide attention in various fields, including optical trapping and manipulation, optical acceleration and nonlinear optics [29–31]. It appears to be a dramatic alternative to using Gaussian beams potentially in some particular scenarios [32]. Some studies have been carried out on the RPF exerted on simple shaped particles by Bessel beams using the Rayleigh model [33], the geometrical optics [34, 35] or the rigorous electromagnetic theory [36–38]. Nevertheless, these studies mainly focused on the cases of spherical particles. Trapping of non-spherical particles is both very different and much less mastered than that of spheres. To the best of our knowledge, the RPF exerted on large arbitrary shaped particles by a high-order Bessel vortex beam (HOBVB) has not been studied yet.

This paper is devoted to the computation of RPF exerted on large arbitrary shaped homogeneous particles by HOBVB. The arbitrarily shaped particles are modeled by using surface triangular patches and the SIEs are discretized by following procedure of the method of moments (MoM). The MLFMA is utilized to speed up matrix-vector multiplication in the iterative solution step of the final matrix equation system. Then the radiation force is computed by integrating the dot product of the outwardly directed normal unit vector and the Maxwell stress tensor over a spherical surface tightly enclosed the particle. We use the accurately computed near region electromagnetic fields and the accurate vector expressions of HOBVB electromagnetic field components instead of the widely used far-field approximation. Some numerical results are given to show accuracy and capability of the proposed approach.

The body of this paper is organized as follows: Section 2 gives in detail the vector description of the high-order Bessel vortex beama. Section 3 gives a brief introduction to the MLFMA enhanced SIE method. Section 4 describes how to compute RPF using near fields. In Section 5, the numerical accuracy of MLFMA is validated and the RPF exerted on large non-spherical particles by a HOBVB are computed and analyzed. The last section is devoted to the conclusions.

2. Description of the high-order Bessel vortex beam

We first give the mathematical expressions of the incident HOBVB. We introduce here two cartesian coordinate systems: one is the particle coordinate systems (O′uvw), the other is the beam coordinate systems (Oxyz). A Bessel beam with its electric field polarized in the u direction is assumed to propagate along the w axis in its own Cartesian coordinate system (O′uvw). To provide satisfactory results for the case when the size of the central spot of the beam is much larger than the wavelength, we adopt the vector expressions presented by Mitri to describe the incident Bessel beam in the beam coordinate system (O′uvw) [39], which take into account the vector nature of electromagnetic waves. The electromagnetic field components in O′uvw are:

\begin{array}{l} E_{u}^{i} = \frac{1}{2} E_{0} [(1 + \frac{k_{w}}{k} - \frac{k_{r}^{2} u^{2}}{k^{2} r^{2}} + \frac{n (n - 1) {(u - j v)}^{2}}{k^{2} r^{4}}) J_{n} (k_{r} r) \\ - \frac{k_{r} (v^{2} - u^{2} - 2 j n u v)}{k^{2} r^{3}} J_{n + 1} (k_{r} r)] e^{j (k_{w} w + n φ)} \end{array}

\begin{array}{l} E_{v}^{i} = \frac{1}{2} E_{0} u v [\frac{2 k_{r} + 2 i n \frac{v^{2} - u^{2}}{u v}}{k^{2} r^{3}} J_{n + 1} (k_{r} r) \\ - \frac{n (n - 1) (2 + j \frac{u^{2} - v^{2}}{u v}) - k_{r}^{2} r^{2}}{k^{2} r^{4}} J_{n} (k_{r} r)] e^{j (k_{w} w + n φ)} \end{array}

E_{w}^{i} = \frac{1}{2} j E_{0} \frac{u}{k r} (1 + \frac{k_{w}}{k}) [\frac{n (1 - j \frac{v}{u})}{r} J_{n} (k_{r} r) - k_{r} J_{n + 1} (k_{r} r)] e^{j (k_{w} w + n φ)}

\begin{array}{l} H_{u}^{i} = \frac{1}{2} H_{0} u v [- \frac{n (n - 1) (2 + j \frac{u^{2} - v^{2}}{u v}) - k_{r}^{2} r^{2}}{k^{2} r^{4}} J_{n} (k_{r} r) \\ + \frac{2 k_{r} + 2 i n \frac{v^{2} - u^{2}}{u v}}{k^{2} r^{3}} J_{n + 1} (k_{r} r)] e^{j (k_{w} w + n φ)} \end{array}

\begin{array}{l} H_{v}^{i} = \frac{1}{2} H_{0} [(1 + \frac{k_{w}}{k} - \frac{k_{r}^{2} v^{2}}{k^{2} r^{2}} + \frac{n (n - 1) {(v + j u)}^{2}}{k^{2} r^{4}}) J_{n} (k_{r} r) \\ - \frac{k_{r} (u^{2} - v^{2} + 2 j n u v)}{k^{2} r^{3}} J_{n + 1} (k_{r} r)] e^{j (k_{w} w + n φ)} \end{array}

H_{w}^{i} = \frac{1}{2} j H_{0} \frac{v}{k r} (1 + \frac{k_{w}}{k}) [\frac{n (1 + j \frac{u}{v})}{r} J_{n +} (k_{r} r) - k_{r} J_{n + 1} (k_{r} r)] e^{j (k_{w} w + n φ)}

where

r = \sqrt{u^{2} + v^{2}}

, k_w = k cos ϕ with ϕ being the half-cone angle, and φ = tan⁻¹(v/u),

j = \sqrt{- 1}

the imaginary unit. The parameters E₀ and H₀ are the amplitudes of the electric field and magnetic field. The J_n is the first kind cylindrical Bessel function of the nth orders.

To give readers a better understanding, before carrying out numerical experiments, we give in Fig. 1 the magnitude plots for the components of the electric fields for a zero-order (n = 0) Bessel beam (ZOBB), as well as the first-order (n = 1) and second-order vortex Bessel beam. All the beams have the half-cone angle ϕ = 20°.

Fig. 1 Theoretical magnitude cross-sectional profiles for the electric field components of the zero-order Bessel beam and the HOBVB with different orders. The units in the (u, v) plane are in λ. The half-cone angle of the Bessel beam is set to ϕ = 20°. Fig. 1(a)–1(c) Zero-order. Fig. 1(d)–1(f) First-order. Fig. 1(g)–1(i) Second-order.

Download Full Size | PDF

In practice, we need the incident beam expressed in the particle coordinate system (Oxyz). The relation between the two systems can be obtained by a translation and three rotations according to Euler angles [40, 41]. As shown in Fig. 2 the coordinate of the beam center in the particle system are (x₀, y₀, z₀) and the three rotation Euler angles are α, β, γ.

Fig. 2 Schematic of arbitrary shaped homogeneous particle illuminated by a shaped beam and definition of Euler angles.

Download Full Size | PDF

The coordinates of a point in the particle system x, y, z can be expressed as function of the coordinates in the beam system u, v, w:

[\begin{array}{l} x - x_{0} \\ y - y_{0} \\ z - z_{0} \end{array}] = A [\begin{array}{l} u \\ v \\ w \end{array}]

where A is the transformation matrix:

A = [\begin{array}{l} a_{11} & a_{12} & a_{13} \\ a_{21} & a_{22} & a_{23} \\ a_{31} & a_{32} & a_{33} \end{array}]

and its elements are defined by

\begin{array}{l} a_{11} = \cos α \cos β \cos γ - \sin α \sin γ \\ a_{21} = - \cos α \cos β \sin γ - \sin α \cos γ \\ a_{31} = \cos α \sin β \\ a_{21} = \sin α \cos β \cos γ + \cos α \sin γ \\ a_{22} = - \sin α \cos β \sin γ + \cos α \cos γ \\ a_{23} = \sin α \sin β \\ a_{31} = - \sin β \cos γ \\ a_{32} = - \sin β \sin γ \\ a_{33} = \cos β \end{array}

Similarly, the electric and magnetic field components of the incident beam in the two coordinate systems are related by:

\begin{array}{l} [\begin{array}{l} E_{x}^{i} \\ E_{y}^{i} \\ E_{z}^{i} \end{array}] = A [\begin{array}{l} E_{u}^{i} \\ E_{v}^{i} \\ E_{w}^{i} \end{array}], & [\begin{array}{l} H_{x}^{i} \\ H_{y}^{i} \\ H_{z}^{i} \end{array}] = A [\begin{array}{l} H_{u}^{i} \\ H_{v}^{i} \\ H_{w}^{i} \end{array}] \end{array}

3. Combined tangential formulation with MLFMA

We give in this section a brief description of the formulation of the problem with MLFMA. The integral equation approach limits the discretization of the unknowns to the surface of the object. The boundary S of the dielectric object is taken as the equivalent surface, as shown in Fig. 2, with filed components of the incident HOBVB denoted as (Eⁱ, Hⁱ) and the equivalent electric and magnetic current sources as (J, M). Using the equivalence principle or the vector Green’s theorem, we can formulate a set of integral equations to calculate the electric and magnetic fields (E, H) in terms of equivalent electric and magnetic currents (J, M) on the boundary surface S of the dielectric body: the electric field integral equation outside the dielectric body (EFIE-O), the magnetic field integral equation outside the dielectric body (MFIE-O), the electric field integral equation inside the dielectric body (EFIE-I), and the magnetic field integral equation inside the dielectric body (MFIE-I) [24]. Different combinations of the four basic equations construct various forms of combined field integral equations [25, 42], among which the Combined tangential formulation (CTF) achieve best accuracy, especially when dealing with objects with sharp edges, corners or large refractive index. The equation form of CTF can be expressed as follow:

{\begin{matrix} Z_{1}^{- 1} \hat{t} \cdot (EFIE-O) + Z_{2}^{- 1} \hat{t} \cdot (EFIE-I) \\ Z_{1} \hat{t} \cdot (MFIE-O) + Z_{2} \hat{t} \cdot (MFIE-I) \end{matrix}

where

\hat{t}

is the tangential vector at a point on the surface, Z_l = (μ_l/ε_l)^1/2 denote the wave impedances in medium l with l = 1, 2 respectively for medium outside and inside object. Equation (11) can be discretized by first expanding (J, M) with the Rao, Wilton and Glisson(RWG) vector basis functions g_i as [43]:

\begin{matrix} J = \sum_{i = 1}^{N_{s}} g_{i} J_{i} & M = \sum_{i = 1}^{N_{s}} g_{i} M_{i} \end{matrix}

where N_s denotes the total number of edges on S. Applying g_i as the trial functions for the tangential field equations and following procedure of the well known method of moment (MoM), a complete matrix equation system can be obtained. After solving the matrix equation system to obtain (J, M) on the boundary surface, fields in any place can be gotten. The scattered fields E^s and H^s at any point can be obtained by:

E^{s} = Z_{1} L_{1} (J) - K_{1} (M)

H^{s} = 1 / Z_{1} L_{1} (M) - K_{1} (J)

Due to the computation and storage complexities in the order O(N²), with N the number of unknowns, the MoM is limited to dealing with small particles. To circumvent this limitation, various acceleration methods, such as the fast Fourier transform (FFT), the adaptive integral method (AIM) and the MLFMA [23, 44, 45] are adopted. Among these methods, the MLFMA has been well developed to reduce the time and memory complexities to O(N log N) [23,24,46]. In MLFMA, interactions between the basis and testing elements can be divided into two types: near interactions and far interactions. The near interactions are computed in the same way as in the conventional MoM. For far interactions, they are calculated in a group-by-group manner consisting of three stages called aggregation, translation, and disaggregation. The detailed description of MLFMA for solving the CTF is given in [24] and is not repeated here.

4. Computation of radiation pressure force

When an arbitrarily shaped particle is illuminated by a shaped beam, the radiation force exerted on the particle can be determined by integrating the dot product of the surface normal $\hat{n}$ and the Maxwell stress tensor $\overset{\leftrightarrow}{T}$ over a surface enclosing the particle [47]:

F = \int_{S_{v}} \overset{\leftrightarrow}{T} (r) \cdot \hat{n} d s

where

\begin{matrix} \overset{\leftrightarrow}{T} (r) = \frac{1}{2} Re [ε_{1} E (r) E^{*} (r) + μ_{1} H (r) H^{*} (r) \\ - \frac{1}{2} (ε_{1} {| E (r) |}^{2} + μ_{1} {| H (r) |}^{2}) \overset{\leftrightarrow}{I}] \end{matrix}

is the time average Maxwell stress tensor, the asterisk * indicates for conjugate, and E(r) and H(r) are total electromagnetic fields:

E (r) = E^{s} (r) + E^{i} (r)

H (r) = H^{s} (r) + H^{i} (r)

If we choose a virtual sphere of radius r_s tightly enclosing the particle with its center located in the scattering object, Eq. (15) can be written as [47]

\begin{array}{l} F = \frac{1}{4} \int_{0}^{2 π} \int_{0}^{π} Re [ε_{1} ({| E_{r} |}^{2} - {| E_{θ} |}^{2} - {| E_{ϕ} |}^{2}) \\ + μ_{1} ({| H_{r} |}^{2} - {| H_{θ} |}^{2} - {| H_{ϕ} |}^{2}) e_{r} + 2 (ε_{1} E_{r} E_{θ}^{*} + μ_{1} H_{r} H_{θ}^{*}) e_{θ} \\ + 2 (ε_{1} E_{r} E_{ϕ}^{*} + μ_{1} H_{r} H_{ϕ}^{*}) e_{ϕ}] r_{s}^{2} \sin θ d θ d ϕ \end{array}

where the electric and magnetic field components are near-fields evaluated on the surface of the virtual sphere.

5. Numerical results and discussions

In this section, a series of numerical experiments are performed. All the computations are performed on a computer platform having 2 Intel X5650 2.66 GHz CPUs with 6 cores for each CPU, 96 GB memory. Mesh density in the numerical realization is set to about 0.08 − 0.1λ. The generalized minimal residual method (GMRES) iteration solver is employed for solving the final matrix equation system with a restart number of 100. The HOBVB is incident on a hard polystyrene particle of refractive index m = 1.573 surrounded by air (m = 1.0) or a biconcave cell-like particle of refractive index m = 1.4 embedded in water (m = 1.3). In all the numerical examples except the last biconcave cell-like particles, we assumed the wavelength of the beams is 0.785 μm. In all the calculations, we set $I_{0} = E_{0}^{2} / 2 Z_{1}$ to be unit.

We first validate our code by comparison with the results of GLMT. Since the Bessel beam has characteristics of nondiffraction, moving away from the center of particle along the z axis has no influence on the RPF. A spherical particle with radius 1 μm is computed. The particle is illuminated by a first order Bessel beam with a half-cone angle of ϕ = 20°, and the beam center is moving along the y axis of the particle system. In such case, we have α = β = γ = 0° and x₀ = z₀ = 0. Figure 3 shows the radiation pressure force as a function of y₀ calculated by MLFMA and GLMT. It is found that the agreement between the two methods is excellent. We observe a nonzero RPF component in the x direction in Fig. 3(a), while for a zero-order Bessel beam it always remains zero. Such force is caused by the ‘vortex’ field phase characteristic of HOBVB. Due to this force, the annular rotation of particles were able to be realized. It can be canceled out using another counter propagating HOBVB with the same parameters. In Fig. 3(b), at the position y₀ = 0, a positive y₀ causes a positive traverse force along the y direction, and vice versa. It means at the beam center, the particle can be trapped in y direction. Due to its ring structure, other radial equilibrium positions can be observed, showing Bessel beams have the ability of trapping simultaneously high refractive index particles in their bright rings. We can also observe a force pushing the particle in the beam propagation direction (+z), which is the same as a Gaussian beam and not given here.

Fig. 3 Comparison of RPF components versus beam center location along the y axis computed by GLMT and MLFMA for a sphere with radius of 1 μm illuminated by a first order Bessel beam of wavelength 0.785 μm. The half-cone angle of the Bessel beam is ϕ = 20°. The refractive index is m = 1.573. The incident beam propagates along z axis.

Download Full Size | PDF

A zeroth-order Bessel beam is termed as a non-vortex beam and it has a very narrow nondiffracting central core (bright central spot), as in Fig. 1(a)). However, a HOBVB propagates over a characteristic length without spreading (dark central region), as shown in Fig. 1(d) and 1(g). The size of the central spot highly depends on the half-cone angle. The smaller half-cone angle, the larger the size of the central spot (bright for zero order and dark for first order), and vice versa. To study the influence of the half-cone angle, we set ϕ = 5°, ϕ = 10°, ϕ = 20° respectively, and the other parameters remain unchanged. The computed RPF components are shown in Fig. 4. Again nonzero traverse forces along x direction are observed in all the cases, as in Fig. 4(a). In Fig. 4(b), when ϕ = 5° and ϕ = 10°, near the beam axis, a positive y₀ causes a negative traverse force along the y direction. It means that the beam will drive the particle off the beam center, and the particle can not be trapped near the beam axis. In all these three cases, the particle experiences a force pushing the particle in the beam propagation direction (+z), but when ϕ = 5° the RPF curve for the first order beam has two maxima locates almost in the middle of the first bright ring, as shown in Fig. 4(c). This is because when ϕ = 5° radius of the first bright ring is larger than the sphere.

Fig. 4 RPFs versus beam center location along the y axis computed by MLFMA for a sphere with radius of 1 μm illuminated by a first order Bessel beam with half-cone angle ϕ = 5°, ϕ = 10° and ϕ = 20°. Other parameters are the same as in Fig. 3.

Download Full Size | PDF

Figure 5 presents the influence of the order of the beam on the RPF. For this case, we fix ϕ = 20° and change the order of HOBVB to n = 1, 2, 3 respectively. The RPFs versus beam center location along the y axis with different orders are given. In theory, if we fix ϕ as a constant, radius of the first bright ring increases with the increment of the beam order. Hence in Fig. 5(b), when order n = 2 and n = 3, positive y₀ results in negative F_y, and vice versa. Moreover, we can see in Fig. 5(c), when n = 3, the radius of the first bright ring is larger than the sphere. Hence the RPF curve of F_y for the third order beam has two maxima.

Fig. 5 Comparison of RPF versus beam center location along the y axis for a sphere illuminated by HOBVB of wavelength 0.785 μm with different orders. The half-cone angle of the Bessel beam is ϕ = 20°. The refractive index is m = 1.573. Other parameters are the same as in Fig. 3.

Download Full Size | PDF

Next, we study the case of non-spherical particles. We first study the influence of incidence angle for RPF exerted on spheroidal particles (a = 3 μm, b = 3 μm, c = 12 μm) by a HOBVB of wavelength 0.785 μm with different orders. The half-cone angle of the Bessel beam is set to ϕ = 5°. We fix (x₀ = y₀ = z₀ = 0), set α = γ = 0° and let the HOBVB incident on a spheroid polystyrene particle (m = 1.573) at different incidence angles β. In such case, the axis of the incident beam remains always in the xz plane and the center of the HOBVB is always coincident with the center of the particle, as shown in Fig. 6(a). The three components of RPF are depicted in Fig. 6. In this calculation, the whole surface of the particle is discretized into 212415 edges with total of 424830 unknowns. This computation costs about 6.2GB memory and 1600 seconds wall clock time with 12 threads for computing each RPF point. Again, the spheroid experiences nonzero force component F_y along y direction, as shown in Fig. 6(c), while F_y remains always zero for the Gaussian beam incidence case [26]. The maximum value of F_x and F_z decreases with the increment of beam order, as shown in Fig. 6(b) and Fig. 6(d). Moreover, depending on the radius of the first bright ring, the incident angle β where the particle experiences the maximum F_z increases with increment of the beam order.

Fig. 6 RPF on an spheroidal (a = 3 μm, b = 3 μm, c = 12 μm) polystyrene particle (m = 1.573) illuminated by a HOBVB of wavelength 0.785 μm with different orders. The half-cone angle of the Bessel beam is set to ϕ = 5°. The incident beam rotates around the center of the particle (x₀ = y₀ = z₀ = 0, α = γ = 0°).

Download Full Size | PDF

The proposed approach can be applied to compute RPF of any shaped particle. As an example, the RPFs exerted on biconcave cell-like particles in water versus beam center position along x axis by a HOBVB of the first order are calculated. The incident wavelength of the beam in vacuum is 0.5145 μm. The half-cone angle is set to ϕ = 10°. The surface function we employed to describe the biconcave cell-like particle is given as [48]:

r (θ_{A}, ϕ_{A}) = a \sin^{q} θ_{A} + b

We suppose all the biconcave cell-like particles have the same volume of a sphere with a radius of 2.66μm. Three different geometry parameters are taken into consideration: q = 2 (a = 2.1 μm, b = 1.12 μm), q = 5 (a = 3 μm, b = 0.75 μm) and q = 9 (a = 3.8 μm, b = 0.41 μm). The cross sections of biconcave cell-like particles with different q are shown in Fig. 7(a). The RPF calculated with different geometry parameters are shown in Fig. 7. We can find from Fig. 7(b) that, the profile of x components of RPF for q = 2 and q = 5 are quite similar, positive x₀ leads to a positive F_x which push the particle toward the beam axis. However, for q = 9, positive x₀ leads to a negative F_x which push the particle off the beam axis. In consequence, when the particle of q = 2 or q = 5 is off but near the beam axis, the RPF pushes it toward the beam axis but along a spiral trajectory. Therefore, the stable trapping at the center is possible, while the beam axis is not a stable position for a particle of q = 9. Furthermore, the x component of RPF may change the sign when the particle is off axis but not too far (about 2μm for the case in Fig. 7), while the direction of y component changes even when the particle is far from the beam axis. These properties are very different from a Gaussian beam. The F_y is about two orders of magnitude smaller than F_x and F_z because refractive index of the particle (m = 1.4) is very close to the background media (m = 1.3). All the three biconcave particles experience a positive force F_z pushing them along the beam propagating direction. Therefore, to obtain a 3D trap this force should be compensated by the gravitation or another beam.

Fig. 7 RPFs on the three biconcave cell-like particles (m = 1.4) in water (m = 1.3) versus beam center position along x axis. The half-cone angle of the first order HOBVB is ϕ = 10°.

Download Full Size | PDF

6. Conclusion

In this paper, we generalize the MLFMA enhanced SIE method to include the case of illumination by a high-order vortex Bessel beam and apply it to study RPFs exerted on large non-spherical particles. The incident HOBVB is described by vector expressions in terms of the electric and magnetic fields. The whole outer surface of the particle is discretized into triangle patches. The MLFMA is utilized to efficiently solve the resulting matrix equations. Then RPF is computed by vector flux of the Maxwell’s stress tensor over a spherical surface tightly enclosing the particle. Good agreements between the computed RPF exerted on a hard polystyrene particle and those from GLMT validate the accuracy and capability of the proposed approach. The influence of beam order, half-cone angle, beam center position, and the incident angles of HOBVB on the RPF are analyzed. The numerical results of RPFs exerted on a spheroidal particle and a biconcave cell-like particle with different shapes are presented, showing capability of the proposed method for computing RPFs on complex shaped particles. All the numerical results prove that, a particle in a HOBVB experiences a force in the direction perpendicularly to the symmetric plane due to the vortex characteristic of the beam which may push some kinds of particles (biconcave cell-like cell of q = 2 and q = 5 for example) toward the trapping position along a spiral trajectory, but not straightforward as in the Gaussian beam. Moreover, a HOBVB has the ability of trapping particles in other positions in addition to the beam axis due to its ring structure. These results are expected to provide useful insights into the RPFs exerted on complex-shaped particles by a HOBVB.

Funding

National Basic Research Program (973) (61320602, 61327301); 111 Project of China (B14010); National Natural Science Foundation of China (NSFC) (61421001); French National Research Agency (ANR-13-BS090008-01).

Acknowledgments

We would thank Renxian Li for providing the radiation pressure force of a sphere calculated by GLMT.

References and links

1. A. Ashkin, “Acceleration and trapping of particles by radiation pressure,” Phys. Rev. Lett. 24(4), 156 (1970). [CrossRef]

2. A. Ashkin, J. M. Dziedzic, J. E. Bjorkholm, and S. Chu, “Observation of a single-beam gradient force optical trap for dielectric particles,” Opt. Lett. 11(5), 288–290 (1986). [CrossRef] [PubMed]

3. K. Svoboda and S. M. Block, “Biological applications of optical forces,” Annu. Rev. Bioph. Biom. 23(1), 247–285 (1994). [CrossRef]

4. N. K. Metzger, M. Mazilu, L. Kelemen, P. Ormos, and K. Dholakia, “Observation and simulation of an optically driven micromotor,” J. Opt. 13(4), 044018 (2011). [CrossRef]

5. R. M. P. Doornbos, M. Schaeffer, A. G. Hoekstra, P. M. A. Sloot, B. G. de Grooth, and J. Greve, “Elastic light-scattering measurements of single biological cells in an optical trap,” Appl. Opt. 35(4), 729–734 (1996). [CrossRef] [PubMed]

6. Y. Harada and T. Asakura, “Radiation forces on a dielectric sphere in the Rayleigh scattering regime,” Opt. Commun. 124(5–6), 529–541 (1996). [CrossRef]

7. G. Gouesbet, B. Maheu, and G. Grehan, “Light scattering from a sphere arbitrarily located in a Gaussian beam, using a Bromwich formulation,” J. Opt. Soc. Am. A 5(9), 1427–1443 (1988). [CrossRef]

8. K. F. Ren, G. Grehan, and G. Gouesbet, “Radiation pressure forces exerted on a particle arbitrarily located in a Gaussian beam by using the generalized Lorenz-Mie theory, and associated resonance effects,” Opt. Commun. 108(4–6), 343–354 (1994). [CrossRef]

9. H. Polaert, G. Grehan, and G. Gouesbet, “Improved standard beams with application to reverse radiation pressure,” Appl. Opt. 37(12), 2435–2440 (1998). [CrossRef]

10. H. Polaert, G. Grehan, and G. Gouesbet, “Forces and torques exerted on a multilayered spherical particle by a focused Gaussian beam,” Opt. Commun. 155(1–3), 169–179 (1998). [CrossRef]

11. F. Xu, K. F. Ren, G. Gouesbet, X. Cai, and G. Grehan, “Theoretical prediction of radiation pressure force exerted on a spheroid by an arbitrarily shaped beam,” Phys. Rev. E 75(2), 026613 (2007). [CrossRef]

12. V. V. Kotlyar and A. G. Nalimov, “Analytical expression for radiation forces on a dielectric cylinder illuminated by a cylindrical Gaussian beam,” Opt. Express 14(13), 6316–6321 (2006). [CrossRef] [PubMed]

13. B. T. Draine and J. C. Weingartner, “Radiative torques on interstellar grains I. superthermal spin-up,” Astrophys. J. 470, 551 (1996). [CrossRef]

14. A. G. Hoekstra, M. Frijlink, L. B. F. M. Waters, and P. M. A. Sloot, “Radiation forces in the discrete-dipole approximation,” J. Opt. Soc. Am. A 18(8), 1944–1953 (2001). [CrossRef]

15. S. H. Simpson and S. Hanna, “Holographic optical trapping of microrods and nanowires,” J. Opt. Soc. Am. A 27(6), 1255–1264 (2010). [CrossRef]

16. T. A. Nieminen, H. Rubinsztein-Dunlop, N. R. Heckenberg, and A. I. Bishop, “Numerical modelling of optical trapping,” Comput. Phys. Commun. 142(1), 468–471 (2001). [CrossRef]

17. S. H. Simpson and S. Hanna, “Numerical calculation of interparticle forces rising in association with holographic assembly,” J. Opt. Soc. Am. A 23(6), 1419–1431 (2006). [CrossRef]

18. W. L. Collett, C. A. Ventrice, and S. M. Mahajan, “Electromagnetic wave technique to determine radiation torque on micromachines driven by light,” Appl. Phys. Lett. 82(16), 2730–2732 (2003). [CrossRef]

19. M. M. A. Zakharian and J. Moloney, “Radiation pressure on a dielectric wedge,” Opt. Express 13(6), 2064–2074 (2005). [CrossRef] [PubMed]

20. R. C. Gauthier, “Computation of the optical trapping force using an FDTD based technique,” Opt. Express 13(10), 3707–3718 (2005). [CrossRef] [PubMed]

21. T. A. Nieminen, H. Rubinsztein-Dunlop, and N. Heckenberg, “Calculation of the T-matrix: general considerations and application of the point-matching method,” J. Quant. Spectrosc. Radiat. Transf. 79–80(2), 1019–1029 (2003). [CrossRef]

22. M. A. Yurkin and A. G. Hoekstraa, “The discrete dipole approximation: An overview and recent developments,” J. Quant. Spectrosc. Radiat. Transf. 106(1–3), 558–589 (2007). [CrossRef]

23. J. M. Song and W. C. Chew, “Multilevel fast multipole algorithm for solving combined field integral equations of electromagnetic scattering,” Micro. Opt. Tech. Lett. 10(1), 14–19 (1995). [CrossRef]

24. X. Q. Sheng, J. M. Jin, J. M. Song, W. C. Chew, and C. C. Lu, “Solution of combined-field integral equation using multilevel fast multipole algorithm for scattering by homogeneous bodies,” IEEE Trans. Antennas Propag. 46(11), 1718–1726 (1998). [CrossRef]

25. Ö. Ergül and L. Gürel, “Comparison of integral-equation formulations for the fast and accurate solution of scattering problems involving dielectric objects with the multilevel fast multipole algorithm,” IEEE Trans. Antennas Propag. 57(1), 176–187 (2009). [CrossRef]

26. M. L. Yang, K. F. Ren, M. J. Gou, and X. Q. Sheng, “Computation of radiation pressure force on arbitrary shaped homogenous particles by multilevel fast multipole algorithm,” Opt. Lett. 38(11), 1784–1786 (2013). [CrossRef] [PubMed]

27. M. L. Yang, K. F. Ren, Y. Q. Wu, and X. Q. Sheng, “Computation of stress on the surface of a soft homogeneous arbitrarily shaped particle,” Phys. Rev. E 89(4), 043310 (2014). [CrossRef]

28. M. L. Yang, Y. Q. Wu, X. Q. Sheng, and K. F. Ren, “Computation of scattering of a zero-order Bessel beam by large nonspherical homogeneous particles with the multilevel fast multipole algorithm,” Submitted to J. Quant. Spectrosc. Radiat. Transf..

29. J. Durnin, J. Miceli Jr, and J. Eberly, “Diffraction-free beams,” Phys. Rev. Lett. 58(15), 1499 (1987). [CrossRef] [PubMed]

30. V. Garces-Chavez, D. McGloin, H. Melville, W. Sibbett, and K. Dholakia, “Simultaneous micromanipulation in multiple planes using a self-reconstructing light beam,” Nature 419(6903), 145–147 (2002). [CrossRef] [PubMed]

31. D. Grier, “A revolution in optical manipulation,” Nature 424(6950), 810–816 (2003). [CrossRef] [PubMed]

32. D. McGloin and K. Dholakia, “Bessel beams: Diffraction in a new light,” Contemporary Physics 46(1), 15–28 (2005). [CrossRef]

33. A. N. Rubinov, A. A. Afanas’Ev, I. E. Ermolaev, Y. A. Kurochkin, and S. Y. Mikhnevich, “Localization of spherical particles under the action of gradient forces in the field of a zero-order Bessel beam. Rayleigh-gans approximation,” Journal of Applied Spectroscopy 70(4), 565–572 (2003). [CrossRef]

34. L. A. Ambrosio and H. E. Hernandez-Figueroa, “Gradient forces on double-negative particles in optical tweezers using Bessel beams in the ray optics regime,” Opt. Express 18(23), 24287–92 (2010). [CrossRef] [PubMed]

35. L. A. Ambrosio and H. E. Hernandez-Figueroa, “Inversion of gradient forces for high refractive index particles in optical trapping,” Opt. Express 18(6), 5802–5808 (2010). [CrossRef] [PubMed]

36. G. Milne, K. Dholakia, D. Mcgloin, K. Volkesepulveda, and P. Zemanek, “Transverse particle dynamics in a Bessel beam,” Opt. Express 15(21), 13972–87 (2007). [CrossRef] [PubMed]

37. J. Chen, J. Ng, P. Wang, and Z. F. Lin, “Analytical partial wave expansion of vector Bessel beam and its application to optical binding: erratum,” Opt. Lett. 35(10), 1674–1676 (2010). [CrossRef] [PubMed]

38. R. X. Li, K. F. Ren, X. E. Han, Z. S. Wu, L. X. Guo, and S. X. Gong, “Analysis of radiation pressure force exerted on a biological cell induced by high-order Bessel beams using debye series,” J. Quant. Spectrosc. Radiat. Transf. 126(5), 69–77 (2013). [CrossRef]

39. F. G. Mitri, “Electromagnetic wave scattering of a high-order Bessel vortex beam by a dielectric sphere,” IEEE Trans. Antennas Propag. 59(11), 4375–4379 (2011). [CrossRef]

40. A. R. Edmonds, Angular Momentum in Quantum Mechanics (Princeton University, 1957). [CrossRef]

41. G. Gouesbet, J. J. Wang, and Y. P. Han, “Transformations of spherical beam shape coefficients in generalized Lorenz-Mie theories through rotations of coordinate systemsIII. special Euler angles,” Opt. Commun. 283(17), 3235–3243 (2010). [CrossRef]

42. J. R. Mautz and R. F. Harrington, “Electromagnetic scattering from a homogeneous material body of revolution,” Arch. Elektr. Uebertrag. 33, 71–80 (1979).

43. S. M. Rao, D. R. Wilton, and A. W. Glisson, “Electromagnetic scattering by surfaces of arbitrary shape,” IEEE Trans. Antennas Propag. 30(3), 409–418 (1982). [CrossRef]

44. H. Gan and W. C. Chew, “A discrete BCG-FFT algorithm for solving 3-D inhomogeneous scatterer problems,” Journal of Electromagnetic Waves and Applications 9(10), 1339–1357 (1995).

45. E. Bleszynski, M. Bleszynski, and T. Jaroszewicz, “AIM: Adaptive integral method for solving large-scale electromagnetic scattering and radiation problems,” Radio Science 31(5), 1225–1251 (1996). [CrossRef]

46. R. Coifman, V. Rokhlin, and S. Wandzura, “The fast multipole method for the wave equation: A pedestrian prescription,” IEEE Antennas Propagat. Mag. 35(3), 7–12 (1993). [CrossRef]

47. J. P. Barton, D. R. Alexander, and S. A. Schaub, “Theoretical determination of net radiation force and torque for a spherical particle illuminated by a focused laser beam,” J. Appl. Phys. 66(10), 4594–4602 (1989). [CrossRef]

48. J. Q. Lu, P. Yang, and X. H. Hu, “Simulations of light scattering from a biconcave red blood cell using the finite-difference time-domain method,” J. Biomed. Opt. 10(2), 024022 (2005). [CrossRef] [PubMed]

Computation of radiation pressure force exerted on arbitrary shaped homogeneous particles by high-order Bessel vortex beams using MLFMA

Abstract

1. Introduction

2. Description of the high-order Bessel vortex beam

3. Combined tangential formulation with MLFMA

4. Computation of radiation pressure force

5. Numerical results and discussions

6. Conclusion

Funding

Acknowledgments

References and links

Cited By

Figures (7)

Equations (20)

Optics Express