Control of Light by Topological Solitons in Soft Chiral Birefringent Media

In practically all branches of physics, different types of solitons, with a number of them enjoying topological protection, are found. Here we explore how one- and two-dimensional topological solitons formed by spatially localized continuous orientational patterns of optical axis in uniaxial birefringent media interact with light. These solitons, in the forms of one-dimensional twist walls and two-dimensional skyrmions, are controllably generated in thin films of cholesteric liquid crystals to introduce spatially localized patterns of effective refractive index. Laser light interacts with these solitons as quasiparticles or extended interfaces of different effective refractive indices seen by ordinary and extraordinary waves propagating within the liquid-crystal medium. Despite our system ’ s complex nature, our findings can be paralleled with the familiar phenomena of total reflection and refraction at interfaces of optically distinct media, albeit these behaviors arise in a medium with homogeneous density and chemical composition but with spatial variations of molecular and optical-axis orientations. By exploiting the facile response of liquid crystals to external stimuli, we show that the twist walls and skyrmions can be used to steer laser beams and to act as lenses and other optical elements, which can be reconfigured by low-voltage fields and other means. Analytical and numerical modeling, with the latter based on free-energy-minimizing configurations of the topological solitons, closely reproduce our experimental findings. The fundamental insights provided by this work potentially can be extended also to three-dimensional solitons, such as Hopfions, and may lead to technological applications of optical-axis topological solitons in telecommunications, nanophotonics, electro-optics, and so on.


I. INTRODUCTION
Multidimensional topological solitons-once introduced by Skyrme to describe subatomic particles with different baryon numbers-found their way into many disciplines of physics ranging from condensed matter to cosmology [1][2][3][4][5][6] and optics [7][8][9]. These topologically protected but continuous fields cannot be unwound by smooth deformations [10] and can exist as energy-minimizing field configurations in theories that possess chiral or high-order nonlinear free-energy terms [11,12] needed to overcome the stability constraints imposed by the Derrick theorem [13]. Recently, there is great interest in studies of multidimensional topological solitons in condensed matter, especially in magnets and in liquid crystals (LCs) [8,[14][15][16][17][18][19][20][21][22][23]. The interest in topological solitons in magnetic solids is partially driven by the solitons' potential for spintronics applications including racetrack magnetic memory devices [24,25]. Moreover, topological solitons in LCs have been used to form various photonic crystal lattices and diffraction gratings [26][27][28]. Both fundamental and applied interests in studies of topological solitons arise largely from their particlelike nature and from their ability to be morphed and reconfigured by weak external stimuli like fields and light [29][30][31]. The effective sizes and physical behaviors of these particlelike structures in chiral LCs can be controlled by confinement, chirality, external fields, laser tweezers, and so on [26,[32][33][34][35][36][37][38][39][40]. Topological solitons in LCs can be created and templated on demand using scanning laser beams [26,41], thus making them good candidates for tunable all-optical devices using the intrinsic optical anisotropy of their LC host media. At the same time, the idea of light propagation control with LCs distinguishes itself in the more general context of advanced LC-based applications in the fields of photonics, telecommunications, and refractive and singular optics [42], which go far beyond the more conventional display applications for which LCs are broadly known. Recent achievements in these fields include biocompatible microresonators and microlasers based on onionlike chiral LC droplets [43,44], light guiding and steering with light-written waveguides in nematic LCs [45], photopatterned large-angle refractive elements [46], multispectral generation of singular optical vortices [47], optical inscription of reconfigurable arrays of microscopic optical vortex generators in chiral LCs [48], and tunable broadband Bragg mirrors stable in a wide range of temperatures [49]. However, these advanced applications of LCs have not utilized topological solitons in the optical axis of the birefringent LC medium so far, to the best of our knowledge.
In this work, we experimentally demonstrate and theoretically model controlled guiding and lensing of laser beams by highly reconfigurable topological solitons in the continuous spatial patterns of the optical axis of uniaxial LCs. The localized nature of topological solitons including one-dimensional (1D) twist domain walls and 2D skyrmions allows for their well-defined interaction with light in the linear optical regime, an interaction similar to that with objects of effective refractive indices that are different from those of the surrounding medium; however, here this interaction emerges in a chemically homogeneous single medium solely due to spatially localized structures of the LC's optical-axis (director) orientation [19,50,51]. The susceptibility of LC topological solitons to external stimuli further allow for reconfiguration of their interactions with light. A set of theoretical tools that we introduce and use allows for an explanation of our experimental findings and can be extended for a more general description of topologyenabled LC optical elements such as beam deflectors, polarization multiplexers, and lenses.
Below (Sec. II) we introduce the structure and topology of the two classes of solitons that we study and describe the experimental and theoretical methods. In Sec. III, we explore linear interactions of light with translationally invariant 1D topological solitons characterized by a twist of the optical axis by π [32,52]. We show that the interaction of an optical beam with a translationally invariant topological soliton may be described by a generalized version of Snell's law, which we introduce in a closed analytical form and complement with numerically calculated Fresnel coefficients. In Sec. IV, we demonstrate controlled deflection and lensing of a beam by a 2D elementary skyrmion [8,52,53], which is consistent with a simple description of the linear optical interactions based on Hamilton's equations for the flow of the Poynting vector in the geometrical optics approximation. We also show robust control over beam deflection and lensing with an electric field. Finally, in Sec. V, we describe the potential impact of our findings on the study of solitons in fields ranging from nuclear physics to optics and cosmology and draw conclusions. Since various solitons are widely studied [54] and since 2D topological solitons have been recently realized in optical fields, our demonstrations of how topological solitons in the optical axis of birefringent media interact with light may create a new research paradigm for exploration and exploitation of interactions between solitons of different types in condensed matter and in optics.

A. Classes of topological solitons
1D and 2D topological solitons, which in one and two spatial dimensions are characterized by elements of the first and second homotopy groups, respectively, can be smoothly embedded in 3D space as localized LC director field configurations within a uniform far-field background [52]. The diverse configurations that embed lower-dimensional topological solitons can be greatly enriched by adjusting the confinement and surface boundary conditions of the LC system [52]. The anisotropy of the LC's elastic constants, cholesteric pitch, external fields, thickness of the gap between confining substrates, and substrate's surfaceanchoring strength can tune solitonic structural features [52]. We emphasize here that whatever the type of embedding, optical solitons cover by definition at least once the order parameter space and therefore always yield a strong effective refractive index contrast seen by a light beam, thus leading to a wide range of strong interaction behaviors. In our experiments and modeling, we focus on two classes of 1D and 2D topological solitons embedded in 3D space.
The first class corresponds to 1D topological solitons associated with a π twist of the director field along one axis orthogonal to the far-field background director n 0 [ Fig. 1(a)], which cannot be realized in vector fields and are specific to LCs because of their nonpolar nature. Such 1D structures can be embedded in 3D space as orientational fields with translational invariance in a plane orthogonal to the axis of twist [52]. Particular confinement and surface boundary conditions are needed to stabilize different possible structures, with the extreme cases being represented by twist walls that are fully continuous [52] or cholesteric fingers of the third type (CF3s) that have the twist wall terminating on two singular disclination lines [55] near the top and bottom confining substrates. For constant elastic anisotropy, cholesteric pitch, and thickness between confining substrates, the transition from twist wall to CF3 or vice versa may occur by adjusting the substrate's anchoring strength and orientation of director at the surface [52]. For strong surface-anchoring strengths, CF3s form with χ disclinations that facilitate matching of the 1D topological soliton with the fixed, homeotropic substrate boundary conditions representing a uniform, topologically trivial state of the field [see Fig. 1(c)]. For weak-to-moderate surface-anchoring strengths, the singular disclination lines of the CF3s escape from the physical volume of the sample, since the homeotropic orientation is not strongly enforced by weak boundary conditions on confining substrates [see Fig. 1(b)]. Note that in this article, we neglect all possible in-plane curvature of CF3 and twist walls and assume that the structures presented in Figs. 1(b) and 1(c) are invariant by translation along the direction normal to the plane of the figure.
The second class of solitonic structures that we study corresponds to 2D topological solitons called skyrmions and also elementary skyrmions, which are low-dimensional topological analogs of Skyrme solitons introduced in highenergy physics [see Fig. 1(d)]. Similar to twist walls and CF3s, these 2D structures possess different structural surface-anchoring-strength-dependent 3D spatial embedding, which in general takes the form of rotationally invariant orientational fields [8,52]. For strong homeotropic anchoring conditions, the 3D structure that embeds the elementary skyrmion is called a toron and in its simplest form is associated with two point defects near the substrates that match the topologically nontrivial soliton's structure to the uniform topologically trivial boundary conditions [8,16,56] [see Fig. 1(f)]. For moderate anchoring strengths, the singular point defects of the toron are pushed outside of the physical sample's dimensions and become virtual [52] [see Fig. 1(e)]. We call this nonsingular 3D structure a baby skyrmion, while the term elementary skyrmion is reserved to describe the underlying 2D topological soliton shown in Fig. 1(d).
All experiments are done with topological solitons in cells with moderately strong anchoring, with which relatively thick LC samples tend to yield torons and CF3s [52]. Nevertheless, twist walls and baby skyrmions are good approximate structural representations of CF3s and torons in the LC cell's midplane. As the interaction of light with the localized patterns of the optical axis mostly takes place in this midplane, we preferentially use baby skyrmions and twist walls-which can be described with simple Ansätze minimizing the system's total free energy-in numerical simulations and theoretical models. With one example, we also illustrate how the precisely simulated director field of a CF3-based on the minimization of the full Landau-de Gennes free energy-allows a theoretical explanation of an effect seen in the experiments (see the Appendix D). This approach allows us to explore different regimes for optical interactions with topological solitons, thus revealing how the choice of a particular 3D embodiment influences the optical properties of the underlying topological soliton.

B. Observation of topological solitons and flow of light
We observe optical interactions between LC topological solitons and linearly polarized Gaussian beams with an inverted optical microscope (Olympus IX-73) that includes a homemade sample stage [57], which enables an insample-plane optical coupling with a monochromatic Gaussian beam and through-sample-plane visualization of the topological solitons with polarized-light optical microscopy (POM). Topological solitons are selectively generated with high-intensity laser pulses according to the method described in Appendix A. The transmitted light through phase-retarding LC solitons between crossed polarizers [polarization axes P and A in Fig. 2(a)] and the light scattered from the in-sample-plane laser beam due to LC orientational fluctuations are collected through an objective and captured by a digital camera revealing the interactions between the beam and solitons. The sample cell consists of two almost parallel but slightly tilted glass plates terminated by a glass coverslip [58], which enforces via tangential anchoring the uniform frustrated cholesteric LC alignment at a common edge, as shown in Fig. 2(a). A coupling objective introduces a Gaussian beam with its focal point near the coverslip's plane. To ensure that the center of the beam interacts with the midplane of the topological solitons, we configure the beam's coupling angle and insertion location so that the beam is coupled into the equidistant center plane defined between the cell's two confining glass plates [ Fig. 2(a)].
In Fig. 2(b), we illustrate a CF3 drawn in the cell's yz plane at a specific azimuthal angle that causes beam deflection. We also selectively create and position torons near the coverslip for interactions with a Gaussian beam (not shown in Fig. 2). The incident beam is subject to lensing and deflection from the torons, whose diameters are adjusted with an electric potential applied across the gap between the glass plates (U RMS ¼ 0-0.8 V at 1 kHz), indicated by the alternating voltage source shown in Fig. 2(a). Our experimental system thus enables the study of the laser beam's interactions with topological solitons under conditions when the beam's intensity is sufficiently low to ensure only linear interactions observed with beam powers <1 mW and exclude nonlinear optical effects, which we observe with powers starting around 20 mW. The incident beam's power, width, and linear polarization are configured with the setup shown in Fig. 2(c). A 532-nm diode laser placed directly before a half-wave (λ=2) plate combined with a polarized beam-splitting cube and an adjustable neutral density filter (NDF) define the incident beam power (< 1 mW) and a fixed linear beam polarization. The insertion of a double-Fresnel-Rhomb prism after the filter tunes the beam's incident linear polarization by rotation about the prism's optical axis as shown with the dashed line in Fig. 2 Fig. 2(c)], where objectives (not shown) are used for the lenses in the telescope, permitted the tuning of the beam waist in the sample. Directly prior to the coupling objective, the beamwidth is also adjusted by an iris (not shown). Inside the cell, the beam's full width at half maximum is measured between 30 and 40 μm. Additional details about the experimental setup, liquidcrystal cell preparation, and laser-assisted generation of topological solitons are presented in Appendix A.

C. Beam propagation and imaging simulations
In this article, we employ theoretical tools based on Maxwell's equations-with an exact eigenmode expansion or geometrical optics approximation-to elucidate the general mechanisms of interaction between light and topological solitons. However, numerical simulations are also used in order to validate and complement theoretical results concerning the flow of light in our system. Full 3D simulations of beam propagation inside our LC sample are performed, with the main axis of propagation parallel to the sample plates. We use a custom-written wideangle-beam propagation method (BPM) that one of us (G. P.) recently designed [59]. Briefly, the wave equation for anisotropic lossless media is recast into a simpler Schrödinger-like equation of the type ∂ u E ⊥ ¼ iLE ⊥ , where by definition the transverse field E ⊥ is orthogonal to the main propagation axis u. The differential operator L includes several contributions-anisotropic diffraction, beam walk-off, and phase operator-but most importantly, it is obtained without any paraxial approximation. The numerical implementation of this Schrödinger equation ensures accurate propagation of beams with deflection angles as high as 20°with respect to the main propagation axis [60]. Note that the more conventional finite-difference time-domain method, which directly solves Maxwell's equations without any approximations, is not a realistic option for our system due to its huge size (typical box size of 40 × 100 × 200 μm).
To simplify the comparison between simulations and experimental images, we also simulate POM micrographs of the studied LC solitons using a simpler version of our BPM code with a paraxial approximation, where the main axis of propagation-identical to the optical path in a microscope (Fig. 2)-is normal to the sample plates [61]. In all simulations, we use the values of the material constants given in Table I.

III. OPTICAL INTERACTIONS WITH 1D TOPOLOGICAL SOLITONS
In this section, we explore the linear optical interactions of Gaussian light beams with π-twisted 1D topological solitons embedded in translationally invariant configurations in 3D space, as we describe above. In order to develop insight into the reflection and transmission of light impinging on such objects, we first provide a generalization of Snell's law and corresponding Fresnel relations for our system. We then complement and verify this procedure by full BPM simulations before its use in the detailed analysis of the experimental results.

A. A generalization of Snell's law for confined modes of light propagation
We show here how to determine the effective modes of propagation of light while a monochromatic beam is reflected or transmitted through a topological soliton invariant by translation in one direction. Specifically, we consider frustrated cholesteric LC samples confined between two plates with homeotropic boundary conditions. By convention, x is the normal to the sample plane, and z is the axis of invariance of an isolated topological soliton inside the LC sample. The background state of the LC far from the topological soliton is assumed to be invariant by translation along both y and z. The Poynting vector of the incident beam is parallel to the sample plates in the yz plane and generates, through interaction with the topological soliton, one or more reflected and transmitted beams, as represented schematically in Fig 3. The situation represented in Fig. 3(b) resembles very much, by design, Snell's law at the interface between two isotropic media. Despite the fact that our physical system includes complex 3D  [62,63]. The value of the refractive index n p is experimentally measured for typical crown glasses [64]. variations of the optical axis, we now demonstrate that the very simple laws of refraction and reflection between isotropic media can be generalized to the case of light beams incident on birefringent-medium topological solitons embedded in a uniform far field. We start with the wave equation for the electric field E in lossless and nonscattering birefringent media [65]: with ϵ the permittivity tensor. Following the formalism presented in Ref. [66], we can use the existence of an axis of invariance (z) to express monochromatic solutions (at angular frequency ω) of the wave equation (1) in the following form: Eðx; y; z; tÞ ¼ where k 0 ≡ ω=c ¼ 2π=λ is the wave vector in empty space. The z component of the wave vector (k 0 p z ) is therefore conserved in the system of Fig. 3 and should be the same for incident, reflected, and transmitted beams. Although the latter statement is quite powerful, it is insufficient to determine the general directions of propagation of possible reflected and transmitted beams since the x and y dependences are still unspecified in Eq. (2). Our approach to express the general form ofẼ ðp z Þ ðx; yÞ in Eq. (2) is based on an eigenmode decomposition of the wave equation. Although we could greatly simplify the problem at hand by ignoring variations along the x axis and assuming that the relevant eigenmodes are simple plane waves-as in the original law of Snell-we choose here to rely on an exact decomposition based on waveguide modes, since we see later that these modes can play an important role in our system. We split the whole sample into a series of slabs parallel to the xz plane and centered on the coordinates y l such that y lþ1 > y l , with l ¼ 1; …; N as the slab number. In each slab, we assume that the permittivity tensor ϵ is independent from y: ϵðx; yÞ ≈ ϵðx; y l Þ in slab l, an assumption that is rigorously exact in the continuous limit (an infinite number of slabs for N → ∞). The first and last slabs of this discretization always correspond to the background domains represented in Fig. 3, while the other slabs correspond to the soliton. We then find for the slab l: In Eq. (3), a subscript ðp z ; l; m; þÞ [resp., ðp z ; l; m; −Þ] corresponds to the mth forward-propagating (resp., backward-propagating) eigenmode in slab l with a z dependence set by p z in Eq. (2). Here, forward propagating (backward propagating) means that the renormalized y component of the wave vector p y is positive (negative), and the actual sorting scheme of the eigenmodes is of little importance for the present discussion. The vectors Ψ ðÁÁÁÞ are defined as the field amplitudes of the eigenmodes, and the coefficients a ðÁÁÁÞ correspond to the weights of each eigenmode for a particular input beam, which can be found recursively by imposing the continuity of the x and z components of the electric and magnetic field at the interface between each slab. Combining Eqs. (2) and (3), we observe that the phase of each eigenmode flows along a vector p ðp z ;l;m;AEÞ ≡ f0; p ðp z ;l;m;AEÞ y ; p z g in the yz plane. Although Eqs. (2) and (3) correspond to a formal analytical solution of the problem at hand, it contains little information since p ðp z ;l;m;þÞ y is unspecified for now. To go further, we therefore need to specify the permittivity tensor of the topological soliton and calculate the eigenmodes.
In the rest of this subsection, we choose to focus on a simple example of a translationally invariant topological soliton-the twist wall already described in Sec. II A. The background state of such a structure corresponds to a uniform director field n 0 parallel to the x axis, while the twist wall itself corresponds to a π rotation of the director field along the y direction [Figs. 3(c) and 3(d)]. In all numerical simulations and comparisons between theory and experiments, we use a simple topological Ansatz for the director field of twist walls, whose mathematical expression can be found in Appendix B.
Since we are interested in the effective modes of propagation for the reflected or transmitted beams far from the soliton, we need to calculate the eigenmodes only inside the background domains of Figs. 3(a) and 3(b). This computation corresponds to the classical calculation of the eigenmodes of an anisotropic slab waveguide (see, for example, Ref. [67]), which can be split into an ordinary (transverse-electric, TE) mode with a polarization orthogonal to n 0 and extraordinary (transverse-magnetic, TM) mode with a polarization parallel to n 0 . In all numerical simulations and experiments, the extraordinary index n e and ordinary index n o always have a greater value than the refractive index of the sample plates n p , so that there always exists confined waveguide modes for both polarizations.
For simplicity, we give the expression of only the field amplitude for each eigenmode inside the LC layer and drop the ðp z ; l; m; AEÞ subscripts in favor of a simpler subscript (α) expressing the polarization state α ¼ e, o (with e for extraordinary and o for ordinary). We express the mode profiles along x as a function of the renormalized coordinate X ¼ k 0 x, and assume that the LC layer exists in the interval x ∈ ½0; h with h the sample thickness; we find for the polarization α ¼ e or o: where Ψ ðαÞ p can be determined from the Maxwelldivergence equation (∇ · E ¼ 0) and is nonzero only when α ¼ e. In the equation above, κ ðoÞ ≡ 1 and κ ðeÞ ≡ n 2 p =ðn e n o Þ are dimensionless factors, and we define the ordinary and extraordinary polarizations u ðαÞ as as well as the refractive index contrasts Δn ðαÞ for ordinary and extraordinary modes as Finally, the renormalized frequency ξ is given by a transcendental equation depending on an integer m: For each value of the integer m, we define ξ ðα;mÞ as the unique solution of this transcendental equation. When the sample is thick with respect to the wavelength, a very good asymptotic expression of ξ ðα;mÞ can be used: The integer m, which we call the waveguide mode number in the following, physically corresponds to the number of intensity lobes in the thickness of the sample. For odd (even) m, the mode is symmetric (antisymmetric) with respect to the center of the cell. Only a finite set of waveguide modes exists for a given p z , since the square roots in Eqs. (4) and (5) become imaginary as m goes to infinity. More specifically, it can be shown that the modes are fully confined in the sample (with exponentially decaying fields in the glass plates) if and only if jξ ðα;mÞ j < Δn ðαÞ .
Combining the previous results for extraordinary and ordinary waveguide modes, we find that the renormalized wave vector p ðα;mÞ along which the phase of the mth waveguide mode with polarization α ¼ e, o flows can always be expressed as Using the conservation of p z demonstrated at the beginning of this subsection, we finally arrive at a generalization of Snell's law for the transformation of an incident beam with incidence angle θ i and effective index n i into a reflected or transmitted eigenmode with effective index n ðα;mÞ (m ∈ N and α ¼ e, o) and angle θ ðα;mÞ : This equation shows that, after interaction, one incident beam can produce as many as four modes per waveguide mode number m: a reflected ordinary mode, reflected extraordinary mode, transmitted ordinary mode, and transmitted extraordinary mode. We emphasize that this result is fully general and is obtained without any approximation. Additionally, the Poynting vector and the wave vector are parallel far from the soliton both for extraordinary and ordinary modes, which means that the angle θ ðα;mÞ in Eq. (12) indeed characterizes the direction of energy flow of the transmitted or reflected modes. Remarkably, our generalization of Snell's law does not depend on the internal structure of the soliton but only on the permittivity profile far from the soliton. Although this remark may seem counterintuitive at first, we recall that Snell's law expresses only the conservation of photon momentum parallel to an interface (p z in our formalism), and therefore does not depend on the structural details of the interface but only on the value of the refractive index far from the interface. However, the soliton's orientational field is extremely important when calculating the coefficients of transmission and reflection of each eigenmode-the so-called Fresnel coefficients. Since a closed form of the eigenmodes inside the soliton cannot be obtained in general, we do not propose an analytical formula for the Fresnel coefficients, which therefore need to be calculated numerically (see Appendix C for calculation details).
In addition to the generalization of Snell's law presented above, we propose a simple estimation of the onset of total internal reflection for the twist wall studied here. With a calculation of the eigenmodes inside an infinitesimal slab parallel to the xz plane and centered on the soliton, we find that all confined waveguide modes (i.e., mode with a real p y ) disappear when jp z j > n o and are replaced by exponentially decaying modes along the y axis. Using once again the conservation of p z , we conclude that an input mode with effective index n i and incidence angle θ i will be totally reflected by a twist wall when θ i > θ c , where θ c is the critical angle of reflection: This result is not exact and relies on the assumption that the optical axis is parallel to z at the center of the soliton, which is not true near the sample plates as can be seen in Fig. 3(c). The validity of this criterion is discussed later in light of experimental and numerical results. We close this theoretical discussion with suggestions on how to use our generalization of Snell's law in practice. Indeed, in experiments and numerical simulations, one never specifies a combination of eigenmodes but rather the transverse field profile and polarization of the input beam. In this article, we always work with input beams with a Gaussian profile and extraordinary (parallel to x) or ordinary (orthogonal to x) polarization. Since the waist of these beams is much bigger than the wavelength (typically 10 times bigger), it can be verified numerically that the effective index of propagation n i in Eq. (12) can be very well approximated by n e (extraordinary polarization) or n o (ordinary polarization). The situation is slightly more complicated for reflected and transmitted beams. As can be seen in Figs. 4(a) and 4(c), where the results of typical BPM simulations for a 40-μm-thick sample are shown, transmitted and reflected beams are diffracted along the sample normal x after interaction with the twist wall [68]-due to variation of the director field along x-with an initially Gaussian profile transformed into a combination of waveguide modes confined between the two plates. In practice, this multiplicity of eigenmodes is not a problem as long as the waveguide mode numbers of the transmission or reflection eigenmodes are small with respect to a critical dimensionless number m c ¼ ðn o k 0 hÞ=π. When m ≪ m c , it can be verified that ξ ðα;mÞ ≪ n o in Eq. (11) and n ðα;mÞ ≈ n α (α ¼ e, o) in Eq. (12).
In other words, for small waveguide mode numbers m, the generalization of Snell's law in Eq. (12) depends on the polarization mode α ¼ e, o only. This approximation is valid for the BPM simulation presented in Fig. 4, since the dashed white lines in Figs. 4(b) and 4(d)-which represent the directions of propagation of reflected or transmitted modes predicted by Eq. (12) assuming n ðα;mÞ ≈ n αperfectly agree with the directions of propagation of the reflected and transmitted beams simulated with BPM. We mention that the approximation n ðα;mÞ ≈ n α is fully equivalent to a simple description based on an effective plane wave associated with a single wave vector p. In Sec. III B, where we discuss our experimental results, we start with this approximation but also present experimental conditions where m can play an important role in our generalization of Snell's law and cannot be assumed to be small-in which case, the use of waveguide modes cannot be avoided.

B. Experimental results and discussion
We present experimental beam reflection (transmission) from (through) a CF3 in whose colors correspond to relative mode powers predicted by numerically calculated Fresnel coefficients for particular incidence angles (see Appendix C for the details of this calculation, which is based on the formalism of transfer matrix for plane waves assuming approximate invariance in the x direction). Note that the curve T o mode in Fig. 5(h) saturates at 90°above an incidence angle of 60°, where transmitted ordinary modes become fully evanescent, as we discuss in more details in Appendix C.
While two reflection and two transmission modes always exist for any incidence angle, relative powers associated with each mode can be negligible or zero, as illustrated by 100 (i)-(k), the sample is illuminated with a backlight to view the birefringent medium's soliton. In (e)-(g) and (l)-(n), the backlight is extinguished to capture the beam's optical interaction with the soliton. Each interaction, as indicated by the incidence angles as shown, is portrayed as a pair of POM images viewed with and without a backlight. In general, as the incidence angle increases, different mode powers become nonzero as shown with a transition from total transmission to a regime with weak transmission and dominant reflection in (h). The scale bar's length is 100 μm.
the mode coloration and corresponding mode power scale in Figs. 5(a) and 5(h). For this reason, experimental data represented by the data points are provided only for regions with sufficiently strong nonzero mode powers. In the case of Fig. 5(a), the reflected e mode (not shown) has a relative power that is identically zero except for a few-degree incidence-angle range close to 90°. The same approach applies to the reflected o mode (not shown) for Fig. 5(h). With experimental mode observation constraints in mind, we conclude from the plots in Figs. 5(a) and 5(h) that the generalization of Snell's law in Eq. (12) with small m predicts reflected and transmitted trajectories of an incident Gaussian beam from a CF3 consistent with experiments, which we use as an illustrative example. However, we emphasize once again that the generalization of Snell's law in Eq. (12) is also applicable to any translationally invariant soliton embedded in a uniform far-field background, as we explain in Sec. III A.
The micrographs for transmitted and reflected beam behavior from incident o-mode beams [Figs. 5(b)-5(g)] portray control over redirection of incident beams using an experimental system fully described by our approach using our generalization of Snell's law assuming n ðα;mÞ ≈ n α (α ¼ e, o) in Eq. (12). Similar to o-mode incidence, incident e-mode beam interactions with CF3s are portrayed with POM images in Figs. 5(i)-5(k) and with no POM illumination in Figs. 5(l)-5(n). Normal incidence interaction angles are indicated directly in the micrographs. For small angles of incidence (typically θ i < 55°), observed micrographs such as the ones in Figs. 5(i) and 5(l) portray an incident e-mode beam interaction that is qualitatively similar to that of its o-mode counterpart, with one or two transmitted modes and no reflected modes as predicted by our numerical calculation of Fresnel coefficients in Fig. 5(h). For angles of incidence above approximately 55°, the situation becomes more complicated, and two interesting effects can be experimentally observed.
First, Figs. 5(k) and 5(n) show beam interaction micrographs associated with an incidence angle θ i above the critical angle for total internal reflection θ c ≈ 60°; although our numerical calculation of Fresnel coefficients predicts that no light can go through the topological soliton when θ i > θ c , a small fraction of light leaks through the experimental CF3 as seen in Figs. 5(k) and 5(n). As rigorously simulated and discussed in Appendix D, this loss is attributable to light slipping above and below the center of the CF3, in the thin space between the CF3's disclinations and the confining plates, a region with no internal reflection constraint. Note that in Fig. 5, the ideal regime of total internal reflection predicted by the numerical calculation of Fresnel coefficients is associated with the red parts of the curves T o mode and T e mode above 60°.
Second, and in stark contrast to the case of o-mode incidence, we experimentally observe that when the angle of incidence is typically above approximately 55°, transmitted and/or reflected beams sometimes split into smaller beams-which we dub beamlets-with slightly different directions of propagation. For example, two very faintly transmitted (resp., brightly reflected) beamlets can be observed in Fig. 5(m) [resp., Fig. 5(n)], with their mean direction of propagation corresponding to the curve T o mode (resp., R e mode ) in Fig. 5(h), which displays all of the measured beamlets' reflection or transmission angles from experiment for each incidence angle as shown. This splitting of transmitted and reflected modes suggests that the waveguide mode number m cannot always be assumed to be small for the incident e mode, and therefore that multiple waveguide mode packets with slightly different directions of propagation can be generated in our system. In the next subsection, we carefully examine the role of the waveguide mode number m in order to get more insight into this splitting effect.

C. Beamlets and higher-order eigenmode corrections
Motivated by the existence of beamlets in our experimental system, we now discuss the conditions under which the waveguide mode number m cannot be neglected in our generalization of Snell's law in Eq. (12). We recall that m becomes relevant in this equation when ξ ðα;mÞ has the same order of magnitude as n o [or equivalently, m=m c ¼ Oð1Þ with m c ¼ n o k 0 h=π], i.e., when eigenmodes with a great number of intensity lobes in the thickness of the sample are generated after transmission or reflection. Since we show in Sec. III A that the maximum value of ξ ðα;mÞ for confined modes is Δn ðαÞ , this condition can be fulfilled only if Δn ðαÞ =n o ¼ Oð1Þ. Using the material constants given in Sec. II, we calculate Δn ðeÞ =n o ¼ 0.502 for extraordinary modes and Δn ðoÞ =n o ¼ 0.125 for ordinary modes. We conclude that the waveguide mode number m can always be neglected in Eq. (12) for ordinary modes-which explains why we do not see beamlets for ordinary input polarization in the experiment-but not for extraordinary modes.
Using BPM simulations, we determine three possible situations in which extraordinary modes with m=m c ¼ Oð1Þ are generated through diffraction in the x direction: (1) Thin samples in which only a limited range of values are accessible for m. Since the maximum value for ξ ðα;mÞ is Δn ðαÞ , the maximum value for m scales with the thickness; see Eq. (9). (2) Misaligned input beams that are not perfectly centered in the midplane of the sample (yz plane equidistant from the two sample plates; see Fig. 3).
In this case, the incident beam can interact with both the birefringent-medium structure of the sample and the confining sample plates during reflection or transmission, thus generating high-frequency modes with high m. (3) Input-beam full width at half maximum that is comparable to the thickness between the LC cell's confining glass plates. In this case, the beam can interact with the glass plates during reflection or transmission and experience a high-frequency maximum m. We simulate with our BPM code the reflection of an extraordinary beam by a twist wall in a thin sample (10 μm assuming situation 1 above) and observe that the reflected modes split into two distinct beams with the same extraordinary polarization detailed in the simulated POM image in Fig. 6(a). For completeness, we also show an experimental image in Fig. 6(c), which is obtained in the same sample as in Sec. III B and shows a similar splitting for reflected modes as well as a weak transmitted mode irrelevant to the present discussion. We perform an eigenmode decomposition of the numerical optical fields along two lines parallel to the x axis and centered on two points P1 and P2 [shown in Fig. 6(a)] associated with the two reflected beams and plot the power spectra of these two decompositions as a function of the waveguide mode number m in Fig. 6(b). From this plot, we deduce that the beam associated with P1 contains a majority of low-frequency modes with a dominant waveguide mode number m ¼ 4, while the beam associated with P2 contains a majority of high-frequency modes with a dominant waveguide mode number m ¼ 14.
Finally, we verify that our generalization of Snell's law predicts the correct directions of propagation for the reflected beams (dashed white lines in Fig. 6) by using Eq. (12) with the dominant waveguide mode numbers of Fig. 6. As can be seen, the agreement between predicted propagation directions and simulated beams is reasonably good.
We emphasize that each beam reflected from a twist wall in fact contains a multiplicity of waveguide modes, which affect the spreading and intensity distribution of the reflected beams emanating from their reflection point. Interestingly, it can be observed that the beam associated with P1 diffracts greatly in the y direction, with a transverse profile skewed toward the soliton. As shown by the associated power spectrum in Fig. 6(b), this reflected beam has non-negligible contributions with moderate waveguide mode number m ∈ ½6; 13, each associated with smaller angle of deflection with respect to the z axis according to Eq. (12). Our numerical study therefore validates our generalization of Snell's law as a powerful tool to study the effective mode of propagation in the system considered here, even in extreme cases where high-frequency eigenmodes are generated through interaction with the soliton and confining plates.

D. Toward topology-enabled optical devices
From our experimental and theoretical studies, we suggest design rules for future research and possible devices that practice the fundamental insights disclosed in this article. We envision at least a few trajectories for future investigations. First, we show that confining effects are quite weak for ordinary modes in the system studied here (Δn ðoÞ =n o ¼ 0.125 ≪ 1), which could become a problem if propagation of ordinary modes with low loss is desirable in the target system. This limitation can be easily mitigated by using confining plates with an index n p less than that of crown glass or alternatively by using a liquid crystal with a higher ordinary index n o . The latter case can be easily realized with negative birefringence LCs (n e < n o ), which could be employed for the total internal reflection of ordinary but not extraordinary modes, if such a feature is wanted [69][70][71][72][73][74][75][76].
Second, future studies could focus on the optimization of truly x-invariant patterns of optical axis to ensure total internal reflection and to mitigate the emergence of beamlets. For twist walls, as discussed in this section, such could be realized experimentally using patterned confining surfaces (see Ref. [77] as an example of a patterning method) or low anchoring strengths for large cell thicknesses h or small h with strong anchoring conditions. Patterning of confining surfaces could also be used to design other birefringent-medium patterns of optical axis to optimize the conversion of e to o modes or vice versa within the pattern and satisfy m=m c ≪ 1 for transmitted and reflected modes, thus yielding simple control over a unique beam without additional beamlets.
Finally, for devices such as all-optical logic gates, the existence of beamlets may be desirable for a multiplicity of information channels, and so efforts would focus on their generation and control. Here, miniaturization of the system would be a very desirable feature, since we expect that the lower number of waveguide modes due to thinner samples would help to control the generation of beamlets. Attention could be given to the splitting of low-and high-order waveguide modes occurring inside the topological solitons of this section, and in particular, examine if multiple internal reflections inside the soliton could be responsible for this splitting in the case of transmitted modes-as suggested by one of the referees. In addition to the pertinent parameters enumerated above for beamlets with e-mode polarizations, the generation or collection of beamlets could be controlled by adjusting-without any change in structure topology-the global orientation of the optical axis inside birefringent-medium structures by, for example, the application of an in-sample-plane (out-of-sample-plane) electric field for LCs with Δϵ > 0 (Δϵ < 0). Thanks to the time invariance of Maxwell's equations, such patterns of optical axis could be used also to collect propagating beamlets and convert them into simple Gaussian beam profiles with e-mode polarizations, as an example.

IV. OPTICAL INTERACTIONS WITH 2D TOPOLOGICAL SOLITONS
In this section, we study experimentally and theoretically the deflection and lensing optical properties of 2D elementary skyrmions embedded in the rotationally invariant structures introduced in Sec. II A. We first use the theoretical framework of geometrical optics to gain insight into the flow of light inside these structures and then use this knowledge to discuss our experimental results.

A. Theory: Ray tracing for skyrmionic structures in birefringent media
First, we specifically focus on the linear interaction between light and baby skyrmions, whose structure is represented schematically in Figs. 7(a) and 7(b). As we explain in Sec. II A, we recall that baby skyrmions can be observed in frustrated cholesteric LC cells treated for homeotropic boundary conditions with moderate anchoring energy and correspond to a two-dimensional vortex of twist in the yz plane, slightly modulated along the x axis by the anchoring potential of the confining sample plates. We assume a simple topological Ansatz for the director field of baby skyrmions, whose mathematical expression can be found in Appendix B. We emphasize that baby skyrmions are topologically equivalent to the torons that we study experimentally, with singular point defects replaced by virtual ones near the confining plates. As long as the optical interaction of a beam with a toron or baby skyrmion happens at the midplane of a sample, we expect identical scattering properties for both types of rotationally invariant topological solitons.
Torons do not have an axis of translational invariance that allows a simple decomposition into propagation eigenmodes as in Sec. III. To provide a general understanding of the effects observed when laser beams are impinging on such objects, we therefore need a different approach. Since the skyrmions that we experimentally study are much bigger than the wavelength of the light (typically 100 times bigger), we choose to rely on the geometrical optics approximation and use the ray-tracing formalism described by Sluijter et al. [78] and Poy and Žumer [79]. The advantage of the approach of Poy and Žumer is that all formulas are given in a coordinate-free covariant form and that it provides energy-conservation laws for the amplitude of the optical fields. Here, we focus only on the field lines of the Poynting vector-called rays-for the e-and o-mode polarizations, which can be traced using the following Hamilton's equations: where α ¼ e, o is the polarization state of the considered ray, r is the spatial position of a virtual observer-which we call a bullet-propagating along the ray, p is the wave vector of the ray renormalized by k 0 ¼ 2π=λ (the wave vector in empty space) and evaluated at the bullet's position,s is the optical length [80], and the Hamiltonians H ðαÞ for extraordinary and ordinary rays are defined as with ϵ ⊥ ≡ n 2 o , ϵ k ≡ n 2 e , and ϵ a ≡ ϵ k − ϵ ⊥ .
Since H ðoÞ does not depend explicitly on r (spatially uniform index of refraction n o ), one may observe that ordinary rays cannot be deflected by LC structures, which is why we focus exclusively on the behavior of extraordinary rays. We assume that these rays initially propagate in the yz midplane of the sample along the z axis defined in Fig. 7. By neglecting the variations of the LC director field nðrÞ along the x axis-an approximation valid as long as the ray is not propagating near the confining plates [see Fig. 7(b)]-and by developing Eqs. (14) and (15) to first order in ϵ a , we find This very simple system of differential equations for the extraordinary ray trajectories clearly shows that the deflection of the photon's momentum along y is controlled by a single quantity g: where z is the initial direction of propagation of light, and y is an axis orthogonal to z and parallel to the sample plates. By examining the sign of dp y =dz in Eq. (19), we find that when g > 0 (g < 0), extraordinary rays are deflected toward negative (positive) y.
Using this simple criterion on the birefringent-medium structure of Figs. 7(a) and 7(b), we can predict the general features of the light flow around baby skyrmions without any numerical calculations. The general shapes of extraordinary rays deflected by the skyrmion are represented in Fig. 7(c) based on the general variation of g [Eq. (20)] represented in Fig. 7(d). We observe that the outer part of the skyrmion acts as a repulsive potential, while the inner core of the skyrmion is analogous to a converging 2D lens. This simple qualitative approach is used to analyze experimental data and is complemented by the full numerical calculations of ray trajectories from the original system of Hamilton's equations (14) and (15), as well as more detailed BPM simulations.
Before closing this theoretical section and discussing our experimental results, we point out that the ray-tracing equations presented here [either the exact Hamilton's equations (14) and (15) or the simplified system (18) and (19)] are invariant by scaling. The deflection paths of the rays incident at a scale-invariant y=R, with R the typical radius of the skyrmion viewed in its yz plane, stay exactly the same while the size of the baby skyrmion is adjusted by a linear scaling operation. This remark is not true in very small systems with characteristic lengths comparable to the wavelength but proves to be reasonably accurate for the systems studied here.

B. Discussion of light deflection in experiments
The 2D ray-tracing model that we introduce above predicts multiple deflection regimes, which allow for a rich assortment of beam-steering and lensing behaviors for an incident e-mode beam, as portrayed in Fig. 8. Near the deformation edge of a toron, beam focusing and deflection are observed, as shown in the simulation of Fig. 8(a) and the corresponding experimental POM image in Fig. 8(d). As in Sec. III, observation of beam interactions are made possible by the scattered light from LC orientational fluctuations away from the local alignment director field interacting with the in-sample-plane laser beam. Focusing and deflection are enabled by the beam's incidence on a conceptual boundary region between g < 0 and g ≈ 0 as shown in Figs. 7(c) and 7(d) above. Rays that interact with the elementary skyrmion structure closer to the skyrmion's center experience a greater angular deflection, while rays that are far from the skyrmion propagate as they would in the absence of a deflection region. The combined effect on In (a)-(c), the white lines correspond to calculated ray trajectories, while the color images are simulated as in Fig. 4 with the beam propagation method. As depicted in Figs. 7(c) and 7(d), the angle of deflection is attributable to the sign and magnitude of the deflection parameter g. Small, negative values of g produce a slightly deflected beam whose simulation is shown in (a) with experimental confirmation in (d). When g switches sharply from a very negative value to 0, defocusing and enhanced deflection occurs, as shown theoretically in (b) and experimentally in (e). By virtue of a skyrmion's rotational invariance, a beam injected at the toron's center experiences symmetric broadening behavior as shown theoretically in (c) and experimentally in (f). The inset in (c) shows an enlarged view of the ray trajectories near the toron's center. Differences in appearance between simulated (a)-(c) and experimental (d),(e) torons mostly result from director perturbation, only near the glass plates, used to pin a toron against diffusive motion within its sample plane and from experimental illumination with a coherent white-light source. The scale bar's length represents 100 μm. a cluster of rays, which approximate a Gaussian beam, shown as white lines in Fig. 8(a) is their simultaneous focusing and deflection, which are experimentally confirmed in Fig. 8(d).
We observe hybrid defocusing and deflection behavior that results from sequential, distinct deflection regions along the beam's trajectory through a baby skyrmion. In the case of Figs. 8(b) and 8(e), rays nearest the skyrmion's center propagate twice through a g < 0 region inside the skyrmion's periphery but also must pass through a g ≈ 0 region near the center. The combined effect results in a spreading of rays nearest the toron center, while rays that propagate through only the peripheral g < 0 region experience uniform redirection as shown by the simulation in Fig. 8(b). We observe concurrent experimental behavior as shown in the POM image of Fig. 8(e).
Injection of an incident beam at the toron's center results in a defocusing effect with symmetric beam broadening originating from a focal point inside the toron. While the beam center propagates through the elementary skyrmion and along the skyrmion-center boundary indicated by the dashed line in Fig. 7(c), each beam half must pass through deflection regions of opposite sign. The resultant broadening is shown in the simulation of Fig. 8(c). Therein, the defocused light appears to form a conelike multitude of distinct beams emanating from the toron's center. It can be shown that this discretization is due to an interference effect at the focal point [81], since a monochromatic beam is used in the simulation. While such a behavior is not distinguishable in the experimental POM image of Fig. 8(f), which demonstrates the defocusing effect from an elementary toron at no field, we experimentally observe this effect from a constricted elementary toron (not shown) subject to an electric field orthogonal to the sample plane. Since electric fields are known to dampen director orientational fluctuations [55], we hypothesize that experiments conducted with an electric field are similar to the ideal situation depicted in BPM simulations with a static director field. Conversely, in Fig. 8(f), orientational fluctuations cause dynamic random deviations from the equilibrium structure employed in the simulations. These deviations result in the destruction of the interference pattern shown in Fig. 8(c) and a blurred output signal as visible in Fig. 8(f). Differences in appearance between simulated Figs. 8(a)-8(c) and experimental Figs. 8(d) and 8(e) torons result from perturbation of the experimental torons' substrate director only near the glass plates, strong scattering and defocusing from the thermal fluctuations in a thick cell, and experimental illumination with a coherent white-light source. We emphasize that the pinning of torons with director perturbations is a surface effect used to prevent their diffusive motion (see Refs. [26,51]) within the sample plane. We have checked that identical deflection behavior is observed for both pinned and unpinned torons and notice that using incoherent illumination with a fully opened condenser aperture averaged out the surface perturbation due to the pinning, as visible in Fig. 9(a), which has demonstrated the expected crossed double-axis symmetry. From these two observations, we deduce that the pinning of the toron does not affect the bulk field alignments and therefore the topology of the baby skyrmion Ansatz.
The beam-center deflection angles from the simulations and experiments depicted in Fig. 8 may be included in a summary of all possible beam-center deflection angles, which are presented in Fig. 9. Each beam-center deflection angle θ may be mapped to the incident beam's normalized 9. Experimental measurements of incident beam deflection and lensing with a theoretical ray-tracing fit for a toron at no electric field. (a) A schematic representation of the coordinate system and parameters used for data measurement. The incident beam k i is deflected by an angle θ from the incidence trajectory, which is depicted as a white dashed line segment. The deflected ray k t may represent the center of a deflected beam profile or its half-intensity-maximum envelope rays. (b) Deflection angles for the beam envelope and beam center rays plotted against the normalized distance δ ¼ y i =R. The incident beam has an extraordinary polarization. The theoretical ray-tracing fit for center angles is plotted as the black curve. Negative angles in (b) correspond to a clockwise deflection of k t as indicated in (a). The scale bar's length represents 100 μm.
distance δ ¼ y i =R, where each parameter is depicted in the coordinate system shown in Fig. 9(a) over a toron's POM image. For a particular δ, theory indicates that a resultant beam-center deflection behavior will generally correspond to the ray paths shown for the deflection regions in Fig. 7(d) whose abscissa also depicts a normalized transverse distance. Solving Hamilton's equations in Sec. IVA for rays incident at various y i yields the solid black curve shown in Fig. 9(b). Corresponding experimental beam-center deflection data are shown as black filled circles in the same figure, with additional gray symbols representing the half-maximum beam envelopes. In general, experimental data are consistent with the theoretical curve except at the theoretical peaks near δ ¼ AE0.1. Though exploring the mechanism responsible for this deviant experimental peakdeflection behavior is outside of the scope of this study, we hypothesize that (1) Director fluctuations combined with the toron's optical thickness diminish sharp deflection behavior in favor of scattering. (2) Off-midplane rays could bias the experimental observation in favor of smaller deflection angles. Rays that interact with the toron outside of the toron's midplane and between the toron's two singularities would still interact with an elementary skyrmion topology; however, as a ray transgresses farther from the toron's midplane, it experiences a deflectionpotential pattern [Figs. 7(c) and 7(d)], wherein each region g has a smaller deflection magnitude. (3) Our simple comparison between the angle of deflection of the beam center and a single ray centered on the input beam neglects the possible shifting and/or scaling of the transverse intensity profile orthogonal to the rays, which can be accurately captured only with a full-wave simulation. We emphasize that such a small discrepancy is not a problem in our opinion, since the goal of our ray-tracing model is only to provide a qualitative description of the flow of light in our system. Should more quantitative data be needed-for example, in the design of specific beamsteering devices-one can always use full BPM simulations with the exact birefringent structures as input instead of Hamilton's equations (14) and (15).

C. Tuning of light-soliton interactions with an electric field
The nearly parallel-plate geometry of the LC cell [see Fig. 2(a)] enables the application of an electric field that is parallel to the uniform far-field director n 0 [along the x axis in Fig. 9(a)] [82]. As the LC has positive dielectric 10. Modulation of beam-center deflection by the application of an electric field along the far-field uniform director n 0 , which is orthogonal to the beam-propagation plane. (a) Incident extraordinary-mode (e-mode) beam-center deflection plotted against the normalized distance δ. Three datasets are collected for three different electric potentials as indicated. The solid curve is the behavior predicted by our simple ray-tracing model. For the beam-center deflection peaks at δ ¼ AE0.5, larger beam-center deflection is attributable to an asymmetric broadening and lensing of the beamwidth for different toron diameters, which map to applied potentials as indicated in (b). (b) Measured toron diameter plotted as a function of the electric potential for the experimental system. (c)-(e) Incident e-mode light deflected by three instances of a toron experiencing different electric potentials. The normalized distance of interaction δ ¼ 0.5 is constant as indicated by the star in (a). The shift in beam-center deflection is attributable to the contracting toron's distortion profile despite a constant incident beamwidth. Negative angles in (a) correspond to the clockwise deflection of the beam in (c)-(e), as defined in Fig. 9(a). The scale bar's length represents 100 μm.
anisotropy Δϵ > 0, the ensuing torque tends to orient the director field parallel to the x axis, causing a shrinking of the toron's effective diameter as the voltage increases. Our Hamilton's equations for ray trajectories predict size-invariant deflection and lensing behavior for rays of light interacting with a baby skyrmion. Like the simulations shown in Figs. 8(a)-8(c), the incident beam's intensity profile may be represented as a ray bundle. As the diameter of the toron shrinks under an applied electric field, the ray bundle width remains constant. At the same time, the distinct beam-deflecting regions depicted in Fig. 7(c) contract proportionately with the diameter of the toron, thus allowing the ray bundle to access multiple beamdeflection regions at large electric potentials but only one region at small or zero potentials. This is demonstrated in Fig. 10(a), wherein beam-center deflection is plotted against the normalized distance δ for the relaxed state (0 V) as well as for two nonzero potentials.
Comparison of beam-center deflection data taken in the relaxed state to data collected at applied potentials indicates sharper deflection at beam-deflection peaks, as indicated by the zero-electric-field theoretical prediction shown as the solid black curve. For the beam-center deflection peaks at δ ¼ AE0.5, larger beam-center deflection is attributable to an asymmetric broadening and lensing of beamwidth for different toron diameters, which map to applied potentials as indicated by Fig. 10(b). Therein, the experimentally measured toron diameter is plotted against applied voltages. At normalized distance δ ¼ 0.5 as indicated by the star in Fig. 10(a), the asymmetric beam-broadening and lensing behavior at different U RMS is shown in Figs. 10(c)-10(e).
Based on our experimental and theoretical studies, we suggest design rules for future research and possible devices that include 2D topological solitons. In our opinion, externally applied fields such as the electric field used in this section are essential to gain tunable optical devices. For enhanced control over beam deflection with reduced beam broadening, a LC with negative dielectric anisotropy Δϵ < 0 may be used to grow the diameters of baby skyrmions and torons with an applied electric field parallel to the symmetry axis of rotational invariance. Further control over deflection and beam lensing behaviors may be realized by patterning of confining substrates or weak anchoring in thick cells for baby skyrmions that are invariant along the x axis. Alternatively, with strong anchoring conditions, thin cells may be used to realize baby skyrmions [29]. Note that the in-sample-plane size of the birefringent-medium structures constrains the choice of suitable modeling frameworks; with large sizes, our raytracing formalism still applies, but for sizes comparable to the wavelength, BPM, or even the finite-difference time-domain method-which directly solves Maxwell's equations-must be used. For the case that lensing and broadening behaviors are desirable, a LC with a larger Δn ¼ n e − n o may be used to enhance lensing, broadening, and deflection according to our theoretical framework discussed in this section.

V. CONCLUSION
By exploiting the facile response of liquid crystals to external stimuli, we show that 1D and 2D topological solitons can be used to steer laser beams and to act as lenses and other optical elements, which can be reconfigured by laser tweezers and applied weak electric potentials. Analytical and numerical modeling, with the latter based on free-energy-minimizing configurations of the topological solitons, closely reproduce our experimental findings. First, we explore optical interactions with 1D topological solitons. We generalize the description of the observed optical transmission and reflection behavior due to such 1D solitons into a modification of Snell's law. We also explore and explain optical interactions with 2D topological solitons and demonstrate controlled deflection and lensing from a rotationally invariant toron. By the application of an electric field, we show tunable modulation of a toron's effective size and of enhanced deflection and lensing control by real-time tuning of the effective index that an incident beam experiences. The fundamental insights provided by our studies of 1D and 2D topological solitons potentially can be extended also to fully 3D solitons, such as Hopfions.
The existence of topologically nontrivial quasiparticles hosted in birefringent media that are accessible at the mesoscale allow optical interactions with them to be described either analytically or numerically, depending on the number of symmetries for the studied structure. Moreover, the optically interactive structures presented herein have no differences in density, chemical composition, permittivity, or permeability from the medium in which they are embedded. Their stability is assisted by topological protection. Because of this unique suite of characteristics, we envision potential technological applications.
While magnetic topological solitons (including 1D walls and 2D skyrmions) can be used as robust spintronic information carriers for applications like energy-efficient and high-density memory [21,83,84] their LC counterparts can potentially attract comparable technological interest because of their interaction with light described here, along with electric reconfiguration [26,38,85], spontaneous self-ordering into crystalline lattices [38,84], coordinated movement in schools [83,86], and diversity of states described with higher-order topological numbers [17,18,53,86,87]. One can also envisage optically writable and readable analogs to an information infrastructure realized by skyrmions and twist walls in magnetic systems. For example, by polymerization of a predesigned configuration of topological solitons, scattering from liquid-crystal director fluctuations can be significantly reduced with the formation of homogeneous or nearly homogeneous polymer domains [88]. This low-signal-loss environment would be ideal for applications like optical read-only memory even in the mesoscopic regime. Moreover, with sufficient system miniaturization, we hypothesize that both twist walls and skyrmions may be reconfigured in real time and using weak stimuli for alloptical transistor analogs that could form all-optical logic gates, circuits, information processors, reconfigurable memory, and information busses [89][90][91][92][93]. Finally, we posit that our fundamental studies may enable novel applications like evanescent beam steerers whose localized effective refractive index is controlled by the location of topological solitons and that may be useful as lidar systems [94]. Since our work provides a fundamental physical understanding of the interaction of light with topological solitons in optical-axis patterns of birefringent media, we envision that it will lay a foundation for a new subfield of informatics realized with topological structures thought previously inaccessible for facile configuration and manipulation of optical signals.

Observation of optical interactions
With reference to Fig. 2(a), a white-light source illuminates the back aperture of a condenser whose back-aperture iris is closed to a pinhole for coherent plane-wave illumination incident upon the sample. Interactions between topological solitons and incident beams are captured with a 4× objective and a digital camera (FLIR BFS-U3-51S5C-C). Meanwhile, a 2× coupling objective (∞ corrected, NA ¼ 0.06) is used to couple the in-sample-plane laser beam into the LC cell. We probe spatial configurations of topological solitons by moving the cell along the y axis of Fig. 2(b) to access linear segments of CF3s at multiple azimuthal orientations spanning 0°to 90°with respect to the y axis.

Liquid-crystal cell preparation
We first wash indium-tin-oxide-coated (ITO) glass, as described previously [95]. By a spin-coating and baking procedure, we apply a homeotropic polyimide substrate (SE-1211, Nissan Chemicals, Inc.) atop the ITO layer. Briefly, 50 μL of SE-1211 is applied to the ITO surface of the glass for spin coating. First, the plate is ramped to 700 rpm over 1 s. Immediately thereafter, it is ramped to 3100 rpm over 5 s and held at that speed for 30 s. The glass is then removed from the spin coater and immediately placed on a hot plate at 90°C for 2 min. It is then placed into an oven at room temperature. The oven is then ramped to 185°C over 30 min and subsequently held at the same temperature for 1 h.
To define the thickness between the glass pieces, we place ultraviolet-(UV) curable glue dots (NOA-65, Norland Optical) near the corners of the cell volume, containing silica spacer spheres (ThermoFisher), as shown in Fig. 2(b). Before stacking the second glass subsection, we affix 50-μm-thick UV glue dots at the corners near the coverslip and 60-μm-thick UV glue dots at the corners opposite those near the coverslip. We build the cell so that the SE-1211 substrates face inward toward each other. With weak UV exposure, the glue cures until hard (30 min). In a different embodiment of the cell, we lap the glass pieces to provide exposed surface area with an ITO substrate. We then solder electrical leads to the exposed ITO substrates.
We treat a washed microscope slide glass [95] with a nylon solution (Elvamide, DuPont) for tangential alignment via a different spin-coating and alignment procedure such that the tangential alignment is parallel to the x axis in Fig. 2(b). The coverslip made in this way imposes uniform alignment on the frustrated cholesteric LC director field at the beam-coupling interface so that the incident polarization suffers no phase retardation and thus maintains its initial state until interaction with a topological soliton. To prepare the nylon spin-coating solution, 0.15-wt % Elvamide pellets are added to pure methanol and vigorously stirred until dissolution, which occurs usually after 2 h at room temperature. An aliquot of 100 μL is deposited on the surface of a 25 × 75 mm 2 microscope slide. Immediately thereafter, the slide is spin coated by ramping to 4500 rpm over 1.5 s and then held at the same speed for 58 s. To establish tangential planar boundary conditions, the Elvamide surface is then slowly rubbed five times using a velvet cloth with an applied force of about 1 N. The slide glass is broken into approximately 3 × 25 mm 2 coverslips. We attach the coverslip with UV glue as shown in Fig. 2(b) and allow the assembly to cure with weak UV exposure until hard (approximately 30 min).
By using capillary forces, we fill an empty cell with a cholesteric LC composed of E7 (Shijiazhuang Chengzhi Yonghua Display Material Co., Ltd.) and a chiral dopant (cholesteryl pelargonate, SigmaAldrich) to define a chiral pitch according to the relationship C ¼ 1=ðHTPpÞ, where C is the weight fraction of the chiral dopant, HTP is the helical twisting power of the chiral dopant in E7, and p is the desired pitch of the cholesteric LC [56]. We adjust the LC's pitch so that the cell has a thickness-to-pitch (h=p) ratio of approximately unity near the beam-injecting coverslip. These particular cell specifications allow a uniform frustrated cholesteric LC background director field n 0 between the coverslip and selectively generate topological solitons farther away from the sample boundary. With reflection-mode POM, we confirm the uniformity of n 0 at the coverslip's interior surface.

Laser-assisted generation of topological solitons
We use near-infrared (1064 nm) laser tweezers described previously [26,51,96,97] to selectively generate topological solitons embedded in a uniform frustrated cholesteric LC background. At a nominal laser power of approximately 1 W a 10× objective (NA ¼ 0.30) focuses laser radiation on either the top or bottom substrate to deform the alignment of the substrate as anchoring points. We position two of these anchoring points at hundreds of micrometers apart to define a line.
With an in-sample-plane power of 30-160 mW, we translationally draw and anchor a CF3 along a line segment from one anchoring point to the other by locally heating LCs to the isotropic phase and quenching back to the nematic phase with the focused laser beam. These anchoring points fix the upper and lower π-twist-disclination termini in a CF3, which will otherwise shrink due to the tension in the singular disclination lines. Torons are generated with an approximately 1-s burst of nominally 190-mW laser radiation at any lateral position in the cell's bulk. When needed, we erase topological solitons by the application of an electric field parallel to the far-field n 0 of the LC's director. Because of the strong aligning effect of the electric field when applied to our LC with positive dielectric anisotropy [98], complete erasure of solitons occurs with an applied potential of U RMS ¼ 2 V at a frequency of 1 kHz. See Refs. [26,51] for a thorough characterization of the laser powers used to generate topological solitons within frustrated cholesteric LC systems.

APPENDIX B: THEORETICAL ANSÄTZE FOR THE TOPOLOGICAL SOLITONS
We present here the theoretical director fields that we use for the topological solitons studied in this article. These director fields are obtained by choosing simple topological Ansätze and minimizing the LC free energy (including elastic and anchoring energy) along the subspace spanned by the Ansätze. Note that the Ansätze presented in this section can be interpreted as a simple generalization of the theory in Ref. [99]. In this reference, the authors assume strong boundary conditions and derive simple Ansätze in isotropic elasticity, whereas here we go a step further and include the finite anchoring energy and the elastic anisotropy.
Based on Ref. [55], the total free energy of cholesteric samples treated for homeotropic anchoring can be written as where V corresponds to the volume of LC, S corresponds to the interfaces between the LC layer and confining plates, and the bulk (Frank-Oseen) and surface (Rapini-Popoular) free-energy densities are defined as ðn · ∇ × nÞ 2 þ K 2 qðn · ∇ × nÞ; ðB2Þ with K 1;2;3 representing the splay, twist, and bend elastic constants of Table I, q the cholesteric spontaneous twist, and W a the anchoring energy. We have checked that the Gauss elastic constant K 4 does not contribute to the elastic deformations of the Ansätze of this section. To simplify the mathematical expressions in the next subsections, we immediately introduce the rescaled Frank constants κ 21 ¼ K 2 =K 1 and κ 31 ¼ K 3 =K 1 , as well as the anchoring length ξ s ¼ K 1 =W a .

Ansatz for the twist wall
We approximate the director field of a twist-wall invariant by translation along the z axis with the following Ansatz: θðx; yÞ ¼ arctan sinh y ρ w ðxÞ ; ðB5Þ where x ∈ ½0; h is the direction normal to the sample plate, h is the sample thickness, and ρ w is a function to be determined through minimization of the total free energy. By inserting Eq. (B4) into Eq. (B1) and minimizing the resulting free energy, we find where we define the following dimensionless factors:

Ansatz for the baby skyrmion
We approximate the director field of a baby skyrmion with the following Ansatz, using cylindrical coordinates ρ ≡ ffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi y 2 þ z 2 p and ϕ (azimuthal angle such that cos ϕ ¼ y, sin ϕ ¼ z) around the rotational invariance axis x: where similar to the twist-wall Ansatz, x ∈ ½0; h is the direction normal to the sample plate, h is the sample thickness, and ρ s is a function to be determined through minimization of the total free energy. By inserting Eq. (B10) into Eq. (B1) and minimizing the resulting free energy, we find where we define the following dimensionless factors: with G Catalan's constant, and ζ the Riemann zeta function.

APPENDIX C: NUMERICAL CALCULATION OF FRESNEL COEFFICIENTS
The numerical calculation of the Fresnel coefficients necessitates the precomputation of the eigenmodes Ψ ðp z ;l;m;AEÞ introduced in Sec. III A and the recursive application of the continuity of the transverse-electric and -magnetic fields at the interfaces between the virtual slabs delineating the topological soliton. In this Appendix, we propose a simple way of estimating the Fresnel coefficients based on the assumption that the permittivity tensor does not depend on the x coordinate (direction normal to the sample): ϵðx; yÞ ≈ ϵðh=2; yÞ, where h=2 corresponds to the midplane of the LC layer. This approach is evidently quite simplistic since we show in the main text that waveguiding effects-caused by a discontinuity of the permittivity tensor in the x direction-can play an important role. Nevertheless, it has the advantage of giving the leading contribution to the transmitted and reflected beam powers due to the midplane permittivity profile at a low computational cost.
Based on this assumption, the eigenmodes inside each slab simply correspond to the usual extraordinary and ordinary plane waves of uniaxial media with uniform optical axis. By replacing the eigenmode index m with the polarization state α ¼ e, o, we find in each slab four eigenmodes with electric polarization vectors Ψ ðp z ;l;e;þÞ (forward extraordinary mode), Ψ ðp z ;l;o;þÞ (forward ordinary mode), Ψ ðp z ;l;e;−Þ (backward extraordinary mode), and Ψ ðp z ;l;o;−Þ (backward ordinary mode). The expressions of these vectors are given in a previous article [79] under a covariant form as a function of the optical axis n and the renormalized wave vector p ðp z ;l;α;AEÞ , which can be found in each slab using the eigenvalue equations also given in Ref. [79], provided the value of the conserved quantity p z is imposed.
Using Maxwell's equations to find the magnetic field as a function of the electric field, we then rewrite Eq. (3) of the main text under the following equivalent form for the electric and magnetic field in the lth slab: where χ ðp z ;l;α;AEÞ ≡ p ðp z ;l;α;AEÞ × Ψ ðp z ;l;α;AEÞ is the magnetic polarization vector of the mode ðp z ; l; α; AEÞ and y Ã l ¼ y l − h l =2 is the starting point of the lth slab along the y axis, with h l the thickness of the lth slab. The continuity of the transverse-electric and -magnetic fields at the interface between the slabs l and l þ 1 is then written as a ðp z ;lþ1Þ ¼ ½M ðp z ;lþ1Þ −1 M ðp z ;lÞ D ðp z ;lÞ a ðp z ;lÞ ; ðC3Þ By recursively using Eq. (C3), one can find a linear relation between the eigenmode amplitudes in the first (l ¼ 1) and last (l ¼ N) layers. However, one must be careful to take into account possible evanescent modes (associated with a complex p y ) which can exponentially amplify numerical noise. To avoid this numerical instability, we use the very reliable algorithm presented in Ref. [100], which recursively updates a modified transfer matrix T based on Eq. (C3) and the matrices D and M. The main output of this algorithm is the following relation: In plain words, the right-hand side of Eq. (C7) contains the amplitudes of the modes propagating toward the soliton, and the left-hand side of the same equation contains the amplitude of the modes propagating away from the soliton. This equation is causal-in the sense that only input modes are on the right-hand side-and therefore, evanescent modes (if present at all) can be only exponentially decaying.
Assuming a single extraordinary or ordinary input mode of amplitude 1 incident on the first layer, Eq. (C7) can be used to compute the amplitudes of transmitted or reflected extraordinary or ordinary modes, i.e., the so-called Fresnel amplitude coefficients. Alternatively, one can also compute the Fresnel power coefficients by using the y component (in-sample-plane direction normal to the soliton) of the Poynting vector of each mode. Up to a constant multiplicative factor, the Poynting vectors of the output modes can be written as S ðp z ;α;þÞ ≡ ½jaj 2 ðΨ × χ Ã Þ ðp z ;N;α;þÞ ; α ¼ e; o; ðC8Þ for transmitted modes, and S ðp z ;α;−Þ ≡ ½jaj 2 ðΨ × χ Ã Þ ðp z ;1;α;−Þ ; α ¼ e; o; ðC9Þ for reflected modes. In the last two equations, the index ðp z ; Á Á ÁÞ applies to a, Ψ, and χ . We then obtain the rescaled transmitted and reflected powers as P ðp z ;α;AEÞ ≡ Re½S ðp z ;α;AEÞ · e y Re½S i · e y ; where S i is the Poynting vector of the incident mode. We remark that for exponentially decaying evanescent modes in the y directions, the real part of S · e y is exactly zeroas expected from the conservation of energy. However, it may also be possible to have extremely small (but nonzero) transmitted power if the input propagating mode transforms into an evanescent mode inside the topological soliton (which is possible for grazing angles of incidence) before becoming again a propagating mode on the other side of the soliton. In practice, we find in all numerical experiments that this evanescent coupling between each side of the solitons is negligible if the thickness of the soliton is much bigger than the wavelength (typically more than 10 μm).
We use the formalism of this section to compute the reflected and transmitted powers assuming a single ordinary or extraordinary mode incident on a CF3 in a 40-μm-thick sample, i.e., the same topological soliton as in Sec. III B. The director field of this structure is obtained numerically (see Appendix D for more details), but we use only the permittivity profile in the midplane of the sample as we explain above (which is topologically equivalent to the permittivity profile of the twist wall in Fig. 3), which is discretized in a series of 40-nmthick slabs (yielding a total of 4800 slabs for the whole profile). The result of this calculation is shown in Fig. 11, and is also visible in Fig. 5 as color-coded variations of the transmitted and reflected powers. As can be seen, there is a clear range of incidence angles [incidence angle greater than approximately 60°, in agreement with the theoretical formula (13) presented in Sec. III A] for which an incident extraordinary mode suffers total internal reflection.

APPENDIX D: IMPERFECT TOTAL INTERNAL REFLECTION IN CF3
In this Appendix, we show how the total internal reflection regime-which is theoretically predicted and numerically observed for twist walls when the incidence angle is ≥ 60°-can be broken in CF3. To understand this, we must recall that the orientational field of a CF3 is topologically equivalent to the midplane director field of a twist wall, with two additional singular defects near the surface plates due to the strong homeotropic anchoring. We numerically minimize the Landau-de Gennes free energy, which takes into account the variations of the scalar order parameter in the defects' cores, using the method described in Ref. [101] and assuming a sample thickness of 40 μm and using the material constants given in Sec. II. The result of this calculation is shown in Fig. 12(a).
We simulate with BPM the propagation of a laser-beam incident on the CF3 at an angle of 70°with an extraordinary polarization. The associated simulated polarizing optical micrograph including the POM signal of the CF3 and x-averaged green intensity of the laser beam is represented in Fig. 12(b). Although the incidence angle is above the critical angle θ c introduced in Eq. (13) of Sec. III A, the reflection is not perfect here [contrary to twist walls; see Fig. 4(b)], and some light goes through the CF3. This effect can be easily understood by taking a look at the intensity profile of incident and transmitted beams along the sample normal x, as represented in Fig. 12(c): Since the director field is homeotropic inside a thin layer comprised between the bottom sample plate (white line) and bottom defect lines (red ellipse), the extraordinary light mode can propagate in this layer without suffering total internal reflection. Conversely, light is fully reflected at the midplane of the sample, as expected from the simple criterion of Sec. III A since the permittivity profile of a CF3 is topologically equivalent to the one of a twist wall at the center of the sample. We emphasize that the actual fraction of light which can slip through below the defects depends on the alignment of the input beam. In Fig. 12(c), we show a situation where the beam is slightly below the midplane of the sample. We have varied the vertical position of the beam in other (not shown) numerical experiments and have found that there is always some light transmitted below the defect, with an increasing (decreasing) fraction of transmitted light for lower (higher) input position between the two confining plates. The fraction of transmitted light can also be affected by other parameters such as the waist of the beam and a possible tilt of the beam with respect to the plane of the sample.  12. (a) Side view of the director field of a CF3 in a 40-μm-thick sample. The director field is obtained from a numerical minimization of the total free energy of the system, assuming strong homeotropic anchoring on the confining plates. CF3s are generally associated with two defect lines parallel to the soliton invariance axis (here, z). These topological defects are represented by two black spots near the sample plate in this figure. (b) Same simulated microscope image as in Fig. 4(b) but with a CF3 instead of a fully nonsingular twist wall. Contrary to Fig. 4(b), the incident extraordinary beam is not totally reflected, and some part of the beam is transmitted through the soliton, as also observed in experiments. (c) Side view of the beam intensity along the dashed black line (b), showing the intensity profile of the incident and transmitted beam along x (sample normal) and k i (incident wave vector). The defect lines of the CF3 are represented with two red ellipses, and the isolines n z ¼ −0.7 are represented with dashed white lines.