Analytical calculation of optical forces on spherical particles in optical tweezers: tutorial

Antonio Alvaro Ranha Neves; Carlos Lenz Cesar

doi:10.1364/JOSAB.36.001525

1. OPTICAL TRAPPING

Optical trapping is a modern tool that employs a structured laser beam to manipulate and investigate cellular-scale physio-chemical phenomena. A few of the advantages of the trapping technique are that it exerts no physical contact and exhibits a very low stiffness, which makes it an ideal tool for soft-matter research. Moreover, the trapping beam can be split or multiplexed into several traps. During the last couple of decades, optical trapping has been greatly perfected due to instrumental development and calibration techniques. Among the different optical trapping techniques, the single beam gradient optical trap is the prevailing one. This trap is also referred to as optical tweezers, a term coined by the technique’s pioneer, Arthur Ashkin [1,2], and acknowledged with a shared Noble Prize in Physics in 2018. Ashkin’s book [3], a collection of early articles with commentaries, illustrates the early years of optical tweezers research along with some of the related modeling briefly mentioned in this tutorial.

The physical principle behind optical trapping lies in the momentum of the electromagnetic field of the trapping beam as stated by Maxwell. When this beam is incident upon an object, it is scattered elastically upon their encounter, resulting in a net change to the rate of linear momentum. Because of the conservation of linear momentum, a force is exerted on the object by the light beam. To obtain an intuitive idea of the force that can be exerted by light, the light’s momentum is related to the de Broglie relation from modern physics, which indicates a minuscule momentum per photon resulting in a minuscule force comparable to the weight of cellular-sized objects. For example, a 5 μm water droplet weight imparts just as much force as the force imparted on a mirror by a 1 mW laser pointer.

To achieve a significant change in the linear momentum, a highly focused beam is required. In standard optical tweezers, this trapping beam is a single highly focused Gaussian profile laser beam. The trapped object is generally a dielectric microsphere, for which forces can be calibrated and employed as a force transducer to determine optical forces in other non-spherical samples within the same medium.

A. Optical Force Models

Early theoretical modeling for optical forces was conducted as early as 1909 by Debye, showing that the force on a sphere can be described in terms of extinction and scattering cross sections [4,5]. Debye considered a linearly polarized plane wave incident on a homogeneous spherical particle.

To theoretically determine these optical forces, initial assumptions are set based on the microsphere size with respect to the beam wavelength. For a small object that is much smaller than the beam wavelength, the object behaves like an electric dipole in agreement with Rayleigh scattering. Hence, a simple description is established for optical forces acting on an electromagnetic dipole, and this is known as the Rayleigh regime [6]. However, for objects greater than the wavelength, the trapping beam is treated as rays. The optical force is thus described in terms of geometrical optics, considering multiple reflections of the rays in the trapped microsphere to determine the net change in momentum [7,8]. However, in the case of most experiments, the trapped object diameter is within the order of the laser wavelength, and the above-mentioned approximations fail. A complete classical electromagnetic (EM) wave theory is needed to model optical forces in this intermediary regime, known as the Mie regime, although the Mie theory is applicable to spheres of any size. This analytical EM approach has been avoided, mainly owing to its troublesome calculations; however, after some algebraic derivations, simple expressions for accurate optical forces, valid for any object size, can be obtained.

The complete EM approach allows us to determine the morphology-dependent resonance, or whispering gallery modes, which have been observed in optical traps and are absent in the approximate size regime descriptions [9–11]. The early studies that considered the complete electromagnetic approach are those of Mie [12] and Lorenz [13]. Historical reviews of the initial developments and contributions are provided by Logan [14] and Kerker [5]. The Lorenz–Mie theory is an analytical method designed for plane-wave scattering by a spherical particle, which is not applicable to optical tweezers that require highly focused beams. Therefore, a generalization for the scattering of generic beams by Mie-sized particles led to the development of the generalized Lorenz–Mie theory (GLMT) [15–17]. In the GLMT framework, the incident and scattered beams are expressed in terms of vector spherical wave functions, whose expansion coefficients describe the applied beam configuration.

Once the vector scattered fields are known, we can therefore proceed with the Maxwell stress tensor (MST) approach to our problem and determine the net force or torque the beam exerts on the microsphere. This approach was first performed by Kim and Lee [18], where they presented a vector theory for optical forces acting on a sphere by a Gaussian beam described from a complex-source-point method. The treatment relied on expanding the fields in vector spherical harmonics, much like GLMT. Subsequently and in a similar manner, Marston and Crichton [19] presented an ab initio calculation of the radiation torque on an isotropic sphere illuminated by a circularly polarized plane electromagnetic wave. The same approach was soon thereafter extended to optical forces and torques by Barton et al. [20], considering a focused beam model from Davis [21,22]. The Davis beam model depicts a paraxial Gaussian beam with higher-order corrections.

From the initial development of GLMT, most applications required knowledge of the scattered field or beam-shape coefficients, which are none other than the multipole expansion coefficients. Several strategies have been developed for the determination of the beam-shape coefficients (BSCs) in the most efficient and/or accurate way.

Apart from the several ways to determine the BSCs, at the heart of this development, the quest for the appropriate incident beam model was also highlighted. The suggested and employed beam models were not a solution to Maxwell’s equation, and most of them were approximations that broke down for a highly focused beam.

Many researchers contributed to the understanding of optical forces by ways of theory and computational modeling [17,23–30]. In this tutorial, we wish to present what we believe is a didactic approach toward an analytical solution, considering a physical and realistic trapping beam along with the exact determination of its beam-shape coefficients. With the MST approach, an exact model for optical forces in a typical optical tweezers is obtained [11]. Therefore, we merge the main contributions from some of the studies highlighted in this work into what we believe is a robust and rigorous treatment for optical forces in optical tweezers. Similar approaches using BSC calculation via plane wave spectrum or for Gaussian beams independent of the radial coordinate have been published [31–34], but were unaware of the proposed method in Ref. [11].

A standard and modern notation for the equations was adopted throughout this work. The tutorial is divided into four parts. The first part introduces the EM fields along with the multipole coefficients. The second part explores how the BSCs are determined, and this is illustrated for the simple case of a plane wave. The third part presents the highly focused beam description along with the BSCs. Finally, the optical force through the Maxwell stress tensor is presented in terms of the BSCs. We hope to have addressed all the important points and developments described above.

2. MODELING

Hence, how do we start modeling the optical forces? Primarily, we rely on the knowledge of the exact EM fields, which travel as a wave with linear momentum. Because the trapped object is a microsphere, it seems reasonable to express the EM fields in terms of spherical harmonic functions, similarly to the standard textbook approach [35]. The encountered problem is essentially the elastic scattering of monochromatic coherent light by the homogeneous and isotropic sphere of radius ( $a$ ).

A. Description of the EM Fields

Let us describe all the required electromagnetic fields, incident (inc), scattered (sca), and internal (int) as time harmonics in terms of spherical vector expansion, also known as multipole expansion, as our scattering particle has spherical symmetry [36]. The time dependence $e^{- i ω t}$ is omitted for clarity. Our coordinate system has its origin at the center of our spherical scatterer, which scatters the incoming light. A similar approach with other notations is present in several textbooks [5,37].

E_{inc} (r) = E_{0} \sum_{l = 1}^{\infty} \sum_{m = - l}^{l} (G_{l m}^{TE} M_{l m}^{(1)} (r) + G_{l m}^{TM} N_{l m}^{(1)} (r)),

Z H_{inc} (r) = E_{0} \sum_{l = 1}^{\infty} \sum_{m = - l}^{l} (G_{l m}^{TM} M_{l m}^{(1)} (r) - G_{l m}^{TE} N_{l m}^{(1)} (r)),

E_{sca} (r) = E_{0} \sum_{l = 1}^{\infty} \sum_{m = - l}^{l} (b_{l m} M_{l m}^{(3)} (r) + a_{l m} N_{l m}^{(3)} (r)),

Z H_{sca} (r) = E_{0} \sum_{l = 1}^{\infty} \sum_{m = - l}^{l} (a_{l m} M_{l m}^{(3)} (r) - b_{l m} N_{l m}^{(3)} (r)),

E_{int} (r) = E_{0} \sum_{l = 1}^{\infty} \sum_{m = - l}^{l} (d_{l m} M_{l m}^{(1)} (r) + c_{l m} N_{l m}^{(1)} (r)),

Z H_{int} (r) = E_{0} \sum_{l = 1}^{\infty} \sum_{m = - l}^{l} (c_{l m} M_{l m}^{(1)} (r) - d_{l m} N_{l m}^{(1)} (r)),

where

l

and

m

are the radial and azimuthal mode index, respectively. The electric field dimension constant

E_{0}

is used to make the expansion coefficients non-dimensional. The mode amplitudes

a_{l m}, b_{l m}, c_{l m}, d_{l m}, G_{l m}^{TM}

, and

G_{l m}^{TE}

are also called expansion coefficients or multipole coefficients; however, the last pair is special, and we will refer to them as beam-shape coefficients (BSCs). The vector spherical harmonics are

k N_{l m} (r) = i \nabla \times M_{l m} (r)

,

M_{l m} (r) = z_{l} (k r) X_{l m} (\hat{r})

, and

X (\hat{r}) = L Y_{l m} (\hat{r}) / \sqrt{l (l + 1)}

[38], where

z_{l} (k r)

is the appropriate spherical Bessel function,

Y_{l m} (\hat{r})

are the scalar spherical harmonics,

Z = \sqrt{μ / ε}

,

k = ω \sqrt{μ ε}

, and

L = - i r \times \nabla

is the angular momentum operator in direct space. The above vector functions are solutions to the vector wave equation, first proposed by Hansen [39], and any field can be expanded in these series given the appropriate coefficients.

The spatial function $z_{l}^{(j)} (k r)$ is a basic solution to the spherical Bessel equation, which for $j = 1, 2, 3, 4$ corresponds to the ordinary spherical Bessel, spherical Neumann, and spherical Hankel functions of the first and second kind, respectively. The choice of the function depends on the particular spatial bounds, either for the fields inside or outside the scatterer. Because the fields must be regular at the origin, only $z_{l}^{(1)} (k r) = j_{l} (k_{1} r)$ can appear at the interior radial function and also for the incident fields that exist in the origin, independent of the scatterer. In the region outside the scatterer, we can write in terms of spherical Hankel functions of the outgoing type based on the physical ground of the scattered wave $z_{l}^{(3)} (k r) = h_{l}^{(1)} (k_{1} r)$ .

B. Scattering Coefficients

To solve the electromagnetic boundary problem on the scatterer surface, the continuity of the tangential components of the electric and magnetic fields is required. The following equations need to be satisfied:

[E_{inc} (r) + E_{sca} (r) - E_{int} (r)] \times \hat{r} |_{r = a} = 0,

[H_{inc} (r) + H_{sca} (r) - H_{int} (r)] \times \hat{r} |_{r = a} = 0 .

To solve each of these tangential components, let us introduce two new abbreviations,

x = k_{2} a = \frac{2 π n_{2} a}{λ}, M x = \frac{n_{1}}{n_{2}} k_{2} a = k_{1} a,

where

x

is the size parameter,

M

is the relative refractive index,

n_{1}

and

n_{2}

are the refractive index of the particle and the surrounding medium, respectively. Note that the refractive index of the particle

n_{1}

is allowed to be complex to treat the effects of absorption such as the ones in metallic particles. For the electric TE modes, the solution is

\sum_{l = 1}^{\infty} \sum_{m = - l}^{l} (G_{l m}^{TE} j_{l} (x) + b_{l m} h_{l}^{(1)} (x) - d_{l m} j_{l} (M x)) X_{l m} (\hat{r}) \times \hat{r} = 0 .

This is valid for all angles; therefore,

G_{l m}^{TE} j_{l} (x) + b_{l m} h_{l}^{(1)} (x) - d_{l m} j_{l} (M x) = 0 .

Analogously, for the magnetic TE modes, the condition valid for all angles is as follows:

G_{l m}^{TM} j_{l} (x) + a_{l m} h_{l}^{(1)} (x) - \frac{Z_{1}}{Z_{2}} c_{l m} j_{l} (M x) = 0 .

The terms proportional to

N_{l m} (r)

simplify to the following:

N_{l m} (r) \times \hat{r} = \frac{i}{k} (\nabla \times (z_{l} (k r) X_{l m} (\hat{r}))) \times \hat{r} = \frac{i}{k r} \frac{\partial (k r z_{l} (k r))}{\partial k r} .

For the TM modes, we obtain

G_{l m}^{TM} \frac{\partial (k_{1} r j_{l} (k_{1} r))}{\partial k_{1} r} + a_{l m} \frac{\partial (k_{1} r h_{l}^{(1)})}{\partial k_{1} r} - c_{l m} \frac{1}{M} \frac{\partial (k_{2} r j_{l} (k_{2} r))}{\partial k_{2} r} = 0 .

Moreover, for the magnetic equation, we have

G_{l m}^{TE} \frac{\partial (k_{1} r j_{l} (k_{1} r))}{\partial k_{1} r} + b_{l m} \frac{\partial (k_{1} r h_{l}^{(1)})}{\partial k_{1} r} - \frac{Z_{1}}{Z_{2}} d_{l m} \frac{1}{M} \frac{\partial (k_{2} r j_{l} (k_{2} r))}{\partial k_{2} r} = 0 .

To further simplify the expressions above, we introduce the Riccati–Bessel functions

ψ_{l} (y) = y j_{l} (y), ξ_{l} (y) = y h_{l}^{(1)} (y) .

These lead to

G_{l m}^{TM} ψ_{l}^{'} (x) + a_{l m} ξ_{l}^{'} (x) - c_{l m} \frac{1}{M} ψ_{l}^{'} (M x) = 0

and

G_{l m}^{TE} ψ_{l}^{'} (x) + b_{l m} ξ_{l}^{'} (x) - \frac{Z_{1}}{Z_{2}} \frac{1}{M} d_{l m} ψ_{l}^{'} (M x) = 0 .

Considering the relations between the EM impedance and the waves at optical frequencies

μ_{1} \approx μ_{2}

,

\frac{Z_{1}}{Z_{2}} \frac{1}{M} = \sqrt{\frac{μ_{1} ϵ_{2}}{ϵ_{1} μ_{2}}} \frac{n_{1}}{n_{2}} = \sqrt{\frac{μ_{1} ϵ_{2}}{ϵ_{1} μ_{2}}} \sqrt{\frac{μ_{1} ϵ_{1}}{μ_{2} ϵ_{2}}} \approx 1 .

The resulting four equations relating all multipole coefficients are

G_{l m}^{TE} ψ_{l} (x) + b_{l m} ξ_{l} (x) - \frac{1}{M} d_{l m} ψ_{l} (M x) = 0,

G_{l m}^{TM} ψ_{l} (x) + a_{l m} ξ_{l} (x) - c_{l m} ψ_{l} (M x) = 0,

G_{l m}^{TM} ψ_{l}^{'} (x) + a_{l m} ξ_{l}^{'} (x) - \frac{1}{M} c_{l m} ψ_{l}^{'} (M x) = 0,

G_{l m}^{TE} ψ_{l}^{'} (x) + b_{l m} ξ_{l}^{'} (x) - d_{l m} ψ_{l}^{'} (M x) = 0 .

Because we are interested in the scattered fields, we can eliminate

c_{l m}

and

d_{l m}

. First, we isolate

c_{l m}

in Eq. (21) and substitute in Eq. (22) as follows:

a_{l m} = - G_{l m}^{TM} [\frac{M ψ_{l}^{'} (x) ψ_{l} (M x) - ψ_{l} (x) ψ_{l}^{'} (M x)}{M ξ_{l}^{'} (x) ψ_{l} (M x) - ξ_{l} (x) ψ_{l}^{'} (M x)}] = - G_{l m}^{TM} a_{l} .

Second, we isolate

d_{l m}

in Eq. (23) and substitute in Eq. (20) as follows:

b_{l m} = - G_{l m}^{TE} [\frac{M ψ_{l} (x) ψ_{l}^{'} (M x) - ψ_{l}^{'} (x) ψ_{l} (M x)}{M ξ_{l} (x) ψ_{l}^{'} (M x) - ξ_{l}^{'} (x) ψ_{l} (M x)}] = - G_{l m}^{TE} b_{l},

where the Mie scattering coefficients are thus defined as in Ref. [37] as

a_{l} = \frac{M ψ_{l}^{'} (x) ψ_{l} (M x) - ψ_{l} (x) ψ_{l}^{'} (M x)}{M ξ_{l}^{'} (x) ψ_{l} (M x) - ξ_{l} (x) ψ_{l}^{'} (M x)},

b_{l} = \frac{M ψ_{l} (x) ψ_{l}^{'} (M x) - ψ_{l}^{'} (x) ψ_{l} (M x)}{M ξ_{l} (x) ψ_{l}^{'} (M x) - ξ_{l}^{'} (x) ψ_{l} (M x)} .

A complete and efficient description of the scattered EM wave can be obtained if information about the scatterer and the incident field is available. The GLMT framework addresses the fields in all space (near and far fields), allowing for a complicated beam that satisfies Maxwell’s equations. In our case, the scatterer is a homogeneous and isotropic dielectric microsphere, and it is completely characterized by the Mie scattering coefficients (

a_{l}

and

b_{l}

) through the knowledge of the microsphere radius (

a

), the refractive index of the microsphere (

n_{1}

), the refractive index of the surrounding medium (

n_{2}

), and the wavelength of the incident wave (

λ

). For the incident beam, this information is present in the beam-shape coefficients (

G_{l m}^{TM}

and

G_{l m}^{TE}

), which have been studied extensively by Gouesbet, Gréhan, and Lock [17,40] (and references therein). The method by which we determine the beam-shape coefficients is presented in the following section.

3. BEAM-SHAPE COEFFICIENTS

The description of BSCs appears in the GLMT framework, which was initially pioneered when studying plasma velocity by laser Doppler velocimetry [17]. The results for the scattering of any beam are known once the beam-shape coefficient is determined. In fact, GLMT can be applied to analytically solve scattering from other objects besides the homogeneous sphere, including the stratified sphere, spheroids, and infinite cylinders of any size. Moreover, any type of beam can be used, the solution not being restricted to Gaussian beams. The main problem in many studies and the subject of much debate is the calculation of the BSCs, mainly due to the non-cancellation of the radial spherical Bessel functions in both sides of the expansion. Primarily, they are defined as follows.

In Eqs. (1)–(6) we described all electromagnetic fields in terms of spherical vector expansion. In the previous section, we explicitly determined the complex scattering coefficients in terms of the BSCs $G_{l m}^{TM}$ and $G_{l m}^{TE}$ . We can, therefore, express any incident EM field by determining these BSCs. Note that $M_{l m} (r)$ is perpendicular to the radial vector $\hat{r}$ ; we can use this and the normalization relation of spherical harmonics

\int Y_{l m}^{*} Y_{p q} d Ω = δ_{l p} δ_{m q},

where

Y_{l m} = \sqrt{\frac{2 l + 1}{4 π} \frac{(l - m)!}{(l + m)!}} P_{l}^{m} (\cos (θ)) e^{i m ϕ} .

To determine the BSCs

G_{l m}^{TM}

and

G_{l m}^{TE}

, we first multiply by the radial vector

r

as follows:

r \cdot E_{inc} (r) = E_{0} \sum_{l = 1}^{\infty} \sum_{m = - l}^{l} G_{l m}^{TM} r \cdot N_{l m}^{(1)} (r) .

The last vector product in the equation above is determined by expanding in terms of spherical vector harmonics and later in terms of spherical harmonics:

r \cdot N_{l m}^{(1)} (r) = - \frac{j_{l} (k r)}{k} \sqrt{l (l + 1)} Y_{l m} .

Therefore, Eq. (30) becomes

r \cdot E_{inc} (r) = - \frac{E_{0}}{k} \sum_{l = 1}^{\infty} \sum_{m = - l}^{l} G_{l m}^{TM} j_{l} (k r) \sqrt{l (l + 1)} Y_{l m} .

To eliminate the double sum, we will multiply both sides by

Y_{p q}^{*}

and integrate over the solid angle as follows:

\int Y_{p q}^{*} r \cdot E_{inc} (r) d Ω = - \frac{E_{0}}{k} \sum_{l = 1}^{\infty} \sum_{m = - l}^{l} G_{l m}^{TM} j_{l} (k r) \sqrt{l (l + 1)} \int Y_{p q}^{*} Y_{l m} d Ω .

Using the normalization property of Eq. (28), we obtain

\int Y_{p q}^{*} r \cdot E_{inc} (r) d Ω = - \frac{E_{0}}{k} G_{p q}^{TM} j_{p} (k r) \sqrt{p (p + 1)}

or

j_{p} (k r) G_{p q}^{TM} = - k r \frac{1}{E_{0} \sqrt{p (p + 1)}} \int Y_{p q}^{*} \hat{r} \cdot E_{inc} (r) d Ω .

Analogous for the magnetic field, the corresponding BSC is given by

j_{p} (k r) G_{p q}^{TE} = k r \frac{Z}{E_{0} \sqrt{p (p + 1)}} \int Y_{p q}^{*} \hat{r} \cdot H_{inc} (r) d Ω .

Equations (35) and (36) above define the BSCs for the incident trapping beam from the GLMT framework. The right-hand side of these equations should contain a spherical Bessel function from the integration to cancel with the one on the left-hand side. The major shortcoming of the other methods is that this radial dependent term has never been explicitly determined to produce BSC independent of the radial coordinates.

The radial component of the electric and magnetic fields of the incident beam needs to be determined. Once the exact EM is known or described in an analytical expression, we are able to recover an exact BSC. However, if the incident EM fields is only approximated, we no longer retain an exact BSC. This is the case for the Davis beam, which results in a BSC that depends on the radial coordinate ( $r$ ). By definition, the expansion coefficient should not depend on the radial coordinate, and hence Barton suggested to determine the BSCs at the particle surface ( $r = a$ ) [20]. However, what would happen if $k a$ is a root of the spherical Bessel function ( $j_{l} (k a)$ )? Note that there is an infinite sum in the order $l$ and infinite roots for each order.

A. BSC for a Plane Wave

The simplest incident field is a plane wave. Here, we determine the BSCs of the plane wave by direct integration. Consider an incident coherent monochromatic plane wave propagating in non-absorbing homogeneous media in the $z$ direction, which is linearly polarized in the $x$ axis as

E_{inc} = E_{0} e^{i k z} \hat{x} .

In spherical coordinates, the radial component is

\hat{r} \cdot E_{inc} (r) = E_{0} e^{i k r \cos θ} \sin θ \cos ϕ .

We can now calculate the plane wave BSCs from the definition of Eq. (35) as

G_{l m}^{TM} = - \frac{k r}{j_{l} (k r)} \frac{1}{\sqrt{l (l + 1)}} \int Y_{l m}^{*} e^{i k r \cos θ} \sin θ \cos ϕ d Ω = - \frac{k r}{j_{l} (k r)} \frac{1}{\sqrt{l (l + 1)}} \sqrt{\frac{2 l + 1}{4 π} \frac{(l - m)!}{(l + m)!}} \int_{θ = 0}^{π} P_{l}^{m} (\cos θ) e^{i k r \cos θ} \sin^{2} θ d θ \int_{ϕ = 0}^{2 π} e^{- i m ϕ} \cos ϕ d ϕ .

The solution for the integral in

ϕ

is

\int_{ϕ = 0}^{2 π} e^{- i m ϕ} \cos ϕ d ϕ = π (δ_{m, + 1} + δ_{m, - 1}) .

Therefore, for plane waves only BSCs with

m = \pm 1

remain. Let us determine the BSC for

m = + 1

by

G_{l, + 1}^{TM} = - \frac{k r}{j_{l} (k r)} \frac{π}{l (l + 1)} \sqrt{\frac{2 l + 1}{4 π}} \int_{θ = 0}^{π} P_{l}^{1} (\cos θ) e^{i k r \cos θ} \sin^{2} θ d θ .

The

θ

integral can be solved using the integral representation of the spherical Bessel function by

\int_{θ = 0}^{π} \sin^{2} θ e^{i k r \cos θ} P_{l}^{1} (\cos θ) d θ = 2 i^{l + 1} l (l + 1) \frac{j_{l} (k r)}{k r} .

Equation (41) thus becomes

G_{l, + 1}^{TM} = - \sqrt{π (2 l + 1)} i^{l + 1} .

For

m = - 1

, we have

G_{l, - 1}^{TM} = \sqrt{π (2 l + 1)} i^{l + 1} .

Similarly, for the other BSC,

G_{l m}^{TE} = \frac{k r}{j_{l} (k r)} \frac{1}{\sqrt{l (l + 1)}} \int Y_{p q}^{*} e^{i k r \cos θ} \sin θ \sin ϕ d Ω = \frac{k r}{j_{l} (k r)} \frac{1}{\sqrt{l (l + 1)}} \sqrt{\frac{2 l + 1}{4 π} \frac{(l - m)!}{(l + m)!}} \int_{θ = 0}^{π} P_{l}^{m} (\cos θ) e^{i k r \cos θ} \sin^{2} θ d θ \int_{ϕ = 0}^{2 π} e^{- i m ϕ} \sin ϕ d ϕ .

The

ϕ

integral now becomes

\int_{ϕ = 0}^{2 π} e^{- i m ϕ} \sin ϕ d ϕ = - i π (δ_{m, + 1} - δ_{m, - 1}) .

For

m = \pm 1

,

G_{l, \pm 1}^{TE} = - i \sqrt{π (2 l + 1)} i^{l + 1} .

We can verify that these are the true BSCs by inserting them into the radial component of the incident electric field [Eq. (1)] as follows:

\hat{r} \cdot E_{inc} (r) = E_{0} \sum_{l = 1}^{\infty} \sqrt{π (2 l + 1)} i^{l + 1} \hat{r} \cdot (- N_{l, + 1} (r) + N_{l, - 1} (r)) = E_{0} \sum_{l = 1}^{\infty} \frac{j_{l} (k r)}{k r} \sqrt{l (l + 1)} \sqrt{π (2 l + 1)} i^{l + 1} (- Y_{l, + 1} + Y_{l, - 1}) = - i \frac{E_{0}}{k r} \cos ϕ \sin θ \frac{d}{d (\cos θ)} \sum_{l = 1}^{\infty} j_{l} (k r) i^{l} (2 l + 1) P_{l} (\cos θ) = E_{0} e^{i k r \cos θ} \sin θ \cos ϕ .

The last summation in Eq. (48) is a well-known relation, known as the plane wave expansion or Rayleigh equation, which describes the plane wave as a sum of spherical waves. As a result, we recover the expected electric field of Eq. (38). Having explained how to obtain the plane wave BSC, let us turn to the beam of interest for optical tweezers, namely, the highly focused Gaussian beam.

4. HIGHLY FOCUSED GAUSSIAN BEAM

Another great challenge for a proper beam description is using an expression for the EM field that exactly satisfies Maxwell’s equation. The beam representation must represent a solution of the vector Helmholtz equation, and the paraxial Gaussian beam is not a solution. The use of these pseudo-fields [41] or non-Maxwellian beams would result in an approximate beam description leading to a pair of BSCs [Eqs. (35) and (36)] that depend on the radial coordinate [42].

Early studies modeled a focused Gaussian beam as a paraxial beam perturbed up to the fifth order through a series expansion in terms of a $s$ parameter [ $s = 1 / (k w_{0})$ , where $w_{0}$ is the beam waist at the focus] [22]. This Davis–Barton fifth-order approximation still does not solve Maxwell’s equation exactly and fails to converge to reproduced experimental results under highly focused conditions ( $s \to 1$ ) as is the case in optical trapping [43,44]. The main limitation of beam models is that they do not exactly satisfy Maxwell’s equation and do not represent a good approximation of a real beam when the spot size is too small. Most paraxial beams used in focusing are also non-physical because they do not exhibit a change in spot size due to polarization alone [45]. Moreover, the diffraction effects due to the truncated Gaussian beam by the finite objective lens and the polarization mixing by the high-numerical-aperture (NA) objective are not correctly represented near the focus of a high-NA objective [24,46].

The suggested most efficient way to overcome this issue of an accurate beam model for optical force prediction was to treat the high-NA focusing as a vector diffraction problem [43].

At high-NA focusing, the paraxial scalar theory does not take into account the vector nature of the beam, it does not justify small-angle approximations, and it does not consider the apodization effects of optical systems [47,48]. This has been thoroughly tested, proven, and compared in numerous experiments, such as the 4Pi confocal fluorescence microscope [49], which led to the development of super-resolved fluorescence microscopy. A tutorial for diffraction calculation in optical systems has recently been published [50], discussing diffraction integrals and high-NA focusing.

The solution for the highly focused beam that represents the vector property of the incident beam is the vectorial Debye–Wolf diffraction integral method proposed by Richards and Wolf [51]. This model describes the paraxial Gaussian beam being focused by a highly convergent aplanatic lens as a sum of plane waves arriving from all the angles spanned by the lens’s numerical aperture as an angular spectrum of plane waves [52]. Recently, it was highlighted that this approach to vector diffraction dates to a much earlier publication in 1920 by Ignatovsky [53,54].

A. Angular Spectrum Representation of Focused Gaussian Beams

The angular spectrum representation depicts a mathematical instrument to calculate the electromagnetic fields near a highly focused beam (by a lens of focal length $f$ ). It comprises the decomposition of the EM field into an infinite series expansion in terms of Fourier components. With respect to optical tweezers, we are concerned with the propagating part of the angular spectrum, which can be solved with the method of stationary phase [55].

The resulting electric field in the vicinity of the focus ( $E_{inc} (ρ, ϕ, z)$ ), as lens focus is chosen as the origin of the coordinate system, with respect to the field on a reference sphere ( $E_{\infty} (θ, ϕ)$ ) is given in cylindrical coordinates as follows:

E_{inc} (ρ, ϕ, z) = \frac{i k f e^{i k f}}{2 π} \int_{0}^{θ_{NA}} \int_{0}^{2 π} E_{\infty} (θ, ϕ) e^{i k z \cos θ} e^{i k ρ \sin θ \cos (ϕ - ϕ)} \sin θ d θ d ϕ,

where

θ_{NA} = \sin^{- 1} (NA / n_{a})

is the maximum acceptance angle given by the NA of the objective lens. The field on the reference sphere (

E_{\infty}

) is itself determined from the paraxial Gaussian

{TEM}_{00}

beam as

E_{\infty} = \sqrt{\frac{n_{b}}{n_{a}} \cos θ} ([E_{par} \cdot \hat{ϕ}] \hat{ϕ} + [E_{par} \cdot \hat{ρ}] \hat{θ}),

where

n_{b}

and

n_{a}

refer to the refractive index of the medium before and after the focusing lens, respectively. The electric field in the focal region, Eq. (49), is described in cylindrical coordinates and completely satisfies Maxwell’s equation. The angular spectrum of a Gaussian beam requires a corresponding Gaussian spectrum that includes evanescent components. These non-homogeneous components decay exponentially with distance from the reference sphere. In optical tweezers, the focal region would be far-field from the reference sphere. Therefore, this method allows for a more rigorous description of the focused field in comparison to the previous ones [52].

To determine the BSCs, we need to convert to spherical coordinates. We would also like to displace this beam from its origin (the focal point) to an arbitrary point in space ( $ρ_{0}$ , $ϕ_{0}$ , $z_{0}$ ) (Fig. 1). This is achieved with the following substitutions:

ρ \cos ϕ = r \sin θ \cos ϕ - ρ_{0} \cos ϕ_{0},

ρ \sin ϕ = r \sin θ \sin ϕ - ρ_{0} \sin ϕ_{0},

z = r \cos θ - z_{0} .

The purpose of shifting the EM field with respect to the sphere origin provides a later advantage in describing the scattering, or the corresponding optical forces, when the focused beam is displaced with respect to the sphere center. In general particle scattering theory, such as the transition matrix (T-matrix), is formulated in the coordinate system located at the particle center [56,57]. The beam and its BSCs are described with respect to another origin with the use of the translational-addition theorem [58–60], which adds another double sum operation that needs to be truncated appropriately. In our case, this is avoided by first generalizing the focused beam description (and its corresponding BSCs) with respect to any arbitrary location.

Fig. 1. Illustration of the variables from the incoming beam (red rays), mapped onto a reference sphere and later focused (converging red rays). The focus is shifted from the origin of coordinates located at the center of the spherical scatterer.

Download Full Size | PDF

We now modify the spherical angular coordinates of the focusing system output pupil to avoid confusion with the BSCs’ integration variables from ( $θ$ , $ϕ$ ) to ( $α$ , $β$ ). The focused field becomes

E_{inc} (ρ, ϕ, z) = \frac{i k f e^{i k f}}{2 π} \int_{0}^{α_{NA}} d α \sin α e^{i k r \cos α \cos θ} e^{- i k z_{0} \cos α} \int_{0}^{2 π} d β E_{\infty} (α, β) e^{i k r \sin α \sin θ \cos (β - ϕ)} e^{- i k ρ_{0} \sin α \cos (β - ϕ_{0})} .

The amplitude of the incoming field, before focusing, has a Gaussian distribution and general polarization as follows:

E_{par} = E_{0} e^{- f^{2} \sin^{2} α / ω^{2}} (p_{x} \hat{x} + p_{y} \hat{y}),

where

ω

is the beam waist at the objective aperture, and (

p_{x}

,

p_{y}

) are the complex polarization parameters describing the state of polarization of the beam before entering an aplanatic objective lens. For example, linear polarization in the

x

axis corresponds to (1,0). Consequently, the magnetic field is given as follows:

Z H_{par} = \hat{z} \times E_{par} = E_{0} e^{- f^{2} \sin^{2} α / ω^{2}} (p_{x} \hat{y} - p_{y} \hat{x}) .

That is, to switch to the BSC, which depends on the magnetic field, we have to change

p_{x} \to - p_{y}

and

p_{y} \to p_{x}

. Hence, we can determine Eq. (50) with the above incident field by changing to the new

α

and

β

angle coordinates. For the BSC, we only require the radial component of the incident field, which implies the radial component of the field at the reference sphere:

\hat{r} \cdot E_{\infty} = E_{x, \infty} \sin θ \cos ϕ + E_{y, \infty} \sin θ \sin ϕ + E_{z, \infty} \cos θ .

Thus, the field at the reference sphere and the radial component at the focus are

\hat{r} \cdot E_{\infty} (α, β) = \sqrt{\frac{n_{b}}{n_{a}} \cos α} E_{0} e^{- f^{2} \sin^{2} α / ω^{2}} ([- \sin β [- p_{x} \sin β + p_{y} \cos β] + \cos β \cos α [p_{x} \cos β + p_{y} \sin β]] \sin θ \cos ϕ + [\cos β [- p_{x} \sin β + p_{y} \cos β] + \sin β \cos α [p_{x} \cos β + p_{y} \sin β]] \sin θ \sin ϕ - \sin α [p_{x} \cos β + p_{y} \sin β] \cos θ),

\hat{r} \cdot E_{inc} (ρ, ϕ, z) = \frac{i k f e^{i k f} E_{0}}{2 π} \sqrt{\frac{n_{b}}{n_{a}}} \int_{0}^{α_{NA}} d α \sin α \sqrt{\cos α} e^{- f^{2} \sin^{2} α / ω^{2}} e^{- i k z_{0} \cos α} e^{i k r \cos α \cos θ} \int_{0}^{2 π} d β e^{i k r \sin α \sin θ \cos (β - ϕ)} e^{- i k ρ_{0} \sin α \cos (β - ϕ_{0})} {[- \sin β (- p_{x} \sin β + p_{y} \cos β) + \cos β \cos α (p_{x} \cos β + p_{y} \sin β)] \sin θ \cos ϕ + [\cos β (- p_{x} \sin β + p_{y} \cos β) + \sin β \cos α (p_{x} \cos β + p_{y} \sin β)] \sin θ \sin ϕ - \sin α (p_{x} \cos β + p_{y} \sin β) \cos θ} .

These double integrals are addressed in the following section.

B. Beam-Shape Coefficients of the Focused Gaussian Beam

We can now determine the BSCs using the definition in Eq. (35) with the radial component of the incident field described in Eq. (59) as follows:

G_{l m}^{TM} = - \frac{k r}{j_{l} (k r)} \frac{1}{E_{0} \sqrt{l (l + 1)}} \sqrt{\frac{2 l + 1}{4 π} \frac{(l - m)!}{(l + m)!}} \int_{0}^{π} d θ P_{l}^{m} (\cos θ) \sin θ \int_{0}^{2 π} d ϕ e^{- i m ϕ} \hat{r} \cdot E_{inc} (r) = - \frac{k r}{j_{l} (k r)} \frac{i k f e^{i k f}}{2 π} \sqrt{\frac{n_{b}}{n_{a}}} \sqrt{\frac{2 l + 1}{4 π l (l + 1)} \frac{(l - m)!}{(l + m)!}} \int_{0}^{α_{NA}} d α \sin α \sqrt{\cos α} e^{- f^{2} \sin^{2} α / ω^{2}} e^{- i k z_{0} \cos α} \int_{0}^{2 π} d β e^{- i k ρ_{0} \sin α \cos (β - ϕ_{0})} \int_{0}^{π} d θ P_{l}^{m} (\cos θ) \sin θ e^{i k r \cos α \cos θ} \int_{0}^{2 π} d ϕ e^{- i m ϕ} e^{i k r \sin α \sin θ \cos (β - ϕ)} {[- \sin β (- p_{x} \sin β + p_{y} \cos β) + \cos β \cos α (p_{x} \cos β + p_{y} \sin β)] \sin θ \cos ϕ + [\cos β (- p_{x} \sin β + p_{y} \cos β) + \sin β \cos α (p_{x} \cos β + p_{y} \sin β)] \sin θ \sin ϕ - \sin α [p_{x} \cos β + p_{y} \sin β] \cos θ} .

We proceed by solving three of the four integrals above. The first integral to solve is that of

ϕ

, where we use the following relation:

\int_{0}^{2 π} d ϕ e^{- i m ϕ} e^{i x \cos (β - ϕ)} (\begin{matrix} \cos ϕ \\ \sin ϕ \\ 1 \end{matrix}) = - 2 π i^{m} e^{- i m β} (\begin{matrix} (- \sin β \frac{m J_{m} (x)}{x} + i \cos β \frac{d J_{m} (x)}{d x}) \\ (\cos β \frac{m J_{m} (x)}{x} + i \sin β \frac{d J_{m} (x)}{d x}) \\ - J_{m} (x) \end{matrix}),

where we abbreviated

y = k r \sin α \sin θ

. This solves the first of four integrals, and we will further introduce some dummy variables to shorten the algebra as follows:

A = i^{m + 1} k f e^{i k f} \frac{k r}{j_{l} (k r)} \sqrt{\frac{n_{b}}{n_{a}}} \sqrt{\frac{2 l + 1}{4 π l (l + 1)} \frac{(l - m)!}{(l + m)!}}, B = [- p_{x} \sin β + p_{y} \cos β], C = [p_{x} \cos β + p_{y} \sin β] .

After rearranging the term proportional to

\cos α

, this becomes

G_{l m}^{TM} = A \int_{0}^{α_{NA}} d α \sin α \sqrt{\cos α} e^{- f^{2} \sin^{2} α / ω^{2}} e^{- i k z_{0} \cos α} \int_{0}^{2 π} d β e^{- i k ρ_{0} \sin α \cos (β - ϕ_{0})} e^{- i m β} \int_{0}^{π} d θ P_{l}^{m} (\cos θ) \sin θ e^{i k r \cos α \cos θ} (B \sin θ \frac{m J_{m} (y)}{y} + i C \cos α \sin θ \frac{d J_{m} (y)}{d y} + C \sin α \cos θ J_{m} (y)) .

Both

C

terms can be combined within the

θ

integral, considering it as a result of the derivative with respect to

α

as follows:

G_{l m}^{TM} = A \int_{0}^{α_{NA}} d α \sin α \sqrt{\cos α} e^{- f^{2} \sin^{2} α / ω^{2}} e^{- i k z_{0} \cos α} \frac{1}{k r} \int_{0}^{2 π} d β e^{- i k ρ_{0} \sin α \cos (β - ϕ_{0})} e^{- i q β} [\frac{m B}{\sin α} + i C \frac{d}{d α}] \int_{0}^{π} d θ P_{l}^{m} (\cos θ) \sin θ e^{i k r \cos α \cos θ} J_{m} (k r \sin α \sin θ) .

The result for the

θ

integral is given in Ref. [61] as follows:

\int_{0}^{π} d θ P_{l}^{m} (\cos θ) \sin θ e^{i k r \cos α \cos θ} J_{m} (k r \sin α \sin θ) = 2 i^{l - m} P_{l}^{m} (\cos α) j_{l} (k r) .

This reduces the BSC to only two integrals as follows:

G_{l m}^{TM} = 2 i^{l - m} \frac{j_{l} (k r)}{k r} A \int_{0}^{α_{NA}} d α \sin α \sqrt{\cos α} e^{- f^{2} \sin^{2} α / ω^{2}} e^{- i k z_{0} \cos α} \int_{0}^{2 π} d β e^{- i k ρ_{0} \sin α \cos (β - ϕ_{0})} e^{- i m β} [\frac{m B P_{l}^{m} (\cos α)}{\sin α} + i C \frac{d P_{l}^{m} (\cos α)}{d α}] .

Note that this last

θ

integral has actually cancelled the radial dependence from the spherical Bessel function defined in the BSC expression of Eq. (35). The analytical result for this integral is the key to avoiding most of the numerical approximation in the literature. This integral has been pointed out as a special case of a Gegenbauer integral [62], even though it is not generally well known and not easily available in integral tables [63]. This result has been useful in several applications other than light scattering, including acoustics [64,65], ionization rates for orbitals [66], the Casimir interaction [67], atom tunneling [68], field theory [69], and cosmology [70].

This analytical calculation was not derived for many years, and other strategies have been employed to solve the problem of determining BSCs, such as the use of van de Hulst’s localization principle and building a localized model of the beam. For an incident Gaussian beam on axis to a sphere, described by the Davis model, the transverse and angular (polar) coordinates are mapped to specific values within the localized approximation as $k r \to n + 1 / 2$ and $θ \to π / 2$ , respectively, thereby no longer depending on the radial variable and avoiding the integral over the polar component [71–73]. A historical overview of the different methods for determining the BSC can be found in Gouesbet’s book [17].

To proceed with the determination of BSCs, we integrate (with respect to $β$ ), as we need to eliminate our dummy indexes that are functions of $β$ . We introduce the following two angular functions related to the associated Legendre function, which are common in light scattering [5]:

π_{l}^{m} (α) = \frac{m P_{l}^{m} (\cos α)}{\sin α}, τ_{l}^{m} (α) = \frac{d P_{l}^{m} (\cos α)}{d α} .

After some algebra, we obtain

G_{l m}^{TM} = i^{l + 1} k f e^{i k f} \sqrt{\frac{n_{b}}{n_{a}}} \sqrt{\frac{2 l + 1}{π l (l + 1)} \frac{(l - m)!}{(l + m)!}} \int_{0}^{α_{NA}} d α \sin α \sqrt{\cos α} e^{- f^{2} \sin^{2} α / ω^{2}} e^{- i k z_{0} \cos α} \int_{0}^{2 π} d β e^{- i k ρ_{0} \sin α \cos (β - ϕ_{0})} e^{- i m β} [(- p_{x} π_{l}^{m} (α) + i p_{y} τ_{l}^{m} (α)) \sin β + (p_{y} π_{l}^{m} (α) + i p_{x} τ_{l}^{m} (α)) \cos β] .

The integral (with respect to

β

) is equal to

\int_{0}^{2 π} d β e^{- i x \cos (β - ϕ_{0})} e^{- i m β} (\begin{matrix} \cos β \\ \sin β \end{matrix}) = 2 π i^{- m} e^{- i m ϕ_{0}} (\begin{matrix} - \sin ϕ_{0} \frac{m J_{m} (x)}{x} + i \cos ϕ_{0} \frac{d J_{m} (x)}{d x} \\ \cos ϕ_{0} \frac{m J_{m} (x)}{x} + i \sin ϕ_{0} \frac{d J_{m} (x)}{d x} \end{matrix}),

where

x = k ρ_{0} \sin α

. Solving this last analytical integral, we obtain the BSC for the highly focused incident beam as follows:

G_{l m}^{TM} = i^{l - m + 1} k f e^{i k f} e^{- i m ϕ_{0}} \sqrt{\frac{n_{b}}{n_{a}}} \sqrt{\frac{4 π (2 l + 1)}{l (l + 1)} \frac{(l - m)!}{(l + m)!}} \int_{0}^{α_{NA}} d α \sin α \sqrt{\cos α} e^{- f^{2} \sin^{2} α / ω^{2}} e^{- i k z_{0} \cos α} {[- p_{x} π_{l}^{m} (α) + i p_{y} τ_{l}^{m} (α)] [\cos ϕ_{0} \frac{m J_{m} (x)}{x} + i \sin ϕ_{0} \frac{d J_{m} (x)}{d x}] + [p_{y} π_{l}^{m} (α) + i p_{x} τ_{l}^{m} (α)] [- \sin ϕ_{0} \frac{m J_{m} (x)}{x} + i \cos ϕ_{0} \frac{d J_{m} (x)}{d x}]} .

The analogous TE BSC is

G_{l m}^{TE} = - i^{l - m + 1} k f e^{i k f} e^{- i m ϕ_{0}} \sqrt{\frac{n_{b}}{n_{a}}} \sqrt{\frac{4 π (2 l + 1)}{l (l + 1)} \frac{(l - m)!}{(l + m)!}} \int_{0}^{α_{NA}} d α \sin α \sqrt{\cos α} e^{- f^{2} \sin^{2} α / ω^{2}} e^{- i k z_{0} \cos α} {[p_{y} π_{l}^{m} (α) + i p_{x} τ_{l}^{m} (α)] [\cos ϕ_{0} \frac{m J_{m} (x)}{x} + i \sin ϕ_{0} \frac{d J_{m} (x)}{d x}] + [p_{x} π_{l}^{m} (α) - i p_{y} τ_{l}^{m} (α)] [- \sin ϕ_{0} \frac{m J_{m} (x)}{x} + i \cos ϕ_{0} \frac{d J_{m} (x)}{d x}]} .

This numerical integration is quite fast once the integrating functions are known. For the calculation of optical forces with the beam positioned in space (

ρ_{0}

,

ϕ_{0}

,

z_{0}

) with respect to the sphere center, only the BSCs need to be recalculated, which is an advantage with respect to other methods that do not share this feature. These BSCs can be expressed in a more compact form, uncoupling the azimuthal coordinate (

ϕ_{0}

) and polarization (

p_{x}

,

p_{y}

) [74], as

G_{l m}^{TM, TE} = - i^{l - m + 1} k f e^{i k f} e^{- i m ϕ_{0}} \sqrt{\frac{n_{b}}{n_{a}}} \sqrt{\frac{4 π (2 l + 1)}{l (l + 1)} \frac{(l - m)!}{(l + m)!}} \int_{0}^{α_{NA}} d α \sin α \sqrt{\cos α} e^{- f^{2} \sin^{2} α / ω^{2}} e^{- i k z_{0} \cos α} (\begin{matrix} [\frac{m J_{m} (k ρ_{0} \sin α)}{k ρ_{0} \sin α} π_{l}^{m} (α) + J_{m}^{'} (k ρ_{0} \sin α) τ_{l}^{m} (α)], & - i [\frac{m J_{m} (k ρ_{0} \sin α)}{k ρ_{0} \sin α} τ_{l}^{m} (α) + J_{m}^{'} (k ρ_{0} \sin α) π_{l}^{m} (α)] \end{matrix}) \cdot (\begin{matrix} \cos ϕ_{0} & \sin ϕ_{0} \\ - \sin ϕ_{0} & \cos ϕ_{0} \end{matrix}) \cdot (\begin{matrix} p_{x}, p_{y} \\ p_{y}, - p_{x} \end{matrix}) .

The BSCs

G_{l m}^{TM, TE}

correspond to the transverse-magnetic (TM) and transverse-electric (TE) multipoles, which correspond to the left and right polarizations of the last vector in Eq. (72), respectively, whose last lines involve a (

1 \times 2

) (

2 \times 2

) (

2 \times 1

) matrix operation. The efficiency of this calculation can be higher by exploiting the symmetry of some special cases as highlighted in the following sections.

C. Case Study: Coincident Microsphere and Beam Axes

This case is a common condition under trapping equilibrium, resulting in the microsphere axis coincident with the beam axis. To determine these special BSCs, we make $x \to 0$ , which leads to the following conditions:

\lim_{x \to 0} \frac{d J_{m} (x)}{d x} = \frac{1}{2} (δ_{m, + 1} - δ_{m, - 1}), \lim_{x \to 0} \frac{m J_{m} (x)}{x} = \frac{1}{2} (δ_{m, + 1} + δ_{m, - 1}) .

For

m = \pm 1

we have

G_{l, \pm 1}^{TM} = \frac{i^{l} (k f e^{i k f})}{l (l + 1)} \sqrt{\frac{n_{b}}{n_{a}}} \sqrt{π (2 l + 1)} (\mp p_{x} + i p_{y}) \int_{0}^{α_{NA}} d α \sin α \sqrt{\cos α} e^{- f^{2} \sin^{2} α / ω^{2}} e^{- i k z_{0} \cos α} (π_{l}^{1} (α) + τ_{l}^{1} (α)),

while for

G^{TE}

we obtain

G_{l, \pm 1}^{TE} = - \frac{i^{l} (k f e^{i k f})}{l (l + 1)} \sqrt{\frac{n_{b}}{n_{a}}} \sqrt{π (2 p + 1)} (p_{y} \pm i p_{x}) \int_{0}^{α_{NA}} d α \sin α \sqrt{\cos α} e^{- f^{2} \sin^{2} α / ω^{2}} e^{- i k z_{0} \cos α} (π_{l}^{1} (α) + τ_{l}^{1} (α)),

which can be summed up as follows:

G_{l, \pm 1}^{TM} = (\mp p_{x} + i p_{y}) G_{l}, G_{l, \pm 1}^{TE} = - (p_{y} \pm i p_{x}) G_{l},

where

G_{l} = \frac{i^{l} (k f e^{i k f})}{l (l + 1)} \sqrt{\frac{n_{b}}{n_{a}}} \sqrt{π (2 l + 1)} \int_{0}^{α_{NA}} d α \sin α \sqrt{\cos α} e^{- f^{2} \sin^{2} α / ω^{2}} e^{- i k z_{0} \cos α} (π_{l}^{1} (α) + τ_{l}^{1} (α)) .

5. OPTICAL FORCES

The optical force is derived from the conservation law of electromagnetic wave momentum. For a trapped spherical scatterer, the total electromagnetic force may be determined by applying Maxwell’s stress tensor (MST) in the space surrounding the scatterer. A step-by-step application of the MST for an arbitrary incident wave and a dielectric sphere is presented. With the EM field described in terms of spherical harmonics, the observable effects of the force components that an EM wave exerts on the microsphere is given by the time average of the Maxwell stress tensor ( $T_{i j}$ ), which for steady-state optical conditions is given in the traditional Minkowski form, defined below as in Ref. [35]:

F_{i} = \oint T_{i j} n_{j} d A = \frac{1}{2} R \oint [ϵ E_{i} E_{j}^{*} + μ H_{i} H_{j}^{*} - \frac{1}{2} (ϵ E \cdot E^{*} + μ H \cdot H^{*}) δ_{i j}] n_{j} d A,

where the EM fields are the ones outside the scatterer (incoming and scattered fields). We integrate over a spherical surface surrounding the scatterer, which can be chosen arbitrarily; hence, we choose one whose radius

k r \to \infty

. Where

n_{j} = \hat{r}

, and because

d A = r^{2} d Ω

, only certain terms of the radial function expansion are needed. At infinity, the electromagnetic fields only contain transverse components, and therefore the first two terms of the right-hand side of Eq. (78) do not contribute, resulting in

F = - \frac{r^{2}}{4} R \int (ϵ E \cdot E^{*} + μ H \cdot H^{*}) \hat{r} d Ω .

The EM fields involve spherical vector terms such as

\nabla \times (z_{l} (k r) X_{l m} (\hat{r})) = k z_{l}^{'} (k r) (\hat{r} \times X_{l m} (\hat{r})) + z_{l} (k r) \nabla \times X_{l m} (\hat{r}) .

For

k r \to \infty

, this would result in

\nabla \times (z_{l} (k r) X_{l m} (\hat{r})) = k z_{l}^{'} (k r) (\hat{r} \times X_{l m} (\hat{r})) .

The asymptotic expansion of the spherical radial functions and their derivatives in this limit are

j_{l} (k r) = {(- i)}^{l + 1} \frac{e^{i k r}}{2 k r} + i^{l + 1} \frac{e^{- i k r}}{2 k r},

j_{l}^{'} (k r) = {(- i)}^{l} \frac{e^{i k r}}{2 k r} + i^{l} \frac{e^{- i k r}}{2 k r},

h_{l}^{(1)} (k r) = {(- i)}^{l + 1} \frac{e^{i k r}}{k r},

h_{l}^{' (1)} (k r) = {(- i)}^{l} \frac{e^{i k r}}{k r} .

The total EM field outside the sphere and the sum of the incident and scattered fields are known in terms of the spherical vector harmonic. Inserting these fields into Eq. (79) and taking into account that the incident field contribution is removed since it is independent of the sphere position, the contribution of both fields to the force amounts to zero, i.e., all terms that do not contain

a_{l}

and

b_{l}

. The vectorial spherical harmonic products can be re-written as follows:

(\hat{r} \times X_{l m}) \cdot (\hat{r} \times X_{p q}^{*}) = \hat{r} \cdot [X_{p q}^{*} \times (\hat{r} \times X_{l m})] = \hat{r} \cdot [\hat{r} X_{l m} \cdot X_{p q}^{*} - X_{l m} (\hat{r} \cdot X_{p q}^{*})] = X_{l m} \cdot X_{p q}^{*}

and

X_{p q}^{*} \cdot (\hat{r} \times X_{l m}) = - X_{p q}^{*} \cdot (X_{l m} \times \hat{r}) = - X_{l m} \cdot (\hat{r} \times X_{p q}^{*}) .

These simplify the electric and magnetic terms to two components each. In this way, we can sum the contributions from the electric and magnetic terms, resulting in a force

F = \frac{ϵ {| E_{0} |}^{2}}{4 k^{2}} R \sum_{l = 1}^{\infty} \sum_{m = - l}^{l} \sum_{p = 1}^{\infty} \sum_{q = - p}^{p} i^{p - l} {[G_{l m}^{TM} G_{p q}^{TM *} (a_{l} + a_{p}^{*} - 2 a_{l} a_{p}^{*}) + G_{l m}^{TE} G_{p q}^{TE *} (b_{l} + b_{p}^{*} - 2 b_{l} b_{p}^{*})] \int X_{l m} \cdot X_{p q}^{*} \hat{r} d Ω + [G_{l m}^{TM} G_{p q}^{TE *} (a_{l} + b_{p}^{*} - 2 a_{l} b_{p}^{*}) - G_{l m}^{TE} G_{p q}^{TM *} (b_{l} + a_{p}^{*} - 2 b_{l} a_{p}^{*})] \int X_{l m} \cdot (\hat{r} \times X_{p q}^{*}) \hat{r} d Ω} .

The force is generally split into components such that the explicit integrals over the solid angle are given in terms of

\int \sin θ X_{l m} \cdot X_{p q}^{*} e^{i ϕ} d Ω = \frac{δ_{q, m + 1} δ_{l, p + 1}}{p + 1} \sqrt{\frac{p (p + 2) (p - m + 1) (p - m)}{(2 p + 3) (2 p + 1)}} - \frac{δ_{q, m + 1} δ_{l + 1, p}}{l + 1} \sqrt{\frac{l (l + 2) (l + m + 2) (l + m + 1)}{(2 l + 3) (2 l + 1)}},

\int \sin θ X_{l m} \cdot (\hat{r} \times X_{p q}^{*}) e^{i ϕ} d Ω = δ_{q, m + 1} δ_{l, p} \frac{i \sqrt{(l - m) (l + m + 1)}}{l (l + 1)},

\int \cos θ X_{l m} \cdot X_{p q}^{*} d Ω = \frac{δ_{q, m} δ_{l, p + 1}}{p + 1} \sqrt{\frac{p (p + 2) (p + m + 1) (p - m + 1)}{(2 p + 3) (2 p + 1)}} + \frac{δ_{q, m} δ_{l + 1, p}}{l + 1} \sqrt{\frac{l (l + 2) (l + m + 1) (l - m + 1)}{(2 l + 3) (2 l + 1)}},

\int \cos θ X_{l m} \cdot (\hat{r} \times X_{p q}^{*}) d Ω = δ_{q, m} δ_{l, p} \frac{i m}{l (l + 1)} .

After substituting these integrals in Eq. (88), we obtain after some simplification the transverse and longitudinal force components as follows:

[\begin{matrix} F_{x} \\ F_{y} \end{matrix}] = \frac{- ϵ {| E_{0} |}^{2}}{{(2 k)}^{2}} [\begin{matrix} R \\ I \end{matrix}] \sum_{p = 1}^{\infty} \sum_{q = - p}^{p} \frac{i}{p + 1} {\sqrt{\frac{p (p + 2) (p + q + 2) (p + q + 1)}{(2 p + 3) (2 p + 1)}} \times [A_{p} G_{p + 1, - (q + 1)}^{TM} G_{p, - q}^{TM *} + A_{p}^{*} G_{p q}^{TM} G_{p + 1, q + 1}^{TM *} + B_{p} G_{p + 1, - (q + 1)}^{TE} G_{p, - q}^{TE *} + B_{p}^{*} G_{p q}^{TE} G_{p + 1, q + 1}^{TE *}] - \frac{\sqrt{(p - q) (p + q + 1)}}{p} [C_{p} G_{p q}^{TM} G_{p, q + 1}^{TE *} - C_{p}^{*} G_{p q}^{TE} G_{p, q + 1}^{TM *}]},

F_{z} = - \frac{ϵ {| E_{0} |}^{2}}{2 k^{2}} R \sum_{p = 1}^{\infty} \sum_{q = - p}^{p} \frac{i}{p + 1} {\sqrt{\frac{p (p + 2) (p + q + 1) (p - q + 1)}{(2 p + 1) (2 p + 3)}} \times [A_{p} G_{p + 1, q}^{TM} G_{p q}^{TM *} + B_{p} G_{p + 1, q}^{TE} G_{p q}^{TE *}] - \frac{q}{p} C_{p} G_{p q}^{TM} G_{p q}^{TE *}},

where the combination of Mie coefficients has been defined as the following:

A_{p} = a_{p + 1} + a_{p}^{*} - 2 a_{p + 1} a_{p}^{*}, B_{p} = b_{p + 1} + b_{p}^{*} - 2 b_{p + 1} b_{p}^{*}, C_{p} = a_{p} + b_{p}^{*} - 2 a_{p} b_{p}^{*} .

The force calculations outlined step-by-step above can easily be generalized for optical torques using the same stress tensor and guidelines presented in this tutorial. This methodology was used to calculate the force in an optical force spectroscopy experiment [75]. In Ref. [42], the figures show a comparison between the expansion method and the real beam. This model also explained the results from a double tweezers setup to perform an ultra-sensitive force measurement of the full optical force profile as a function of radial position and wavelength for both TE and transverse magnetic TM mode excitation [11]. Aberrations can also be considered, where the model successfully explained the ripple structure in the optical forces profile due to spherical aberrations due to refractive index mismatching [74].

6. FINAL REMARKS

The infinite summation for the optical force may be of concern. Computationally, we truncate it to a particular value ( $N_{\max}$ ) because all other terms have a negligible contribution. This maximum summation term scales with the scatterer size (size parameter, $x$ ), just like in the determination of the scattering cross section, which uses the Wiscombe criterion [76]. By using a different floating-point variable for the numerical computation, a suitable truncation number can be chosen [77] by taking the closest integer of

N_{\max} = 0.76 ϵ^{2 / 3} x^{1 / 3} - 4.1,

where

10^{- ϵ}

is the relative truncation error and

x

is the size parameter defined in Eq. (9). This is correct only for the total cross section, as has been emphasized recently [78], or for EM fields within a region of space delimited by

k_{2} r < N_{\max}

. However, for specific fields instead of cross sections, a similar rule is applicable in which the maximum number of terms should be greater than the largest

k_{2} r

in the plot region.

For most applications, when the bead is near its equilibrium position in the optical trap, the trap stiffness ( $κ$ ) may suffice, which can be directly calculated from the optical force from Hooke’s law as

κ_{q} = {- (\frac{\partial F_{q}}{\partial q})}_{r = r_{eq}},

where

q = x, y

, or

z

and is evaluated at the equilibrium position (

r_{eq}

). Many refinements have been made to render the approach described here more complete. One can take into account the spherical aberrations produced by the refractive index mismatch at the coverslip interface [79], whose consequences on optical forces are known [80–82]. The appropriate aberrated beam coefficients can be determined in a similar way [74], making the optical force more realistic when trapping close to an interface.

One of the main efficiencies in determining optical forces with this procedure is the analytically determined BSCs at an arbitrary position with respect to the microsphere, without the need of introducing computational error from translation operations. The present method allows us to obtain the resultant force on a spherical dielectric by displacing the particle with respect to beam focus, and hence a full force profile can be obtained.

The precise determination of the BSC for a highly focused beam has been validated experimentally, making it a robust strategy for BSC determination. One of these validations was the measurement of a complete force profile for microspheres of different sizes and the effects of coupling the beam polarization to different TM and TE modes separately [11]. Recently, this force profile or force profile mapping measurement and comparison have been performed with an improved setup [83–85]. The analytical BSCs were also validated independently in a single particle photo-thermal experiment [86].

Another method to analytically determine the BSCs was derived by Moreira et. al. using the Fourier transform of the angular momentum operator [87,88]. In this work, the BSCs are for the first time shown to be independent of the radial coordinate if the beam model exactly satisfies Maxwell’s wave equations. Several analytical BSCs for different beam modes, such as Bessel beams and circular and rectangular wave-guided beams were likewise derived. This new method also applies to plane waves and highly focused beams.

As a final contribution of this study, the computer codes generated as a result of all of the work presented in this paper are released and are freely available, also allowing for contributions from the community [89].

7. SUMMARY

We presented an overview of the approach to model optical forces occurring in optical tweezers. The optical forces are determined by the Maxwell stress tensor approach and expressed in terms of the beam-shape coefficients from the generalized Lorenz–Mie theory. A robust beam model is employed for the tightly focused Gaussian beam with a high-NA objective lens in terms of the diffraction integrals presented by Richards and Wolf and expressed by the angular spectrum representation. From this beam, the beam-shape coefficients are determined by analytically integrating under the solid angle.

This work allows for accurate and precise computational modeling to quantitatively compare the theoretical results with experimental outcomes. Using the analytical expressions presented in this work requires low computational complexity. We hope that this work can be of help to young researchers in the field and provide them with a more solid background in optical forces theory, stimulating further development for trapping applications.

Funding

Fundação de Amparo à Pesquisa do Estado de São Paulo (FAPESP) (2017/13543-3).

REFERENCES

1. A. Ashkin, “Acceleration and trapping of particles by radiation pressure,” Phys. Rev. Lett. 24, 156–159 (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, 288–290 (1986). [CrossRef]

3. A. Ashkin, Optical Trapping and Manipulation of Neutral Particles Using Lasers (World Scientific, 2006).

4. P. Debye, “Das Verhalten von Lichtwellen in der Nähe eines Brennpunktes oder einer Brennlinie,” Ann. Phys. 335, 755–776 (1909). [CrossRef]

5. M. Kerker and E. M. Loebl, eds., The Scattering of Light and Other Electromagnetic Radiation (Academic, 1969).

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

7. A. Ashkin, “Forces of a single-beam gradient laser trap on a dielectric sphere in the ray optics regime,” Biophys. J. 61, 569–582 (1992). [CrossRef]

8. R. Gussgard, T. Lindmo, and I. Brevik, “Calculation of the trapping force in a strongly focused laser beam,” J. Opt. Soc. Am. B 9, 1922–1930 (1992). [CrossRef]

9. A. Ashkin and J. M. Dziedzic, “Observation of resonances in the radiation pressure on dielectric spheres,” Phys. Rev. Lett. 38, 1351–1354 (1977). [CrossRef]

10. A. Fontes, K. Ajito, A. A. R. Neves, W. L. Moreira, A. A. de Thomaz, L. C. Barbosa, A. M. de Paula, and C. L. Cesar, “Raman, hyper-Raman, hyper-Rayleigh, two-photon luminescence and morphology-dependent resonance modes in a single optical tweezers system,” Phys. Rev. E 72, 012903 (2005). [CrossRef]

11. A. A. R. Neves, A. Fontes, L. de Y. Pozzo, A. A. de Thomaz, E. Chillce, E. Rodriguez, L. C. Barbosa, and C. L. Cesar, “Electromagnetic forces for an arbitrary optical trapping of a spherical dielectric,” Opt. Express 14, 13101–13106 (2006). [CrossRef]

12. G. Mie, “Beiträge zur Optik trüber Medien, speziell kolloidaler Metallösungen,” Ann. Phys. 330, 377–445 (1908). [CrossRef]

13. L. Lorenz, Lysbevægelsen i og uden for en af plane Lysbølger belyst Kugle, Kongelige Danske Videnskabernes Selskabs Skrifter (Naturvidenskabelig og Mathematisk Afdeling, 1890).

14. N. A. Logan, “Survey of some early studies of the scattering of plane waves by a sphere,” Proc. IEEE 53, 773–785 (1965). [CrossRef]

15. G. Gouesbet and G. Gréhan, “Sur la généralisation de la théorie de Lorenz-Mie,” J. Opt. (Paris) 13, 97–103 (1982). [CrossRef]

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

17. G. Gouesbet and G. Gréhan, Generalized Lorenz-Mie Theories (Springer, 2017).

18. J. S. Kim and S. S. Lee, “Scattering of laser beams and the optical potential well for a homogeneous sphere,” J. Opt. Soc. Am. 73, 303–312 (1983). [CrossRef]

19. P. L. Marston and J. H. Crichton, “Radiation torque on a sphere caused by a circularly-polarized electromagnetic wave,” Phys. Rev. A 30, 2508–2516 (1984). [CrossRef]

20. 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, 4594–4602 (1989). [CrossRef]

21. L. W. Davis, “Theory of electromagnetic beams,” Phys. Rev. A 19, 1177–1179 (1979). [CrossRef]

22. J. P. Barton and D. R. Alexander, “Fifth-order corrected electromagnetic field components for a fundamental Gaussian beam,” J. Appl. Phys. 66, 2800–2802 (1989). [CrossRef]

23. P. A. M. Neto and H. M. Nussenzveig, “Theory of optical tweezers,” Europhys. Lett. 50, 702–708 (2000). [CrossRef]

24. A. Rohrbach and E. H. K. Stelzer, “Optical trapping of dielectric particles in arbitrary fields,” J. Opt. Soc. Am. A 18, 839–853 (2001). [CrossRef]

25. J. A. Lock, “Calculation of the radiation trapping force for laser tweezers by use of generalized Lorenz-Mie theory. I. Localized model description of an on-axis tightly focused laser beam with spherical aberration,” Appl. Opt. 43, 2532–2544 (2004). [CrossRef]

26. J. A. Lock, “Calculation of the radiation trapping force for laser tweezers by use of generalized Lorenz-Mie theory. II. On-axis trapping force,” Appl. Opt. 43, 2545–2554 (2004). [CrossRef]

27. T. A. Nieminen, V. L. Y. Loke, A. B. Stilgoe, G. Knöner, A. M. Brańczyk, N. R. Heckenberg, and H. Rubinsztein-Dunlop, “Optical tweezers computational toolbox,” J. Opt. A 9, S196–S203 (2007). [CrossRef]

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

29. T. A. Nieminen, N. D. Preez-Wilkinson, A. B. Stilgoe, V. L. Loke, A. A. Bui, and H. Rubinsztein-Dunlop, “Optical tweezers: theory and modelling,” J. Quant. Spectrosc. Radiat. Transfer 146, 59–80 (2014). [CrossRef]

30. P. Jones, O. Marago, and G. Volpe, Optical Tweezers: Principles and Applications (Cambridge, 2015).

31. J. M. Taylor and G. D. Love, “Multipole expansion of Bessel and Gaussian beams for Mie scattering calculations,” J. Opt. Soc. Am. A 26, 278–282 (2009). [CrossRef]

32. A. Salandrino, S. Fardad, and D. N. Christodoulides, “Generalized Mie theory of optical forces,” J. Opt. Soc. Am. B 29, 855–866(2012). [CrossRef]

33. M. Šiler and P. Zemànek, “Optical trapping in secondary maxima of focused laser beam,” J. Quant. Spectrosc. Radiat. Transfer 162, 114–121 (2015). [CrossRef]

34. W. Lu, H. Chen, S. Liu, and Z. Lin, “Rigorous full-wave calculation of optical forces on dielectric and metallic microparticles immersed in a vector airy beam,” Opt. Express 25, 23238–23253 (2017). [CrossRef]

35. J. D. Jackson, Classical Electrodynamics, 3rd ed. (Wiley, 1999).

36. A. J. Devaney and E. Wolf, “Multipole expansions and plane wave representations of the electromagnetic field,” J. Math. Phys. 15, 234–244 (1974). [CrossRef]

37. C. Bohren and D. Huffman, Absorption and Scattering of Light by Small Particles, Wiley Science Series (Wiley, 2008).

38. E. L. Hill, “The theory of vector spherical harmonics,” Am. J. Phys. 22, 211–214 (1954). [CrossRef]

39. W. W. Hansen, “A new type of expansion in radiation problems,” Phys. Rev. 47, 139–143 (1935). [CrossRef]

40. G. Gouesbet, J. A. Lock, and G. Gréhan, “Partial-wave representations of laser beams for use in light-scattering calculations,” Appl. Opt. 34, 2133–2143 (1995). [CrossRef]

41. T. Nieminen, H. Rubinsztein-Dunlop, and N. Heckenberg, “Multipole expansion of strongly focused laser beams,” J. Quant. Spectrosc. Radiat. Transfer 79-80, 1005–1017 (2003). [CrossRef]

42. A. A. R. Neves, A. Fontes, L. A. Padilha, E. Rodriguez, C. H. de Brito Cruz, L. C. Barbosa, and C. L. Cesar, “Exact partial wave expansion of optical beams with respect to an arbitrary origin,” Opt. Lett. 31, 2477–2479 (2006). [CrossRef]

43. W. H. Wright, G. J. Sonek, and M. W. Berns, “Radiation trapping forces on microspheres with optical tweezers,” Appl. Phys. Lett. 63, 715–717 (1993). [CrossRef]

44. L. Novotny, E. J. Sánchez, and X. S. Xie, “Near-field optical imaging using metal tips illuminated by higher-order Hermite-Gaussian beams,” Ultramicroscopy 71, 21–29 (1998). [CrossRef]

45. R. Dorn, S. Quabis, and G. Leuchs, “Sharper focus for a radially polarized light beam,” Phys. Rev. Lett. 91, 233901 (2003). [CrossRef]

46. D. Ganic, X. Gan, and M. Gu, “Exact radiation trapping force calculation based on vectorial diffraction theory,” Opt. Express 12, 2670–2675 (2004). [CrossRef]

47. C. W. McCutchen, “Generalized aperture and the three-dimensional diffraction image,” J. Opt. Soc. Am. 54, 240–244 (1964). [CrossRef]

48. C. J. R. Sheppard and H. J. Matthews, “Imaging in high-aperture optical systems,” J. Opt. Soc. Am. A 4, 1354–1360 (1987). [CrossRef]

49. S. Hell and E. H. K. Stelzer, “Properties of a 4Pi confocal fluorescence microscope,” J. Opt. Soc. Am. A 9, 2159–2166 (1992). [CrossRef]

50. J. Kim, Y. Wang, and X. Zhang, “Calculation of vectorial diffraction in optical systems,” J. Opt. Soc. Am. A 35, 526–535 (2018). [CrossRef]

51. B. Richards and E. Wolf, “Electromagnetic diffraction in optical systems, II. structure of the image field in an aplanatic system,” Proc. R. Soc. London Ser. A 253, 358–379 (1959). [CrossRef]

52. L. Novotny and B. Hecht, Principles of Nano-Optics (Cambridge, 2012).

53. J. Peatross, M. Berrondo, D. Smith, and M. Ware, “Vector fields in a tight laser focus: comparison of models,” Opt. Express 25, 13990–14007 (2017). [CrossRef]

54. V. Ignatowsky, “Diffraction by a lens having arbitrary opening,” Trans. Opt. Inst. Petrograd1 (1919).

55. L. Mandel and E. Wolf, Optical Coherence and Quantum Optics (Cambridge, 1995).

56. D. W. Mackowski and M. I. Mishchenko, “Calculation of the T matrix and the scattering matrix for ensembles of spheres,” J. Opt. Soc. Am. A 13, 2266–2278 (1996). [CrossRef]

57. F. Borghese, P. Denti, R. Saija, and M. A. Iatì, “Optical trapping of nonspherical particles in the T-matrix formalism,” Opt. Express 15, 11984–11998 (2007). [CrossRef]

58. O. R. Cruzan, “Translational addition theorems for spherical vector wave functions,” Q. Appl. Math. 20, 33–40 (1962). [CrossRef]

59. A. Doicu and T. Wriedt, “Computation of the beam-shape coefficients in the generalized Lorenz-Mie theory by using the translational addition theorem for spherical vector wave functions,” Appl. Opt. 36, 2971–2978 (1997). [CrossRef]

60. O. Moine and B. Stout, “Optical force calculations in arbitrary beams by use of the vector addition theorem,” J. Opt. Soc. Am. B 22, 1620–1631 (2005). [CrossRef]

61. A. A. R. Neves, L. A. Padilha, A. Fontes, E. Rodriguez, C. H. B. Cruz, L. C. Barbosa, and C. L. Cesar, “Analytical results for a Bessel function times Legendre polynomials class integrals,” J. Phys. A 39, L293–L296 (2006). [CrossRef]

62. S. Koumandos, “On a class of integrals involving a Bessel function times Gegenbauer polynomials,” Int. J. Math. Math. Sci. 2007, 1–5 (2007). [CrossRef]

63. P. J. Cregg and P. Svedlindh, “Comment on analytical results for a Bessel function times Legendre polynomials class integrals’,” J. Phys. A 40, 14029–14031 (2007). [CrossRef]

64. S. Koyama, K. Furuya, K. Wakayama, S. Shimauchi, and H. Saruwatari, “Analytical approach to transforming filter design for sound field recording and reproduction using circular arrays with a spherical baffle,” J. Acoust. Soc. Am. 139, 1024–1036 (2016). [CrossRef]

65. Z. Gong, P. L. Marston, W. Li, and Y. Chai, “Multipole expansion of acoustical Bessel beams with arbitrary order and location,” J. Acoust. Soc. Am. 141, EL574–EL578 (2017). [CrossRef]

66. I. Barth and O. Smirnova, “Nonadiabatic tunneling in circularly polarized laser fields: physical picture and calculations,” Phys. Rev. A 84, 063415 (2011). [CrossRef]

67. R. Messina, P. A. Maia Neto, B. Guizal, and M. Antezza, “Casimir interaction between a sphere and a grating,” Phys. Rev. A 92, 062504 (2015). [CrossRef]

68. D. Braun and J. Martin, “Spontaneous emission from a two-level atom tunneling in a double-well potential,” Phys. Rev. A 77, 032102 (2008). [CrossRef]

69. A. P. Polychronakos and K. Sfetsos, “High spin limits and non-abelian T-duality,” Nucl. Phys. 843, 344–361 (2011). [CrossRef]

70. J. J. Blanco-Pillado and M. P. Salem, “Observable effects of anisotropic bubble nucleation,” J. Cosmol. Astropart. Phys. 2010, 007 (2010). [CrossRef]

71. K. F. Ren, G. Gouesbet, and G. Gréhan, “Integral localized approximation in generalized Lorenz-Mie theory,” Appl. Opt. 37, 4218–4225 (1998). [CrossRef]

72. G. Gouesbet, “Validity of the localized approximation for arbitrary shaped beams in the generalized Lorenz-Mie theory for spheres,” J. Opt. Soc. Am. A 16, 1641–1650 (1999). [CrossRef]

73. J. A. Lock, “Angular spectrum and localized model of Davis-type beam,” J. Opt. Soc. Am. A 30, 489–500 (2013). [CrossRef]

74. A. A. R. Neves, A. Fontes, C. L. Cesar, A. Camposeo, R. Cingolani, and D. Pisignano, “Axial optical trapping efficiency through a dielectric interface,” Phys. Rev. E 76, 061917 (2007). [CrossRef]

75. A. A. R. Neves, “Força óptica em pinças ópticas: estudo teórico e experimental,” Ph.D. thesis (Instituto de Física Gleb Wataghin, Universidade Estadual de Campinas (2006).

76. W. J. Wiscombe, “Improved Mie scattering algorithms,” Appl. Opt. 19, 1505–1509 (1980). [CrossRef]

77. A. A. R. Neves and D. Pisignano, “Effect of finite terms on the truncation error of Mie series,” Opt. Lett. 37, 2418–2420 (2012). [CrossRef]

78. J. R. Allardice and E. C. L. Ru, “Convergence of Mie theory series: criteria for far-field and near-field properties,” Appl. Opt. 53, 7224–7229 (2014). [CrossRef]

79. P. Török, P. Varga, Z. Laczik, and G. R. Booker, “Electromagnetic diffraction of light focused through a planar interface between materials of mismatched refractive indices: an integral representation,” J. Opt. Soc. Am. A 12, 325–332 (1995). [CrossRef]

80. A. Rohrbach and E. H. K. Stelzer, “Trapping forces, force constants, and potential depths for dielectric spheres in the presence of spherical aberrations,” Appl. Opt. 41, 2494–2507 (2002). [CrossRef]

81. E. Fällman and O. Axner, “Influence of a glass-water interface on the on-axis trapping of micrometer-sized spherical objects by optical tweezers,” Appl. Opt. 42, 3915–3926 (2003). [CrossRef]

82. K. C. Vermeulen, G. J. L. Wuite, G. J. M. Stienen, and C. F. Schmidt, “Optical trap stiffness in the presence and absence of spherical aberrations,” Appl. Opt. 45, 1812–1819 (2006). [CrossRef]

83. M. Jahnel, M. Behrndt, A. Jannasch, E. Schäffer, and S. W. Grill, “Measuring the complete force field of an optical trap,” Opt. Lett. 36, 1260–1262 (2011). [CrossRef]

84. I. A. Martínez and D. Petrov, “Force mapping of an optical trap using an acousto-optical deflector in a time-sharing regime,” Appl. Opt. 51, 5522–5526 (2012). [CrossRef]

85. I. A. Martnez, E. Roldán, J. M. R. Parrondo, and D. Petrov, “Effective heating to several thousand kelvins of an optically trapped sphere in a liquid,” Phys. Rev. E 87, 032159 (2013). [CrossRef]

86. M. Selmke, M. Braun, and F. Cichos, “Photothermal single-particle microscopy: detection of a nanolens,” ACS Nano 6, 2741–2749 (2012). [CrossRef]

87. W. Moreira, A. Neves, M. Garbos, T. Euser, P. S. J. Russell, and C. Cesar, “Expansion of arbitrary electromagnetic fields in terms of vector spherical wave functions,” arXiv:1003.2392 (2010).

88. W. L. Moreira, A. A. R. Neves, M. K. Garbos, T. G. Euser, and C. L. Cesar, “Expansion of arbitrary electromagnetic fields in terms of vector spherical wave functions,” Opt. Express 24, 2370–2382 (2016). [CrossRef]

89. A. A. R. Neves, “Optical forces in optical tweezers,” https://github.com/aneves76/OpticalForcesInOpticalTweezers (2019).

Analytical calculation of optical forces on spherical particles in optical tweezers: tutorial

Abstract

1. OPTICAL TRAPPING

A. Optical Force Models

2. MODELING

A. Description of the EM Fields

B. Scattering Coefficients

3. BEAM-SHAPE COEFFICIENTS

A. BSC for a Plane Wave

4. HIGHLY FOCUSED GAUSSIAN BEAM

A. Angular Spectrum Representation of Focused Gaussian Beams

B. Beam-Shape Coefficients of the Focused Gaussian Beam

C. Case Study: Coincident Microsphere and Beam Axes

5. OPTICAL FORCES

6. FINAL REMARKS

7. SUMMARY

Funding

REFERENCES

Cited By

Figures (1)

Equations (97)

Journal of the Optical Society of America B