## Abstract

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

## INTRODUCTION

In recent years, there has been an increased interest in topological aspects of optical phenomena, including topological order in photonic structures (*1*–*15*), optical vortices (*16*–*18*), and other robust phenomena in optics [for example, optical bound states in the continuum (*19*–*21*)]. 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 (*22*–*26*), 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*, *29*–*32*), and specially shaped beams give rise to optical transport (*33*) and negative radiation pressure force (*34*–*36*). 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 (*37*–*40*).

## RESULTS

A plane, linearly polarized light wave of the form 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 (

**) 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 (**

*M**41*–

*44*), an effect most frequently observed in high–aspect ratio metallic nanowires, whiskers, and spheroids where the induced dipole moment

**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**

*p**U*= −

*p***⋅**

**). 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**

*E***) makes with respect to the coordinate system (Fig. 1B). Angles θ,φ correspond to the polar (with respect to**

*P**ẑ*) and the azimuthal angle (with respect to ), 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.

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* + 90^{o}) for *i* = 0°, 90°, 180°, and 270°. The regions are related because of the symmetries of the incident field (). 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 (

*P*

_{y}= 0). Because and

*k*_{0}∥

*ẑ*, 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

**Figure 1D shows the value of cos(α) =**

*P*.

*M***⋅**

**/|**

*P***||**

*M***| 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.**

*P*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 (*45*–*47*)(1)(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

*M*_{opt}(

**, λ) and force**

*P*

*F*_{opt}(

**, λ) depend on the particle orientation**

*P***and the wavelength of light λ. The viscous drag is given by**

*P**c*

_{t}= 6πη

*R*(for translation) and

*c*

_{rot}= 8πη

*R*

^{3}(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 and and

**ξ**(

*t*) is a set of independent Gaussian variables with zero mean and covariance ξ

_{i}(

*t*)ξ

_{j}(

*t*′) = δ

_{ij}δ(

*t*−

*t*′) [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

*d*

**(**

*P**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**| ≈

*c*

_{rot}ω), 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}= |

*M*_{att}|/

*c*

_{rot}(approximately when normalized to the intensity of the incident beam;

*M*_{att}is the value of the torque at the attractor orientation).

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/2*D*_{rot}, where *D*_{rot} = *k*_{bT}/*c*_{rot} 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 ~2*I*_{0} or larger, the attractor is very stable to perturbation by rotational diffusion (*I*_{0} = ). 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 10*I*_{0} and 10^{3}*I*_{0}, 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 τ = ρ*R*^{2}/η 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**c*

_{rot}

**ω**(from Eq. 1), so we approximate the evolution of the particle orientation as and proceed to plot the streamlines of the vector field

**≡**

*N***×**

*M***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**

*P*.**in the (θ,φ) space. We can characterize these vortices by their topological charge (known in nonlinear dynamics as an “index”), defined by (**

*N**50*)(3)Here,

**is the spherical line element, γ(**

*l***) = arg[**

*l**N*

_{φ}(

**) +**

*l**iN*

_{θ}(

**)] is the angle of the vector**

*l***, and**

*N**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, ∑

_{i}

*q*

_{i}= 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).

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 (φ ∈ [0^{o}, 90^{o}]) 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.

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.

## DISCUSSION

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.

## MATERIALS AND METHODS

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*)where the stress tensor is defined as , 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 , 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

**(Eq. 3). Naturally, for any closed contour, we have Δθ = 2π and Δγ = 2π**

*N**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 MATERIALS

Supplementary material for this article is available at http://advances.sciencemag.org/cgi/content/full/3/6/e1602738/DC1

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.

## REFERENCES AND NOTES

**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.

- Copyright © 2017 The Authors, some rights reserved; exclusive licensee American Association for the Advancement of Science. No claim to original U.S. Government Works. Distributed under a Creative Commons Attribution NonCommercial License 4.0 (CC BY-NC).