Spontaneous Formations of Dynamical Steady States in Polariton Condensates

We present a numerical analysis on dynamical steady states of polariton Bose–Einstein condensates (BECs) in an incoherent exciton reservoir driven by a ring-shaped optical pump. The balance between the loss and gain of polariton BEC induces a variety of steady states with different configurations, including approximate Gaussian distribution and topological defects, such as vortex–antivortex pairs, vortices with a winding number, and solitons. Besides, the system becomes unstable under fast decay rates and small pumping ring, where BECs can no longer exist in the long-time limit. We also confirm the soliton is dynamically stable in this system, with a steady polariton current induced by the repulsive polariton–polariton and polariton–exciton interactions.


INTRODUCTION
Polariton is a strong light-matter-coupled quasiparticle formed in semiconductor microcavities. Sufficiently strong optical pump can stimulate a spontaneous transition of bosonic polaritons to Bose-Einstein condensates (BECs) [1,2]. Many features of atomic BEC have been successfully demonstrated in polariton BEC, such as long-rang coherence [3,4], superfluid flow [4][5][6], and quantized vortices [7][8][9][10]. The combination of photon component with an extremely light effective mass and exciton component with strong interaction is essential for the high transition temperature [11][12][13] and fast thermalization in polariton BECs. Recent experiments have shown that polariton interactions can be enhanced by orders of magnitude through injecting an itinerant electron gas into a monolayer transition metal dichalcogenides [14,15]. Moreover, the nonlinear interaction between polaritons and reservoir excitons causes a variety of reconfigurable potential barrier for polariton BECs [16], such as parabolic potential [17], chiral polaritonic lenses [18], and ring-shaped potential [19,20]. Besides, rich vortex multistability [21][22][23] and Thomas-Fermi regime [24] revealed in the ring-shaped pumping potential provide opportunities to explore topological defects in nonequilibrium systems.
The polariton BEC is a dynamical steady state instead of a thermal equilibrium due to the rapid radiative decay of polariton and optical pump maintaining excitons. It provides a platform to study nonequilibrium quantum system with a macroscopic coherence. The leakage of photons from cavity and the decay of excitons via radiative and nonradiative processes set a finite lifetime of polaritons. The development of new structures can extend the cavity photon lifetime from 20 to 30 ps [25] to 200 ps or more [26,27]. Since the leaking photons carry the density, momentum, energy, spin, and phase information of polaritons, it is possible to characterize the nonequilibrium dynamics of polariton BECs via photon luminescence of the escaped photons [28,29].
The aim of this study was to investigate dynamical steady states of polariton BECs in incoherent exciton reservoir driven by a ring-shaped optical pump. Owing to the nonequilibrium nature of polariton BEC, the ring-shaped pump and decay can induce considerable effect on the steady states. First, the phase diagram of the steady states of polariton BEC by varying pumping and decay rates is depicted. The phase diagram consists of regions where the polariton BEC acquires an approximate Gaussian distribution, a vortex ring structure, a soliton structure, and a ring mode configuration. Second, different types of steady states in the phase diagram are, respectively, investigated by the density and phase distributions. Some exotic topological defects, such as higher order vortices and vortex-antivortex pairs, are then characterized and distinguished. And the effects of noise of initial condition and particle number on the dynamical evolution processes are discussed. Finally, the stability of solitons in the two-dimensional nonequilibrium system is analyzed. We find that several steady states can form spontaneously inside the pumping ring due to particle current induced by the nonlinear polariton-polariton and polariton-exciton interaction.
The remainder of this article is organized as follows. In Section 2, we present the mean-field formalism to incorporate the reservoir exciton. The phase diagram of steady states of polariton BEC is depicted and analyzed in Section 3. The features of different phases in the phase diagram are investigated in Section 4. Finally, we summarize the main results in Section 5.

THEORETICAL MODEL
At the mean-field level, the condensate dynamics can be described in terms of the Gross-Pitaevskii equation (GPE), which can be generalized to nonequilibrium systems by considering pumping and loss of polaritons with stimulated scattering [30,31]. As a result, an open-dissipative GPE (ODGPE) with complex terms for condensate wave function φ(r, t) can be introduced, together with a dynamical equation for the density n(r, t) of exciton reservoir. A coupled equation set is used to describe the complex coupling system of polariton BEC and exciton reservoir [32] as follows: Here, r (x, y) denotes the spatial coordinate and P(r) represents the pumping rate of exciton reservoir. In the following discussion, we choose a ring-shaped pump which can be written as follows: where J 6 is the first-kind Bessel function of the 6th order and the maximal pumping rate is fixed at P max 35.57 by the pre-factor. The first node of J 6 is located at radius R by tuning the parameter α. A small α corresponds to a large radius R of the pumping ring.
In an experiment, the pumping ring can be achieved by shaping the pulsed laser beam into a ring [24].
The exciton reservoir provides a source for creating condensed polaritons via stimulated scattering processes. The dynamics of polariton BEC is controlled by the following parameters: the decay rates of coherent polaritons c c and incoherent exciton reservoir c R , the stimulated scattering rate R R from the exciton reservoir into the polariton BEC, the nonlinear polariton-polariton interaction g c and polariton-exciton interaction g R , and an additionally pump-induced shift G [33]. In Eq. 1, the variables and coefficients are dimensionless by using the time scale τ 1/Γ c and the length scale of the diffraction length of polariton matter wave L Z/(m p Γ c ) , where m p is the effective mass of polariton and Γ c is the decay rate of coherent polaritons in physical units [22]. According to the typical parameters of a GaAs-based microcavity, we use the dimensionless parameters c R 2, R R 0.1, g c 0.001, g R 0.01, and G 0.01 [32,34] in Eq. 1 to numerically simulate the dynamics of polariton BEC by a time-splitting Fourier pseudospectral method. Such an algorithm is unconditionally stable, and of spectral accuracy in space and second-order accuracy in time [35], and has been widely used to study dynamics of vortices and solitons in BECs [36][37][38]. In the following discussion, all results are tested for convergence and stability.

PHASE DIAGRAM
Previous investigations on polariton BEC have revealed various types of spontaneous formations of steady states [39,40]. However, the boundaries between different steady states are ambiguous in the presence of balance between polariton loss and gain from reservoir. In this section, we discuss the steadystate solutions of polariton BEC by numerically evolving Eq. 1 with a different pumping radius R and the ratio between polaritons and reservoir decays c c /c R . In the seeds of the numerical evolution, we keep the ratio between initial numbers of excitons and polaritons more than 10 4 : 1 and induce random noises in both the exciton reservoir and polariton to simulate fluctuations. In the presence of a nonresonant ring-shaped pumping laser, more polaritons condense in this two-dimensional system through stimulated scattering of excitons. In the long-time limit, the total numbers of polaritons and reservoir excitons become constants upon time evolution, and their corresponding distributions are considered as the dynamical steady state.
In Figure 1, we show the phase diagram of steady states with tuning the pumping radius R and decay rate ratio c c /c R . According to the structures of density and phase distributions, the phase diagram can be separated into different regions, including approximate Gaussian distribution (AGD, red), vortices (yellow), solitons (green), ring mode (blue), and unstable region (gray). Here, the red circles label the phase transition points obtained by numerical simulation and the blue dotted lines are drawn by connecting neighboring circles to show the phase boundaries approximately.
We find that faster loss of polaritons (large c c /c R ) and slower replenishment from reservoir excitons (small pumping radius R) drive the system into an unstable region where no BEC can exist in the long-time limit. When the loss of polaritons is slow (c c /c R ≲ 0.6), steady states of BEC can form with increasing pumping radius. Specifically, the combination between loss and gain of polariton in the system induces a special region of vortex formation (yellow). Such a vortex state is under a subtle competition with the soliton phase (green) since the energies of the two states are very close in a broad region around the phase boundary, and the time evolution results are highly sensitive on the initial noise introduced to polaritons and reservoir excitons.

DYNAMICAL EVOLUTION
The density and phase distributions of polariton BEC in different steady states are investigated, respectively, in Figures 2, 3, 5. For polariton BECs with fixed decay ratio c c /c R 0.1, different steady states appear with increasing pumping ring radius R as in Figures  2A-D, corresponding to the black triangles in the phase diagram Figure 1. When pumping radius R is very small (R 1 in Figure 2A), the system presents no BEC in the long-time limit since the polariton experiences more severe loss than creation. As the pumping radius increases, more reservoir excitons supplement polaritons due to stimulated scattering, and a steady state of BEC starts to emerge as shown in Figure 2B. The density distribution inside the pumping ring tends to an approximated Gaussian profile as in the Thomas-Fermi regime. Similar phenomena have been used to measure polariton-polariton interaction strength [24]. By further increasing R, the systems transform to a vortex ring structure ( Figure 2C). A remarkable finding is that in this rotationally symmetric system, the vortex ring acquires a winding number |m| 2 as shown by the phase distribution of Figure 2C. Such a vortex solution spontaneously breaks rotation symmetry, where polaritons are squeezed to the edge of pumping ring, and phase singularity is at the ring center. When the pumping radius is increased to R 8 ( Figure 2D), the steady state is a soliton solution, where most density are localized in the center of pumping ring. Soliton and vortex are spatially localized solutions in nonlinear system, which can also present in nonequilibrium system.
High-order vortices are unstable in conservative systems, which will evolve into vortices with winding numbers m ± 1. However, it is possible to stabilize high-order vortices in a nonequilibrium system since the reservoir excitons continuously replenish polaritons [23]. To  characterize such exotic topological defects in detail, we fix the pumping radius at a relatively small value of R 5 in Figure 3 and choose four typical values of decay rate ratio c c /c R to obtain the density and phase distributions of steady states. The choices are labeled by black dots in the phase diagram in Figure 1. By increasing the decay rate ratio c c /c R , the polariton BEC transforms from AGD ( Figure 3A) to vortex ring structures with winding number up to |m| 3 ( Figures 3B-D). In the high-density region of vortex, the phase forms a spiral geometry. We also find that the vortex with larger winding number has a larger vortex core. Besides, notice that the winding numbers in Figures 3B,C are positive and the vortices rotate counterclockwise, while the winding number in Figure 3D is negative with clockwise rotating vortex. This parity symmetry is broken by the fluctuation of initial noises of the polariton and reservoir excitons. We also note that high-order vortex can be easily found in numerical simulations when the decay rate ratio c c /c R is large and the pumping radius R is relative small in the vortex region (yellow area) of Figure 1, and the maximal winding number we have found is |m| 5 for the parameters chosen in this study.
To further characterize the effect of fluctuations in the initial condition, we study the time evolution of the total number of polaritons N p (t) φ(r, t)| 2 dr (solid lines), and reservoir excitons N R (t) n(r, t)dr (dash lines) are shown in Figure 4, using the same parameters as in Figures 3B-D. We find that the initial state only influences the transient dynamics of the system but has no effect on the steady states in the long-time limit. In fact, the dynamical steady state is quite robust against the initial conditions of the polaritons and reservoir excitons, and it is solely determined by the pumping radius and decay rates. Comparing the results with different decay ratios c c /c R 0.03, 0.08, 0.3 in Figure 4, polariton population N P first decreases very rapidly and gradually accumulate because of the stimulated scattering process from the exciton reservoir. By contrast, the overall behavior of the number of reservoir excitons N R is to decrease upon evolution. When dynamical steady states have been achieved with the balance between gain and loss of polariton BECs, a larger decay rate ratio c c /c R leads to a smaller polariton population and a lower density distribution φ(r, t)| 2 around the border of the pumping ring. As a result, the stimulated scattering process is reduced and more reservoir excitons are retained in the long-time limit.
If the exciton reservoir is pumped by a large ring-shaped laser with radius R 25, the steady states show quite distinct structures by varying the decay rate ratio c c /c R as illustrated in Figures 5A-D, which correspond to the solid squares marked in Figure 1. For the case of small decay rate ratio c c /c R 0.014 in Figure 5A, a polariton BEC with AGD appears within the pumping ring. A random distribution of polariton exists outside of the pumping ring, owing to the low decay rate of polariton and the outward particle flow induced by repulsive polariton-polariton and polariton-exciton interactions. By increasing the polariton decay rate, vortex-antivortex pairs start to emerge as shown in Figure 5B. In two dimensions, the spontaneous formation of vortex pairs is usually the result of interaction of the condensed particles with thermal phase  fluctuations [41]. We notice that thermal effects are not strong enough to excite vortex pairs in this mean-field system. But the interplay between disorder potential induced by polariton-exciton interaction and the nonequilibrium nature can excite vortex-antivortex pairs [8]. The number of vortex-antivortex pairs varies with decay rates, and the pairs can move instead of having a pinned position with evolution. As the polariton decay rate further increases, the soliton solution spontaneously appears with c c /c R 0.1 as in Figure 5C. In this case, a ring mode structure is present in the density profile, located at the exact position of the phase jump. At even higher c c /c R as in Figure 5D, a ring mode structure is present. Soliton is a self-localized, shape-preserving solution of nonlinear partial differential equations. Solitons in polariton BEC are inherently dissipative and different from the ones in closed systems. To confirm the stability of solitons, we investigate the radial density evolution of soliton as shown in Figure 6A, with R 25 and c c /c R 0.1. This choice of parameters corresponds to the red star labeled in the phase diagram in Figure 1. We can find that a soliton establishes very rapidly from an initial random distribution of polariton BEC and survives for a very long time. This indicates that a dynamical steady soliton can be excited in the nonequilibrium system. In Figure 6B, we choose a typical time of t 1000 and plot the radial density ρ(r) φ(r)| 2 (black solid line) and phase θ(r) (red dotted line) distributions of polariton BEC. Notice that the phase distribution displays sharp gradient at the exact positions of density minima indicated by arrows. The fast variation of phase implies the existence of a large superfluid velocity along the radial direction. In Figure 6C, we show the distribution of superfluid velocity v s (r) zθ(r)/zr (blue solid line) and the polariton current j(r) ρ(r)v s (r) (magenta dotted line). Notice that for r < R, both the superfluid velocity and the polariton current is negative, indicating a current flow toward the center inside the pumping ring. For r > R, both v s (r) and j(r) are positive and the polaritons are flowing to the exterior of the pumping ring. This observation is a direct consequence of the steady solution of soliton, which is induced by the nonlinear repulsive interaction of polariton-polariton g c and polariton-exciton g R .

CONCLUSION
In this work, we investigate the spontaneous formation of steady states of polariton BEC described by ODGPE under ring-shaped pump and incoherent exciton reservoir. By varying the pump radius and decay rates, we systematically map out the phase  diagram via the density and phase distributions of polariton BEC and investigate the features of different phases. Specifically, we distinguish the regions of approximate Gaussian distribution (AGD), spontaneous formed vortices, solitons, and ring mode phases. We also find a region where the polariton BEC is unstable in the long-time limit. We also verify the soliton to be dynamically stable in this two-dimensional nonequilibrium system, where a balance is established in the long-time limit between excited polaritons around the ring-shaped laser and the polariton current flowing away from the boundary.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material, and further inquiries can be directed to the corresponding author.