{'fraction_non_alphanumeric': 0.04284445759952642, 'fraction_numerical': 0.01465147254698831, 'mean_word_length': 4.443474758324382, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
The impact of microcavity wire width on polariton soliton existence and multistability

We have developed a model of the nonlinear polariton dynamics in realistic 3D non-planar microcavity wires in the driven-dissipative regime. We find that the typical microcavity optical bistability evolves into multistability upon variation of the model parameters. The origin of the multistability is discussed in detail. We apply linear perturbation analysis to modulational instabilities, and identify conditions for localisation of composite multi-mode polariton solitons in the triggered parametric oscillator regime. Further, we demonstrate stable polariton soliton propagation in tilted and tapered waveguides, and determine maximum tilt angles for which solitons still exist. Additionally, we study soliton amplitude and velocity dependence on the wire width, with a view to engineering quantum photonic devices.

We have recently developed a driven-dissipative, scalar (spinless) Gross-Pitaveskii model of the nonlinear polariton dynamics in realistic microcavity wires. In addition, we have numerically demonstrated multi-mode polariton solitons and a peculiar type of spatial multistability 15 . Under suitable conditions, different modes within the polariton soliton wave packet interact among themselves in such a way as to give rise to a self-localisation mechanism that prevents the pulse from broadening. Self-localisation phenomena in multi-mode systems arise due to counter-balancing nonlinearity with a combination of dispersive effects: material (polariton) dispersion due to a frequency-dependent dielectric response, cavity modes dispersion, and variation of the group velocity of each mode. In a recent work 16 , using the coupled-mode expansion technique we showed that a new component of the nonlinearity has to be taken into account -namely, the intermodal nonlinear coupling. The latter arises from two distinct nonlinear phenomena: intermodal cross-phase modulation (through Kerr nonlinearity) and degenerate four-wave mixing (through polariton parametric scattering).
In this work, we provide more details of our general mean-field, driven-dissipative Gross-Pitaevskii model. We apply it to the problem of polariton soliton propagation in microcavity wires with different widths -with a view to designing integrated polaritonic devices. Using the model, we study the soliton existence domain, multistability and group velocity dependence on the wire width, and identify optimum geometry parameters for soliton formation. We provide a qualitative explanation of the non-monotonic dependence of the maximum soliton amplitude and group velocity dependence on the wire width. Additionally, we investigate the dependence of the multistability on the cavity-emitter system's dissipation parameter, which includes exciton reservoir dephasing and cavity photon decay. Furthermore, we demonstrate numerically a new coherent propagation phenomenon of a radiating polariton soliton exhibiting periodic collapses and revivals, and identify regimes for existence of this new class of polariton solitons in a realistic microcavity wire. Finally, we simulate polariton soliton propagation in tilted and tapered microcavity wires, and determine the maximum tilt angle for which the soliton persists.

II. THEORETICAL MODEL
We consider a microcavity wire, fabricated by laterally etching a planar GaAs/AlGaAs microcavity, of the type described in 13 . The corresponding geometry is schematically shown in Fig. 1 (a). The wire is pumped by a CW laser linearly polarised in the y direction, such that it supports quasi-TE modes in the wire. The pump is inclined along the x direction at an angle θ contained within the x-z plane, so that it has a non-zero momentum along the wire. Similar to the case of a planar microcavity 12 , solitons propagating along the wire are triggered by a short seed pulse also incident at θ, having the same linear polarisation as the CW pump.
We describe the dynamics of the system by a scalar mean-field, driven-dissipative Gross-Pitaevskii model 15 : where E p is the normalized pump amplitude, g is the strength of exciton-exciton interaction 12 , ω R is the Rabi frequency in a planar homogeneous cavity, time is measured in units of T = 1/ω R , and E and Ψ are the averages of, respectively, the photon and exciton creation or annihilation operators. The normalisation is such that (ω R /g)|E| 2 and (ω R /g)|Ψ| 2 are the photon and exciton numbers per unit area, respectively. The unit length, L = h/(2m c ω R ), is determined by the effective cavity photon mass, m c . δ e , δ c , and the pump frequency, ∆, are detunings from a reference frequency, hω 0 = 1.55 eV (λ 0 = 800 nm), normalised with respect to the Rabi frequency. The reference frequency is blue-shifted with respect to the resonance frequency of the lower polariton branch for a given angle of incidence of the pump; the blue-shift is important because only then is excitation of a bistable mode possible over a finite range of E p 16 . γ c and γ e are cavity and exciton decay rates.
The lateral confinement in the cavity plane (along the y axis) is described by an effective potential U (y) in the photonic component, and a spatially confined normalised coupling, Ω R (y): where w is the dimensionless wire width and U bg is the background confinement potential. As demonstrated in our previous work 15 , the above super-Gaussian potentials give an accurate description of the dispersion of the lowest six guided modes of the wire, which includes the fundamental mode.  Figure 1: (Color online) (a) Scheme of the microcavity polaritonic wire structure and the excitation geometry. The CW pump and the seed pulse are incident along the x direction at angle θ, and have no components in the y direction. (b) Linear polariton dispersion for the fundamental and next five higher-order modes at zero detuning (∆ = 0), wire width w = 3 µm and a pump incidence angle θ = 20 • . Both the lower (LP) and upper (UP) polariton branches are shown. The planar microcavity dispersion is displayed with solid black lines for comparison. Note that ω is normalised with respect to the Rabi frequency. Amplitudes of the corresponding linear (c) photonic (E) and (d) excitonic (Ψ) electric-field modes.

A. Free-polariton mode dispersion and model parameters
The free-polariton dispersion is calculated from our mean-field model (Eqs. (1) and (2)). By setting E p = 0 in Eq. (1) and neglecting the nonlinear term in Eq. (2), using the ansatz E = A(y)e iκx−iωt , Ψ = B(y)e iκx−iωt , the dispersion of the linear modes, ω(κ), is given by the following eigenvalue problem: The resulting linear polariton dispersion for the fundamental and the next 5 modes is shown in Fig. 1 (b), along with the single-mode planar microcavity dispersion (black solid lines) for a microcavity wire width of 3 µm and a pump tilt angle of 20 • at zero pump detuning (∆ = 0). The planar microcavity can be regarded as the limiting case when the width, w → ∞. Thus, it is clear from the figure that, for the pump energy that we are using (blue-shifted from the LP resonance frequency), the planar microcavity dispersion in the reciprocal momentum space is steeper than that of the microcavity wire. The real part of the corresponding linear E-field and Ψ-field fundamental and higher-order mode profiles are shown in Fig. 1 (c) and (d), respectively. Using the finite-element electromagnetic-field solver, COMSOL, we obtain the model free parameters -namely the Rabi frequency ω R , the spatial scale L, and the effective potential depth U bg -by fitting the free-polariton dispersion with the numerically computed dispersion of a realistic microcavity wire structure, as described in 15 .

B. Stationary nonlinear modes
The stationary solutions of the nonlinear system, Eqs. (1) and (2), in the presence of a spatially homogeneous monochromatic pump are obtained by setting all time derivatives to zero, restoring the pump and the nonlinearity, and seeking stationary solutions to the homogeneous system of the form {E, Ψ} = A (y) e iκx , B (y) e iκx . The system then reduces to: After discretising the y co-ordinate, this system transforms into a set of coupled nonlinear algebraic equations, which is solved numerically for E and Ψ, e.g. using Newton-Raphson iterations for each fixed set of pump parameters E p and κ.
In what follows, all calculations are carried out for a microcavity polaritonic wire of width w = 3 µm at pump incidence angle θ = 20 • and model parameters as set out at the end of Sec. III A. We choose to plot the Ψ-field power, |Ψ| 2 dy, integrated in a transverse direction to the channel, in order to avoid any spurious non-monotonic behaviour which arises if the maximum or the channel midpoint value of |E| 2 or |Ψ| 2 are plotted. The integrated power of the Ψ-field, plotted at zero pump amplitude (E p = 0) and no dissipation (γ e = γ c = 0) for different detuning (∆) is shown in Fig. 2 (a). The number of confined modes for the specific potential considered is 5 (including the fundamental mode). Scanning through a range of pump amplitudes ( Fig. 2 (b, d)) shows that, depending on the detuning, a bi-or multistability family of curves is obtained for a range of dissipation parameters, γ. We cover the range 0 ≤ γ ≤ 0.04, as 0.04 (corresponding to Γ e = Γ c = 0.22 meV) is the most realistic value of γ 12 , and 0 is the limiting case of absence of dissipation. Note that, for both detunings considered, the bi-and multistability behaviour becomes increasingly pronounced as the dissipation is reduced. The curve in (d) clearly exhibits three distinct branches and two loops, as opposed to the commonly found single bistability loop observed in (b).
The multiple loops found is a new result specific to 1D microcavity wires, and is markedly distinct from the simple bistability observed in planar microcavities. Here below we analyse the origin of the multistability.
Each point of the bi-or multistability curve corresponds to a polariton (E, Ψ) mode. For the case ∆ = −0.1, only the first two (fundamental and first-order) modes of the family could potentially be excited (see left vertical dashed line in Fig. 2 (a) intersecting two modes). However, since the pump is spatially symmetric, the first-order mode does not arise, and therefore only the fundamental polariton mode (shown in Fig. 2 (c)) is observed in both the lower and upper bistability branches in Fig. 2 (b) -including at the turning point, denoted '1'.
At ∆ = 0, four modes could be excited (see right vertical line in Fig. 2 (a), intersecting four nonlinear modes). Note that, again, due to the even parity of the pump, only even modes survive -in this case, the fundamental and second order. The real and imaginary parts of the modes corresponding to each turning point (denoted '2' and '3') are plotted in Fig. 2 (e) and (f), respectively. By comparing Fig. 2 (d) with the intersection points in (a) for ∆ = 0, one can clearly see that these turning points correspond to a second-order and a fundamental mode, respectively. From this it could be inferred that the appearance of additional turning points may herald the appearance of additional stable states.
As we explain above, the predicted multistability is an interplay between the symmetric Gaussian pump and the number of modes supported by the microcavity wire potential. The number of supported modes can be controlled by selecting the confinement potential depth (U bg ) and width of the wire (w). The exact combination of modes can then be fine-tuned using the pump detuning.

C. Modulational instability
Once a stationary mode {A, B} is found, we can analyze its stability by introducing a small perturbation: where δ and λ are real, α is a scaling parameter, q is the perturbation wave vector, and the subscripts f and b stand for 'forward' and 'backward'-propagating waves, respectively. Assuming small amplitude perturbations ε f,b and p f,b , neglecting the second-order terms O ε 2 f,b and O p 2 the system can be recast as the eigenvalue problem: We find the eigenvectors x and the complex eigenvalues (δ + iλ) of the large sparse 4N × 4M matrix,M , each matrix element of which is a N × M block matrix, N and M being the number of grid points along the x and y directions (typically, N = M = 2048). We designate the names 'branch 1', 'branch 2' and 'branch 3' within the curves shown in Fig. 2 (d) for the consecutive sections with a positive gradient, starting to count from the bottom of the graph. We then plot the maximum of the eigenvalues' imaginary part, λ, against a small range of values of the q wave vector corresponding to each of these branches (Fig. 3), choosing ∆ = 0. In Fig. 4 (a), we plot the equivalent of the multistability curve from Fig. 2 (d) at γ = 0.04 for the integrated E-field power as a function E p . The stability of each point on the curve is worked out as follows. If, for a given E p , the relevant curve in Fig. 3 contains at least one point for which λ > 0 (as long as it is not the one at q = 0), then the corresponding point on the multistability curve is deemed to be 'modulationally unstable'. By inspection of Fig. 3 (a), this is the case for branch 3, and so the manifold of all points belonging to this branch is modulationally unstable (colour-coded in green in Fig. 4 (a)). This type of instability may nevertheless lead to the formation of stable localised entities, such as solitons. If, on the other hand, the curve lies below zero ( Fig. 4 (c)), the corresponding point on the multistability curve is classified as stable (such as is the case for branch 1). If the curve lies entirely above zero, or if at least λ(q = 0) > 0 ( Fig. 4 (b)), the corresponding point is classified as 'unstable' (branch 2). The modes (real and imaginary parts of the E-and Ψ-fields) corresponding to each of the multistability curve branches, are shown in Fig. 4 To verify the above linear stability analysis, we perform a dynamical check by solving Eqs. (1) and (2) in the time domain by the split-step method. First, we propagate the perturbed mode at E p = 0.068 and q = 1.2 as in Fig. 3 (a). Snapshots at t = 0.4 and 20 ps of the time evolution of the perturbed solution are shown in Fig. 5 (a) and (b), respectively. In agreement with the stability analysis above, this perturbed solution is modulationally unstable, and exhibits localisation and soliton-like features at later times ( Fig. 5 (b)). By contrast -and still in agreement with the above analysis -the perturbed solution for E p = 0.087 and q = 0.2 (corresponding to Fig. 3 Figure 2: (Color online) (a) Integrated power of the excitonic (Ψ) field at pump amplitude Ep = 0 and γc = γe = 0 for five nonlinear modes (F, fundamental; 1, first order; 2, second order; 3, third order; 4, fourth order) as a function of ∆, the detuning from the reference pump frequency, ω0. Vertical dashed lines help to determine which modes are present for different detunings. (b) Ψ-field power as a function of Ep for different dissipation parameters, γe = γc = γ, and ∆ = −0.1. A bistability curve is observed. Note that we take the pump to be spatially symmetric. The single folding point (denoted '1') corresponds to the fundamental mode (see panel (a) and text for more details). (c) Profile of the real parts of the photonic (E) electric field and Ψ-field of the fundamental polariton mode, corresponding to point 1 in (b). (d) Ψ-field power as a function of Ep for different dissipation parameters, γe = γc = γ, and ∆ = 0. A multistability curve is observed. The multistability arises from the fact that, for this detuning, two rather than one polaritonic modes are supported (see panel (a) and text for more details). The multistability curve touches the y axis at two folding points, denoted '2' and '3'. The bi-and multistability in (b) and (d) becomes more pronounced as γ is decreased -from 0.04 to 0. (e), (f) Real part of the E-field and Ψ-field of the nonlinear polariton mode profiles at ∆ = 0 at the folding points 2 and 3 in (d), respectively. Point 2 corresponds to a second-order mode, while point 3 corresponds to a fundamental mode.  Fig. 6 (a) and (b), respectively, is unstable. In this case, the time evolution exhibits filamentation, rather than localisation and soliton-like features. The stationary nonlinear mode profiles and perturbation eigenvectors computed from Eqs. 10-11 are given in 19 for both cases considered above.

ps in
The magenta curve above the third (topmost) branch in Fig. 4 (a) is the soliton branch. It is represented by the integrated E-field power along the transverse cross-section through the soliton maximum, which was computed by the 2D Newton method. The shaded area indicates the region in which solitons exist.

D. Dynamical modelling of soliton formation
We choose the pump amplitude within the interval of existence of solitons (see Fig. 4 (a), shaded area) and initialise the system with the stable mode from the lowest branch. The soliton formation is triggered at time t = 0 by a short (τ = 2 ps) Gaussian seed pulse with intensity FWHM of 3 µm. The pulse is collinear with the pump and provides perturbation to the stationary solution. We solve the system of equations, Eqs. (1) and (2), in the time domain by the split-step method. Following a transient process of initial pulse reshaping (data not shown), a soliton is formed. The stable soliton solution at ∆ = 0 and γ = 0.04 is displayed in Fig. 7 (a) for E p = 0.075. Longitudinal cross-sections of the soliton profile in (a) along the microcavity wire (y = 0) are shown for different simulation times in Fig. 7 (b). From these, it is clear that the soliton shape is preserved during propagation along the wire (see also a movie of the soliton propagation in 19 ).
In Fig. 8 (a), a transverse cross-section of the E-field power, |E| 2 , through the soliton maximum and the stationary nonlinear mode of branch 3 for E p = 0.075 are plotted on the same graph. The soliton is clearly more intense in the vicinity of that point than the mode of branch 3. In Fig. 8 (b), a transverse cross-section of the E-field power in the soliton tail, along with the nonlinear mode corresponding to branch 1 with E p = 0.075 are shown. Whereas the soliton and nonlinear mode differ significantly in the peak, they appear quite similar in magnitude in the tail.
It is worth noting that for E p = 0.075, the soliton temporal width (at FWHM), T p = 0.727 ps, is nearly one-half of the soliton temporal width in a planar microcavity (∼ 1.25 ps) 12 , and the soliton is squeezed along the microcavity wire.

E. Radiating solitons: collapses and revivals
In this section, we study multistability and soliton formation in a microcavity wire with the same width and incidence angle as before (i.e., w = 3 µm, θ = 20 • ), whilst keeping ∆ = 0 and selecting γ = 0.01. These conditions correspond to the multistability curve denoted by the red dashed line in Fig. 2 (d).
The multistability curves for the E-field and Ψ-field integrated mode powers, obtained from Eqs. (6) and (7), are shown in Fig. 9 (a) and (b), respectively. After performing a linear stability analysis, we identified the sections of the curve corresponding to stable, modulationally unstable and unstable nonlinear homogeneous solutions. However, in this case we could not find solitons of the type that we found for the case of γ = 0.04. Instead, we found an interesting   1) and (2) in the time domain using the split-step method. Snapshots of the time evolution of the perturbed mode's photonic electric-field intensity at (a) t = 0.4 ps and (b) t = 20 ps, after excitation with a Gaussian-shaped seed pulse of duration 2 ps at time t = 0. The perturbation wave vector is q = 1.2, corresponding to the local maximum of growth rate in Fig. 3 (a) for Ep = 0.068, lying on the upper branch (3) of the multistability curve in Fig. 4 (a). Clearly, the modulational instability on the third branch leads to localisation and concomitant formation of a soliton.  1) and (2) in the time domain using the split-step method. Snapshots of the time evolution of the perturbed mode's photonic electric-field intensity at (a) t = 0.4 ps and (b) t = 28.4 ps, after excitation with a Gaussian-shaped seed pulse of duration 2 ps at time t = 0. The perturbation wave vector is q = 0.2, corresponding to the local maximum of growth rate in Fig. 3 (b) for the pump amplitude, Ep = 0.087, lying on branch 2 of the multistability curve in Fig. 4 (a). Perturbation of the mode on the second branch clearly leads to filamentation, rather than localisation. Therefore, the mode is unstable.  type of radiating soliton in the vicinity of E p = 0.0341, which recovers after emitting radiation in a cyclical manner (see movie of the propagation of the radiating soliton in 19 ).
Regular snapshots of the time evolution of this radiating soliton is displayed in Fig. 9 (c). The resonant radiation, manifested in the appearance of a second peak in the trailing edge of the pulse, is likely due to higher-order dispersion effects, and leads to pulse damping. As a result, the soliton dynamics is described by a series of collapses, followed by revivals. After a pulse break-up and emission (see pulse profile at 107.5 ps, Fig. 9 (c)), the soliton is recovered (pulse profile at 125 ps). The revivals and collapses take place at regular intervals. An estimate for the time interval between soliton collapses, t c , can be derived from comparing the first and the last frames in Fig. 9 (c), giving t c ≈ 10-15 ps. Remarkably, the observed periodic dynamics is strongly dependent on the initial conditions, and in particular on the pump amplitude. Only excitations in a very narrow interval around the pump amplitude (approximately [0.025, 0.035], determined from our simulations) result in this peculiar dynamics.
This periodic behaviour is reminiscent of the collapses and revivals of semiclassical Rabi oscillations occurring when resonant coherent pulses propagate in and interact with a resonantly absorbing/amplifying medium. On the condition that the pulse duration is smaller than the relaxation times in matter, the resonant coherent pulse effectively becomes a polariton soliton, whose pulse area is conserved during propagation owing to the Pulse Area Theorem 17 . In the present case, the Pulse Area Theorem is indeed satisfied, since γ = 0.01 corresponds to a relaxation time T (= T 1 = T 2 ) ≈ 12 ps, which is greater than the pulse duration (2 ps).
It should be noted, however, that the above analogy with Rabi oscillation collapses and revivals is incomplete, since the revival component of the phenomenon cannot be explained classically (i.e., without introducing photon field quantisation) 18 . The semiclassical Gross-Pitaevskii mean-field model we have presented in this paper could only explain destructive interference effects, such as pulse collapses. We therefore find numerically that it must be the interference of the multiple discrete nonlinear modes within the microcavity wire that results in a semiclassical coherent effect mimicking the quantum dynamics of collapses and revivals. This interesting coherent propagation phenomenon deserves further investigation and should be a subject of future studies.

F. Impact of wire width on the linear polariton energy and mode group velocity dispersion
In this section, we investigate the dependence of the linear free polariton mode dispersion on the channel width. The linear energy dispersion, ω(κ), of the fundamental mode is shown in Fig. 10 (a) for different microcavity widths. Similar to Fig. 1 (b), the dispersion becomes steeper with increasing wire width. In agreement with this trend, the group velocity ( Fig. 10 (b), being the slope of the energy dispersion curves at each point) monotonously increases with increasing wire width. The upper bound is the group velocity corresponding to a single-mode planar microcavity.
Let us compare first the soliton velocity with the one for non-interacting polaritons (i.e., absence of a pump and nonlinearity) at a channel width w = 3 µm. The group velocity of the fundamental free polariton mode is calculated as the slope of the fundamental mode linear dispersion curve at pump wave vector κ/L = 2.6862 µm −1 , corresponding to an angle of incidence θ = 20 • . This gives v (0) g ≈ 1.472 µm/ps. This value is very similar to the soliton velocity calculated by the 2D Newton-Raphson method for E p = 0.075 (v g,sol ≈ 1.3 µm/ps), as well as to the one that can be inferred from the dynamical simulation through comparing consecutive time frames (v (0) g ≈ 1.1875 µm/ps). The free polariton group velocity monotonically decreases with the mode number (see slopes for lower polariton branch in Fig. 1 (b)). At pump wave vector κ/L = 2.6862 µm −1 , the first-order mode group velocity is v (1) g = 0.96 µm/ps; for the second-order mode, v (2) g = 0.549 µm/ps, or nearly one-half of the soliton velocity. For the 3rd and 4th-order modes, it decreases even further: v (3) g = 0.313 µm/ps and v (4) g = 0.227 µm/ps, respectively. Thus, we find that the soliton group velocity obtainable with this microcavity wire geometry ( Fig. 1 (a)) is always lower than that of a soliton propagating in a single-mode, planar microcavity (1.68 µm/ps) 12 .

G. Soliton branches at different channel widths
The soliton branches, computed by the 2D Newton method for different channel widths, are shown in Fig. 11. From inspection of Fig. 11 (a), (c) and (e), as well as data we have not shown, the interval of soliton existence is largest for channel widths in the interval [2.5, 2.81] µm. As the width is increased just beyond this range, in the vicinity of w = 2.86 µm a bifurcation point appears, Fig. 11 (c), and the singular soliton splits into two solitons. This point corresponds to the appearance of an additional branch -or equivalently, a second loop (see arrow in Fig. 11 (d))in the bistability curve, hence becoming a multistability curve. Upon increasing the wire width beyond w = 4.8 µm, the multistability curve acquires yet another loop (Fig. 11 (f), curve in magenta). However, note that we found no solitons for w = 6 µm (the 2D Newton-Raphson method did not converge).   In general, we posit that the occurrence of a bifurcation heralds the appearance of an additional stable light state. We will argue this in the ensuing passage.
The group velocity and the maximum amplitude of the soliton as a function of wire width are shown in Fig. 12. Our model predicts a non-monotonous behaviour of both quantities. From Fig. 12 (a), it is clear that the slowest soliton is obtained for the narrowest wire; however, we observe further local maxima and minima in the curve. Here we will attempt to give a qualitative explanation of the non-monotonous behaviour we have predicted.
Starting with w = 1.5 µm, as the width is increased, the soliton group velocity increases, reaching a maximum at w = 2.81 µm. At this point, we speculate that a second-order mode appears, whose group velocity is necessarily much lower than that of the fundamental mode (nearly a 2/3 reduction; see end of Sec. III F). This results in a decrease of the overall soliton group velocity. The soliton can then be viewed as a composite multi-mode entity, consisting of a fundamental and second-order modes. Further increase of the wire width leads to the bifurcation point at w = 2.86 µm in Fig. 11 (c), where we suppose the fourth-order mode appears. With a group velocity of only approximately 1/5th of the fundamental mode group velocity, composite soliton group velocity is decreased even further, reaching a minimum at w = 3 µm. As the wire width is increased beyond this point, the soliton group velocity increases again, as its constituent mode group velocities increase with wire width. This continues until another maximum is reached at around w = 4.4 µm. Beyond this point, we expect that a sixth-order mode appears, and the soliton velocity drops again. Following this line of reasoning, we would expect to encounter another bifurcation just below w = 6 µm, where a third loop appears in the multistability curve (see curve in magenta in Fig. 11 (f)). We expect that the soliton group velocity drops to zero at w = 6 µm, since we were not able to find a soliton for this wire width.
The dependence of the maximum soliton amplitude on the wire width is shown in Fig. 12 (b). It is very similar to the relationship between the soliton group velocity and wire width, and hints at a possible link between the two. Although at this stage we cannot comment on a plausible origin of this, it would be interesting to investigate the nature of this correlation in the future. In order to explain quantitatively this peculiar non-monotonous behaviour, we have recently developed a coupled-mode expansion method 16 . A detailed analysis using this method will be presented in a forthcoming paper.

H. Soliton propagation in tilted microcavity wires
We have performed a number of dynamical simulations with wires oriented at different tilt angles with respect to the x axis, for ∆ = 0 and γ = 0.04. In Fig. 13 (a), a microcavity wire tilted at an angle α = 5 • is shown for w = 3 µm, θ = 20 • and E p = 0.075. We trigger the polariton soliton formation using a seed pulse, as described in the beginning of Sec. III D. Following brief initial reshaping, the soliton travels undistorted along the tilted microcavity wire (Fig. 13  (b) and video in 19 ). For these parameters, we found that the soliton is preserved in amplitude up to tilt angles of approximately 6 • . Larger tilt angles, however, lead to a gradual decline in amplitude and eventual abrupt destruction of the soliton.
Our dynamical results above were checked through the 2D Newton-Raphson method. Using the latter, the depen-   dence of the soliton E-field maximum on the microcavity wire tilt angle is plotted for different pump amplitudes in Fig. 14 (a). The dependence is nonlinear and exhibits typical critical behaviour; the soliton disappears at different critical angles for different E p . The critical angle monotonically increases with increasing E p , reaching a maximum of α max ≈ 12 • at E p = 0.075. The soliton amplitude remains virtually unchanged up to 8 • . In fact, the soliton can persist up to almost 12 • with only a small loss in amplitude (about 2.5%). A transverse cross-section through the soliton maximum, calculated by the 2D Newton-Raphson method for α = 8 • and E p = 0.075, is shown in Fig. 14 (b). The soliton exhibits a pronounced asymmetric shape.
The maximum critical angle obtained by the 2D Newton-Raphson method differs by approximately 2 degrees from the one calculated by the time-dependent split-step method. The two methods are therefore in very good agreement with each other, and give a useful indication of the microcavity wire tilt limitations in designing soliton polaritonic splitters, routers and other functional components of future polaritonic integrated circuits.

I. Soliton propagation in tapered microcavity wires
From Fig. 12, it follows that the soliton amplitude increases with increasing microcavity wire width in the intervals w ∈ [1.5, 2.81] and [2.86, 4.5] µm. This effective soliton amplification effect can be exploited in the design of special, tapered microcavity wires, which could be used for a number of applications.
In Fig. 15 (a) and (b), the confinement potential of a tapered wire with initial and final widths w = 1.00 and 4.02 µm, respectively, is shown. The Rabi frequency potential is constructed in the same way as in Eq. (4), whereby the coupling is 1 within the channel and zero outside (data not shown). A soliton is launched from the left boundary at time t = 0. The time evolution of the photon component of the soliton power, |E| 2 , is shown in Fig. 15 (c) and (d) for two simulation times: t = 48 and 120 ps, near the beginning and end of the wire, respectively. Snapshots of the evolution of the soliton's longitudinal profile over time are shown in Fig. 15 (e). From Fig. 15 (c-e) it is clear that, as the soliton traverses the structure, it expands in the y direction, but holds together along x and does not dissipate. A movie of the soliton's travel through the structure is also provided in 19 . In Fig. 15 (e), the first (leftmost) profile shown is at a point near the left boundary, where the soliton is being excited by the seed pulse. After some initial reshaping, a soliton is formed. Due to the smooth change in the width of the wire, the soliton power increases adiabatically, reaching a local maximum at w ≈ 2.44µm. Subsequently, the power drops down sharply, reaching a minimum at w ≈ 2.7µm, whereupon it gradually increases again. As we shall show below, this non-monotonic behaviour is in very good quantitative agreement with the calculated by the 2D Newton-Raphson method dependence of the soliton amplitude on wire width ( Fig. 12 (b)).
We dynamically simulated the excitation and propagation of a soliton in a tapered wire with initial and final widths, w = 1.00 and w = 4.02 µm, respectively, at E p = 0.075. The soliton's maximum amplitude at each point along the wire is displayed in Fig. 16 (a) Fig. 11 (c)), then subsequently decreases, reaching a minimum, after which it increases again. The maximum amplitude of the soliton is also plotted in Fig. 12 (b). Figure 16: (Color online) (a) Dynamical simulation of the soliton excitation and propagation in a tapered microcavity wire (red curve), along with the dependence of the maximum soliton amplitude on the microcavity wire width, w, from Fig. 14 (b), calculated with the 2D Newton-Raphson method (blue curve). The pump amplitude is Ep = 0.075. Save for the initial soliton reshaping (first few points on the red curve), the dynamic simulation yields the same type of non-monotonic dependence on w as the 2D Newton-Raphson one. The shift to the left of the red curve along w is due to the excitation of a soliton at a point to the left of the first calculated point by the 2D Newton-Raphson method (at w = 1.5 µm). (b) Longitudinal FWHM of the soliton peak as a function of wire width for the same tapered wire as in (a). After formation, the soliton FWHM is not constant, but deviates by no more than 30%.
by 2D Newton-Raphson method (see Fig. 12 (b)) is plotted on the same graph for comparison. Following initial reshaping, the soliton reaches a maximum at a lower wire width than that predicted by the 2D Newton-Raphson method (2.81 µm). In spite of this, the overall trend of the two curves is very similar, and the maximum soliton amplitudes computed by the two methods are in a very good agreement. The lateral shift, in fact, is due to the soliton being excited dynamically at a point along the wire, w < 1.5 µm, preceding the first point calculated by the 2D Newton-Raphson method. Hence, the maximum is reached earlier, and the curve is translated to the left.
For E p = 0.075, there is no soliton at w = 4.8 µm (Fig. 11 (e), green curve). This is why we have plotted the dependence of the maximum soliton amplitude on wire width in Fig. 12 (b) only up to w = 4.5 µm.
The soliton FWHM for the same tapered microcavity wire is plotted in Fig. 16 (b). After the transient soliton formation process, the FWHM of the soliton peak taken in the x direction is not constant, but changes by no more than 30% by the time it reaches the end facet of the wire. The driving force behind this FWHM variation is currently unclear; however, this will be a subject of further investigation using our newly developed coupled-mode expansion technique 16 .

IV. CONCLUSION
We have developed a model of the polariton nonlinear dynamics in non-planar microcavity wires, based on drivendissipative, 2D mean-field Gross-Pitaevskii coupled equations.
For a realistic microcavity wire structure, we found that the typical for planar microcavities conventional bistability of the nonlinear mode power upon variation of the pump amplitude evolves into complex multistability. We have discussed the origin of the multistability in detail, and traced it back to the number of modes arising as a result of the interplay of three distinct factors. Firstly, the number of modes supported by the microcavity confinement potential is determined by the potential depth and width. Secondly, the pump frequency detuning from the reference frequency (blue-shifted with respect to the lower polariton branch resonance frequency for the relevant angle of incidence of the pump) selects a subset of these modes. Thirdly, depending on its parity, the spatial symmetry of the pump further reduces the number of available modes; for example, a symmetric, Gaussian-shaped pump leads to suppression of the odd modes and survival of the even modes.
Further, by employing linear stability analysis and conducting dynamical simulations, we discovered a range of pump amplitudes for which polariton solitons exist in a realistic structure operating in the triggered parametric oscillator regime. We have confirmed these results using the 2D Newton-Raphson method. In contrast to the singlemode polariton solitons in planar microcavities, polariton solitons obtained by our scalar model in microcavity wires exhibit a complex multi-mode structure. This type of spatially localised polariton solitons deserve special attention, since multi-mode solitons are capable of propagating further along an imperfect microcavity than their single-mode counterparts 16 . In view of this fact, we have developed a modal expansion method 16 to investigate the inter-modal nonlinear interactions and mechanisms of localisation, and attempt to find a means to leverage these multi-mode interactions and polariton dynamics for the development of polaritonic integrated circuits.
Bearing in mind recent advances in fabrication of microcavities with ultrahigh polariton lifetimes 7,8 , we investigated the possibility of exciting solitons in an environment with lower dissipation (i.e., when the exciton and cavity decay rates are lower). For these conditions, we numerically demonstrated an interesting type of radiating polariton soliton exhibiting collapses and revivals during propagation along the microcavity wire.
In order to gauge the plausibility of using microcavity wires for building Y-splitters, we examined the effect of tilted microcavity wires on the solitons' ability to propagate. Using dynamical simulations in tandem with 2D Newton-Raphson calculations, we identified pump amplitude-dependent maximum angles for which the maximum of the soliton amplitude was virtually unaltered ( 1% reduction); the limit appeared to be around 8 • . We also established the critical tilt angles above which the solitons vanished; the highest such angle was approximately 12 • .
With a view to creating novel polaritonic devices, we dynamically studied polariton soliton propagation in a tapered microcavity wire. We found that, as the polariton soliton travels over certain ranges of microcavity wire width, both its maximum amplitude and group velocity increase. This microcavity wire geometry could hence be exploited for reshaping and creation of solitons with specific characteristics, or simply to generate soliton amplification. For instance, tapered wire sections can be used as repeaters on a transmission line, compensating for signal loss at junctions or degradation due to structural imperfections. Furthermore, appropriate tapered sections could be selected according to whether fast or slow polariton solitons are required for a particular device application. For instance, fast polariton solitons are necessary for ultrafast polariton switches and modulators with applications in quantum information technologies -e.g., routing entangled single photons. Conversely, similar to the slow-light effect 20 , designing geometries supporting slow polariton solitons should lead to an enhancement of the nonlinear response, thereby decreasing the threshold for nonlinear effects, such as polariton parametric scattering and the Kerr effect. These slow polariton solitons may have applications as controllable delay lines and buffers.
In summary, this work lays the foundations for the simulation of the basic building blocks of future polaritonic integrated circuits, such as X-and Y-splitters, couplers and routers, based on soliton logic. Clearly, theoretical studies like these will be of great importance for the physical realisation of low-consumption, low-intensity threshold polariton integrated circuits of the future.