Research ArticlePHYSICS

Topologically enabled optical nanomotors

See allHide authors and affiliations

Science Advances  30 Jun 2017:
Vol. 3, no. 6, e1602738
DOI: 10.1126/sciadv.1602738


Shaping the topology of light, by way of spin or orbital angular momentum engineering, is a powerful tool to manipulate matter on the nanoscale. Conventionally, such methods focus on shaping the incident beam of light and not the full interaction between the light and the object to be manipulated. We theoretically show that tailoring the topology of the phase space of the light particle interaction is a fundamentally more versatile approach, enabling dynamics that may not be achievable by shaping of the light alone. In this manner, we find that optically asymmetric (Janus) particles can become stable nanoscale motors even in a light field with zero angular momentum. These precessing steady states arise from topologically protected anticrossing behavior of the vortices of the optical torque vector field. Furthermore, by varying the wavelength of the incident light, we can control the number, orientations, and the stability of the spinning states. These results show that the combination of phase-space topology and particle asymmetry can provide a powerful degree of freedom in designing nanoparticles for optimal external manipulation in a range of nano-optomechanical applications.

  • topology
  • nanoparticles
  • light-particle interaction
  • dynamical motion
  • nano-motors
  • optical manipulation
  • opto-mechanics
  • optical angular momentum


In recent years, there has been an increased interest in topological aspects of optical phenomena, including topological order in photonic structures (115), optical vortices (1618), and other robust phenomena in optics [for example, optical bound states in the continuum (1921)]. Shaping the topology of a beam of light can enable unique abilities of light to manipulate matter. For example, light that carries angular momentum [spin or orbital (16, 17)] can cause objects to rotate (2226), in a notion first suggested by Poynting (27) and later experimentally demonstrated by Beth (28). To date, the most common and relevant uses of light to control nanoparticles necessitate manipulation of the phase front of the light beam. For example, tightly focused laser beams are needed to realize optical trapping (17, 2932), and specially shaped beams give rise to optical transport (33) and negative radiation pressure force (3436). However, these approaches can be particularly sensitive to scattering or have small active volume, thus limiting their applications.

Although shaping the topology of a beam of light has opened opportunities for particle manipulation, here, we propose a more fundamental, and powerful, approach. By tailoring the topology of the phase space of the light-particle interaction, we can facilitate new means of manipulation of the position and orientation of a particle. Combining the phase-space topology with the appropriate particle asymmetry further enables dynamics that usually require shaping of the beam of light itself.

This allows us to explore an intriguing idea: When a beam of light illuminates a particle and rotates it into a steady orientation, it is commonly expected that no additional angular momentum transfer would occur, unless the light itself carries some form of angular momentum. On the contrary, here, we show that the dynamical process governing the rotation of particles in optical fields contains rich physics: Specifically, even linearly polarized light that carries zero angular momentum can give rise to a steady state in which a particle is spinning indefinitely around its stable orientation. The topological nature of these stable orientations implies that one should expect their existence in a broad, experimentally relevant, region of the parameter space. This suggests new possibilities of controlling the dynamics of a particle that do not involve shaping of the light. Such an approach may be attractive because it makes the interaction less sensitive to unwanted scattering and relevant to many particles at once.

To achieve this, instead of tailoring the light beam, we turn to particle asymmetry as a degree of freedom. Specifically, we provide examples with asymmetric particles known as Janus particles, which consist of two mirrored (but different) faces that can have independent functionalities (22). Here, we show that a simple spherical Janus particle illuminated with light in the appropriate spectral range can become a stable nanoscale motor. When illuminated by a plane wave, such a Janus particle—consisting of a dielectric core and a thin metallic half-shell—exhibits rich rotational dynamics, which include the existence of precessing steady states in a light field that carries no inherent angular momentum. We find that the dynamics are characterized by a set of vortices that evolve according to basic topological rules. The emergence of spinning steady states is consistent with these rules, arising from topologically protected anticrossing behavior of the vortices of the torque vector field.

Furthermore, we characterize the complete phase space for the interaction of the Janus nanoparticle and identify all stable orientations—the attractors for rotational dynamics. We show that the wavelength of the incident light is closely coupled to the asymmetry and the size of the particle. Accordingly, the wavelength of the light beam allows for external control of the number and the location of these attractors. Moreover, the wavelength also determines the relative sizes of the corresponding attractor basins in phase space. With the right particle design, we can find critical spectral regimes where the dynamics are dominated by a single vortex. This enables the control and prediction of steady-state positions regardless of the initial orientation.

These observations imply that the combination of phase-space topology and asymmetry can provide a powerful degree of freedom in nanoparticle design that enables interaction not achievable by shaping of the light alone and applicable to many particles. These interactions could facilitate optimal external control in applications that include molecular winding, nanoparticle sorting, and in vivo manipulation (3740).


A plane, linearly polarized light wave of the form Embedded Image impinges on a Janus particle consisting of a dielectric (silica) core with one-half coated with a thin layer of gold, as shown in Fig. 1A. In the regime where the size of the particle is comparable to the wavelength of incident light, neither the ray-tracing nor the Rayleigh scattering approximation is applicable. Instead, we take the most general approach and obtain the force (F) and the torque (M) on a particle by integrating the Maxwell stress tensor around a boundary enclosing the particle (see Materials and Methods). In general, a linearly polarized light beam can exert torque on small particles (4144), an effect most frequently observed in high–aspect ratio metallic nanowires, whiskers, and spheroids where the induced dipole moment p tends to align the long axis of the particle with the polarization direction of the incoming beam (the torque arises from the tendency to minimize the dipole energy U = −p E). However, the particle dynamics under this torque are either transient (once the particle is oriented, the torque becomes zero) or restricted in phase space or an optical trap. By contrast, we analyze the entirety of the phase space of a Janus particle interacting with a plane light wave of constant intensity. This is required because scattering of light strongly depends on the particle orientation with respect to the incident beam. In our case, the orientation is determined by the angles that the apex of the gold cap (point P) makes with respect to the coordinate system (Fig. 1B). Angles θ,φ correspond to the polar (with respect to ) and the azimuthal angle (with respect to Embedded Image), respectively. Figure 1 shows an example of a 1-μm particle covered with a coating of gold (60 nm), where the scattered electric and magnetic fields are computed by a finite element solver (see Materials and Methods), assuming the incident wavelength of λ = 1064 nm and water (n = 1.33) as the ambient medium.

Fig. 1 Optical torque on a Janus particle.

(A) A linearly polarized plane wave (Embedded Image can exert force and torque on an asymmetric (Janus) particle (silica sphere with a gold half-cap). (B) Particle orientation (cap apex P) is completely determined by angles θ,φ relative to the fixed light source. Absolute value of the torque (C) and the cosine of the angle between the torque and the cap (D), for all possible orientations of the particle. Here, the wavelength of light is λ = 1064 nm. Circles indicate rotationally stable orientations where the torque is zero (gray) or the particle is spinning (colored). (E) Four rotationally stable points corresponding to a spinning particle (as viewed from the direction of the light source).

We visualize the torque on such a particle by plotting its magnitude and its angle (with respect to particle orientation) in Fig. 1 (C and D, respectively) for different orientations of the particle. The phase space parameterized by (θ,φ) is split into four regions, corresponding to the azimuthal angle φ ranges of [i, i + 90o) for i = 0°, 90°, 180°, and 270°. The regions are related because of the symmetries of the incident field (Embedded Image). The darkest regions in Fig. 1C indicate orientations where the net torque is zero and the particle is in equilibrium without spinning. Among these orientations, there are four (shown as gray circles, for θ ≈ 25° and 170° and φ = 0° and 180°) that correspond to stable rotational equilibria. These stable fixed points are along the high-symmetry direction, where the particle orientation vector P is in the xz plane (Py = 0). Because Embedded Image and k0, we refer to these orientations as belonging to the E-k plane (that is, the particle is left-right symmetric with respect to the beam polarization).

In addition to these high-symmetry equilibria, we also observe four other orientations (colored circles in Fig. 1, C and D) that are also steady states despite the fact that the net torque is not zero. This behavior is explained by examining the angle α between the torque M and the particle orientation vector P. Figure 1D shows the value of cos(α) = M P/|M||P| for all values of particle orientation. It is immediately seen that the four colored directions relate to the cases where the external torque is precisely aligned (or precisely anti-aligned) with the particle orientation. To better elucidate, we show the four respective steady-state orientations in Fig. 1E; these correspond to a phenomenon where a particle spins steadily with a fixed angular momentum while being illuminated by a linearly polarized plane light wave that has zero angular momentum (both spin and orbital). At this point, we briefly remark that this effect originates from an avoided crossing of two topological charges associated with the torque vector field and is closely related to the size of the particle relative to the incident beam wavelength. The fact that these special orientations happen to occupy the central regions of the (reduced) θ-φ map [for example, for the red point in Fig. 1 (C and D), values of θ,φ are close to 90° and 45°, respectively] appears to be accidental: As we vary the wavelength, the positions of these special orientations will emerge, shift, and disappear in accordance with basic topological rules.

We now turn to the dynamics of an asymmetric particle in a plane electromagnetic wave. Because the particle is roughly spherical, we approximate its motion by a set of differential equations for translation and rotation (4547)Embedded Image(1)Embedded Image(2)where x is the position of the particle, ω is its angular velocity, and m and I are the particle mass and moment of inertia, respectively. Time-averaged optical torque Mopt(P, λ) and force Fopt(P, λ) depend on the particle orientation P and the wavelength of light λ. The viscous drag is given by ct = 6πηR (for translation) and crot = 8πηR3 (for rotation), where R is the particle radius and η is the viscosity of the medium (48). The diffusion in the system arising from the Brownian motion is characterized by the Langevin stochastic terms (ξ, ζ), where Embedded Image and Embedded Image and ξ(t) is a set of independent Gaussian variables with zero mean and covariance ξi(tj(t′) = δijδ(tt′) [the same is true for ζ(t)]. Finally, the orientation of the particle (which is fixed in the frame of the particle) is coupled to its angular velocity through dP(t)/dt = ω(P, t) × P(t). This, together with Eqs. 1 and 2, is used to evolve the position and orientation of the particle in time. Note that these equations allow for precession. In addition, they take into account the relevant forces and torques that act on the particles (including Brownian translational and rotational diffusion) and are valid as long as the flows are laminar (where the expressions for the translational and rotational viscous drags are applicable) (48).

Numerically evolving the rotational equation in time (for all possible initial orientations), we observe that in the overdamped regime (|M| ≈ crotω), any initial orientation of the particle will equilibrate into one of eight directions from Fig. 1 (C and D) (marked by four colored points and four gray points). These directions are the attractors of the system, and the corresponding map of attractor basins is shown in Fig. 2A. Here, the entire phase space of our system (all possible initial θ,φ) is partitioned (and colored) according to the corresponding steady-state attractor. The four gray circles correspond to rotationally stable orientations where the net torque vanishes. In contrast, for the four equilibrium orientations where the particle is spinning (colored circles), the terminal angular velocity is given by ωT = |Matt|/crot (approximately Embedded Image when normalized to the intensity of the incident beam; Matt is the value of the torque at the attractor orientation).

Fig. 2 Rotational dynamics in a linearly polarized plane wave.

(A) Attractor map for the rotational dynamics of the system (λ = 1064 nm). Background colors indicate attractor basins—regions of initial (θ,φ) in the phase space that equilibrate into one of the eight orientations (attractors) in the damped regime. (B) Examples of orientation trajectories (solid lines) obtained by numerically evolving Eq. 2 for different parameters (cases 1 to 4). (C) Stability of attractors to rotational diffusion and the dependence of particle orientation on the beam intensity. Averaging over many realizations, we find that the probability of the particle to reach and remain at a stable point increases with the beam intensity.

The dynamic Eqs. 1 and 2 allow us to characterize the stability of the attractors. For weak beam intensities, rotational diffusion will constantly reorient the particle, forcing it from one attractor basin to the next. The characteristic time scale for rotational diffusion is given by τrot = 1/2Drot, where Drot = kbT/crot is the rotational diffusion coefficient. For a 1-μm particle in water at room temperature, such as in our case, τrot ≈ 0.34s. We note that a medium of different viscosity would showcase equivalent dynamics but on a different time scale, as dictated by τrot. To quantify the resistance of a spinning attractor to rotational diffusion, we introduce a variable X as follows: If after t = τrot the orientation is still in the same attractor basin, then X = 1; otherwise, X = 0, as sketched in Fig. 2C (initial orientation is at the attractor). For this definition, the value 〈X〉, averaged over many realizations, is giving the probability of the particle to remain trapped at the attractor without leaving the basin of the attractor during the time τrot. For example, for the red attractor in Fig. 2A (θ ≈ 90° and φ ≈ 45°), the dependence of 〈X〉 on the beam intensity is shown in Fig. 2C (right). We observe that, for intensities of ~2I0 or larger, the attractor is very stable to perturbation by rotational diffusion (I0 = Embedded Image). Naturally, as the intensity is increased, the particle orientation becomes less susceptible to diffusion, with the threshold generally depending on the particle parameters and geometry. Figure 2B shows example trajectories of the particle orientation in time, assuming various starting orientations with no initial angular velocity. The arrows represent time evolution. For sample trajectories 1 and 2 (solid blue and orange lines), the beam intensity is 10I0 and 103I0, respectively, and the ambient medium is water. We observe that, for a wide range of intensities, the particle orientation remains confined to its initial attractor basin and eventually reaches the corresponding attractor. However, a strong enough beam or a medium of lower viscosity may allow for additional dynamics, such as a particle’s orientation precessing about the attractor point (trajectory 3) in a periodic motion (on average) or even becoming unstable and escaping into a neighborhood basin (trajectory 4). Our numerical simulations indicate that the onset of instability occurs when ω0 ~ τ−1, where ω0 is the angular velocity corresponding to the average radiation torque, and τ = ρR2/η is the characteristic viscous time. However, at these high angular velocities, the applicability of Eqs. 1 and 2 is limited, and additional terms may need to be considered (48, 49).

We now discuss the topological properties of the steady-state attractor points. They give an intuitive understanding of the emergence of the attractors that correspond to a steady-state particle spinning in a plane-wave light field with zero angular momentum. Because linearly polarized light comprises both left and right circularly polarized light (of equal magnitude and nonzero angular momentum), a physical explanation for this behavior comes from the difference in coupling of the asymmetric Janus particle to the two circular polarizations. The damped regime is characterized by M(P) ≈ crotω (from Eq. 1), so we approximate the evolution of the particle orientation as Embedded Image and proceed to plot the streamlines of the vector field NM × P. These are shown in Fig. 3 for all possible values of the particle orientation. Naturally, the streamlines in Fig. 3 show the same steady orientations consistent with the attractor map in Fig. 2A and the torque plots from Fig. 1 (C and D). The rotational equilibria are the vortex centers of the vector field N in the (θ,φ) space. We can characterize these vortices by their topological charge (known in nonlinear dynamics as an “index”), defined by (50)Embedded Image(3)Here, l is the spherical line element, γ(l) = arg[Nφ(l) + iNθ(l)] is the angle of the vector N, and C is a simple path on the θ,φ surface that goes around the vortex center in the counterclockwise direction. Because the torque vector returns to itself after a closed loop, the overall change must be an integer multiple of 2π; consequently, q is an integer. Figure 3 shows the position and the charge of all vortices in the system for the wavelength λ = 1064 nm. The topology and the symmetries of our system dictate the properties of the streamlines and the attractors. First, because, according to design, the complete phase space (θ,φ) forms the topology of a sphere with a Euler characteristic of 2, the sum of all topological charges within this phase space must be equal to 2, namely, ∑iqi = 2. Second, a stable attractor must have +1 charge (see Materials and Methods). The statements above imply that the existence of candidates for stable attractors in the experimentally relevant region of the parameter space is topologically protected. Consistent with these statements, in our system we observe q = ±1 vortices: +1 charges correspond to attractors (green circles) and the unstable extremum points at the poles (θ = 0° and 180°, indicated by horizontal blue lines; here, the gold cap faces directly away/toward the light beam); −1 charges are saddle points, shown as red asterisks. Close inspection of the vortex map in Fig. 3 confirms that the sum of all charges in the system is 2 (the Euler characteristic for a sphere).

Fig. 3 Flow of the torque vector field.

Streamlines of the N = M × P vector that governs the evolution of the particle orientation P in the damped regime (see main text). Rotational equilibria are the vortex centers of the vector field N. The vortices are indicated with circles and asterisks. Each vortex is also characterized by its topological charge: −1 (asterisks) or +1 (circles). For the current parameters, the poles (θ = 0° and 180°) carry +1 charge and are shown in blue; +1 charges that are also the attractors are colored green. Here, λ = 1064 nm. Notice that the topological invariant Euler characteristic ensures that the sum of all the topological charges equals 2.

The dynamics of vortex creation and annihilation are best illustrated by varying the wavelength of the incident plane-wave light beam, with the polarization unchanged (but of course, other parameters could be varied just as well, such as the index of refraction of the surrounding liquid, the thickness of the cap, the radius of the particle, or even the shape of the particle). Without loss of generality, we proceed to analyze only one (φ ∈ [0o, 90o]) of the four equivalent regions; this reduced phase space is consistent with the symmetry of our half-coated Janus particle. Figure 4A shows the corresponding vortex map for the vector field N when λ = 532 nm (half of the original wavelength used in Figs. 1 to 3). We observe only two +1 vortices located at the poles (consistent with the Euler characteristic of 2). The top one (θ = 180°) corresponds to the gold cap facing the incident beam and presents a rotationally unstable equilibrium. The bottom one (θ = 0°) indicates that the gold cap is facing away from the incident beam. This vortex is also the only attractor in the system; that is, no matter the initial orientation of the particle, it will always equilibrate such that the gold cap is facing away from the beam. Once at that point, the particle will not spin because the total torque (due to symmetry and the lack of angular momentum of the incident light) must be zero.

Fig. 4 Vortex dynamics for the torque vector field.

Reduced phase space corresponding to θ ∈ [0o, 180o], φ ∈ [0o, 90o]. (A) Attractor map for λ = 532 nm. (B) Increasing the wavelength gives rise to a rich vortex map (B1, B2, and B3), with charge annihilation (B1 → B2, red arrow; see fig. S1) and the emergence of spinning attractors as the result of an avoided crossing of +1 charges (B2 → B3, center, gray arrows). Dashed lines indicate attractor basin boundaries. (C) As the wavelength is further increased, the spinning attractors begin to migrate back to the φ = 0° high-symmetry edge (C1 and C2) and, in an avoided crossing similar to before, finally return to it (C3 and C4, arrows). Another example of charge annihilation (C1, red arrow; see fig. S2) is shown.

As we increase the wavelength, a rich vortex map emerges. For λ = 808 nm (Fig. 4B1), we identify numerous +1 and −1 vortices, including several attractors (green circles). The dashed lines separate the basins that correspond to different attractors (similar to different colors that partition the attractor map in Fig. 2A); consequently, there is only one green circle in each of these basins. To verify that the Euler characteristic is indeed 2, one should account for the multiplicity of vortices along the edges of the reduced coordinate space; for each φ = 0° (90°) vortex, there is an equivalent φ = 180° (270°) vortex. Figure 4 (B1 to B3) shows the dynamics of vortex annihilation and the emergence of spinning attractors. For slightly increased wavelength (λ = 827 nm), we first observe an annihilation process of two +1 vortices and four −1 vortices into two −1 vortices. This occurs near the θ = 180° pole in Fig. 4B1 (red arrow) and is better visualized in fig. S1. We notice two +1 vortices approaching each other (indicated by gray arrows). A further small increase in wavelength first brings these charges closer together and finally repels them away from the high-symmetry edge (φ = 0°). One of the +1 vortices moves inside the shown reduced phase-space region; the other does the same but in the neighboring (reduced) region (dimmed green circle). Because the particle orientation leaves this high-symmetry direction, there is now no symmetry that requires the net torque to be zero, which explains the unexpected rotational dynamics, as we presented in Fig. 1. This gives rise to a spinning attractor—a result of the avoided crossing of two +1 vortices. Increasing the wavelength causes the spinning attractor to traverse the reduced phase space. One example is, of course, the previously analyzed case of λ = 1064 nm and the attractor map shown in Fig. 2A. However, for larger wavelengths, the spinning attractors begin to come back to the edge and, in a reversed way compared to the avoided crossing before, return to the high-symmetry direction. This process is shown in Fig. 4 (C1 to C3) (gray arrows). We also observe another charge annihilation event where two −1 vortices and a single +1 vortex (on the pole) combine to give one −1 vortex (also on the pole). This is shown in Fig. 4C1 (red arrow) and also better visualized in fig. S2.

A particularly interesting situation occurs in Fig. 4C2 (where the wavelength λ = 1562 nm). Here, no symmetric orientation (whether along the edges or on the two poles) represents a stable rotational equilibrium: The only four attractors in the system are spinning attractors, with no obvious symmetry associated with these orientations. Our numerical simulations indicate that the spinning attractors, due to their topological nature, appear in the entire wavelength regime from λ ≈ 833 nm to λ ≈ 1575 nm, guaranteeing the existence of at least four rotationally stable orientations at each wavelength in this range. Developing new analytical techniques to build on these observations could be an important direction of future research. This behavior also illustrates the sensitivity of the attractor locations and the corresponding basin areas to the wavelength of light. Consequently, this can be exploited for external control and particle orientation. In addition, because there is no requirement of shaping or focusing the light beam, many such particles can be manipulated simultaneously.


Our results show that the interaction between a plane-wave light field and a simple metal-dielectric Janus particle forms a phase space that can give rise to a new class of rotationally stable particle orientations. Full-wave simulations demonstrate the existence of steady-state orientations with nonzero net torque, resulting in a spinning particle even in the absence of the angular momentum of the external electromagnetic wave. Using dynamical equations of motion, with stochastic terms describing diffusion, we characterized the orientation space of the particle and established the criteria for the stability of these special orientations to rotational diffusion. These spinning attractors are associated with the vortices of the projected torque vector field and arise as an avoided crossing of two vortices with positive topological charge. Although we distinguish these findings from the effects of topological photonics found in band structures and edge states, they are nevertheless built upon topological concepts and connect the interaction of light with nano- and microparticles to a wide range of physical phenomena, including electromagnetic bound states in the continuum in photonic crystals (19, 20), topological defects (50), and general vortex physics (51). In a broader context, just as topological principles protect edge states in condensed matter systems with nontrivial topology, such as topological insulators or systems exhibiting the quantum Hall effect, here, similar principles are responsible for the behavior of attractors in the phase space as the parameters of the system are varied. In particular, by characterizing the complete phase space of the system, we showed that the wavelength of light can directly control the location of stable orientations and the size of the corresponding attractor basins. This could facilitate the control and prediction of steady-state orientations no matter the initial rotational state of the particle. Although the attractor map is particularly sensitive to the wavelength of incident light, other parameters in the system can also be used to tune the particle dynamics. As an example, in fig. S3, we show the basins of attraction for the case of a particle of different size and also a different ambient refractive index.

We remark that these spinning attractors are characteristic of axially symmetric particles. Although the presented topological apparatus is specific to these particles, we expect similar states in higher dimensions. Moving away from structures with axial symmetry could be an exciting direction of future research. In addition, we note that this kind of rotational dynamics should be demonstrable because the relevant parameters are well within the experimental capabilities of an optical trapping setup, for example. In this case, some care may need to be taken to address the radiation pressure on the particle (the dominant force in the system) as well as any optically induced thermophoretic force that could arise because of localized asymmetric heating. To this end, tailored microfluidic channels (52) could provide the necessary physical confinement without affecting the optical properties of the surrounding environment. Finally, our case of a simple spherical Janus particle in a plane light field is already sufficient to enable rich and controllable rotational dynamics. More complex scenarios could include multiple asymmetric particles with a tailored range of spectral responses for operation at several wavelengths (53). Ultimately, our work demonstrates that shaping the topology of the phase space of the light-particle interaction through particle asymmetry can provide a powerful degree of freedom in designing nanoparticles for optical manipulation.


The optical forces and torques acting on the Janus particle were evaluated numerically from the Maxwell stress tensor T. Their expressions are given by (54)Embedded Imagewhere the stress tensor is defined as Embedded Image, Embedded Image is the unit vector normal to the enclosing surface, and εm is the permittivity of the medium surrounding the particle. The stress tensor is derived from the scattered electromagnetic fields obtained using a finite-element-method solver COMSOL Multiphysics. To calculate the torques, the stress tensor was integrated on a closed surface surrounding the particle. Permittivity of gold used in the finite element simulation is obtained from the literature (55). For the stochastic equations of motion, the viscosity of water was assumed to be 8.9 × 10−4 Pa·s. For all presented results, the electromagnetic source was a linearly polarized plane wave of the form Embedded Image, and the particle size was 1 μm (the cap was 60 nm thick). The orientation of the particle was uniquely determined by angles θ,φ. In most of our finite-element simulations, we calculated the stress tensor and the torques for all possible values of θ,φ in steps of 5°; to better resolve the situations corresponding to charge annihilation and anticrossing, the angular step size was 2°. These discrete values of torque obtained from COMSOL simulations were used in a separate numerical stochastic simulation of Eqs. 1 and 2. For orientations corresponding to intermediate angles, the torque components were obtained by bilinear interpolation.

A vortex is considered stable (that is, an attractor) if it is enclosed by a contour where the vector N is directed “inward” at every point of the contour. This implies that, under the action of the optical torque, the orientation of the particle will not escape the enclosing region. Mathematically, this stability condition can be written as γ(l) − θ(l) ∈ [0, π] along the contour l, where θ denotes the angle of the tangent of the contour and γ is the angle of the vector N (Eq. 3). Naturally, for any closed contour, we have Δθ = 2π and Δγ = 2πq, implying that Δ(γ − θ) = 2π(q − 1). This, combined with the stability condition γ(l) − θ(l) ∈ [0, π], indicates that only q = +1 vortices can be stable.


Supplementary material for this article is available at

fig. S1. Vortex annihilation and conservation of topological charge.

fig. S2. Vortex annihilation and conservation of topological charge.

fig. S3. Attractor map for a larger particle/air as the ambient medium.

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.


Acknowledgments: We thank L. Lu for helpful discussions. Funding: This work was partially supported by the Army Research Office through the Institute for Soldier Nanotechnologies under contract no. W911NF-13-D-0001 and by the Materials Research Science and Engineering Centers program of the NSF under award no. DMR-1419807. The research of I.K. was partially supported by the Seventh Framework Programme of the European Research Council (FP7–Marie Curie International Outgoing Fellowships) under grant agreement number 328853—MC–BSiCS. O.D.M. was supported by the Air Force Office of Scientific Research under award number FA9550-17-1-0093. H.B. acknowledges support from the QuantiXLie Center of Excellence. Author contributions: All authors contributed extensively to this work, including the analysis, discussion, and interpretation of the results. O.I. conceived the idea and performed the main calculations. M.S. supervised the project. All authors discussed, revised, and approved the manuscript. Competing interests: I.K., M.S., and O.I. are authors on a nonprovisional patent application related to this work filed with the United States Patent and Trademark Office (serial no. 15/245,641, filed on 24 August 2016; International Application no. PCT/US2016/048315). All other authors declare that they have no competing interests. Data and materials availability: All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. Additional data related to this paper may be requested from the authors.

Stay Connected to Science Advances

Navigate This Article