Controlling microresonator solitons with the counter-propagating pump

Considering a bidirectionally pumped ring microresonator, we provide a concise derivation of the model equations allowing us to eliminate the repetition rate terms and reduce the nonlinear interaction between the counter-propagating waves to the power-dependent shifts of the resonance frequencies. We present the simulation results of the soliton control by swiping the frequency of the counter-propagating wave in the forward and backward directions and with the soliton-blockade effect either present or not. We highlight the non-reciprocity of the forward and backward scans. Furthermore, we report the soliton crystals and breathers existing in the vicinity of the blockade interval.


Introduction
Frequency comb generation in the high-Q ring microresonators and the associated dissipative Kerr solitons have reached unprecedented heights of practical relevance [1]. The complexity of the microresonator soliton properties is hard to exaggerate. A challenge has emerged after a series of experiments with the bidirectionally pumped microresonators, where combs and solitons [2][3][4] and symmetry breaking [5,6] have been observed in counterpropagating waves. Studies into the gyroscope applications of these devices are also becoming increasingly important [7][8][9][10][11].
The problem of the mathematical modeling of the multimode regimes of operation of the bidirectional resonators has been addressed and to a large degree resolved only recently by demonstrating the equivalence of the coupled-mode model derived from the ab initio Maxwell system to the envelope equations where the nonlinear interaction (cross-phase modulation, XPM) between the counter-propagating waves is reduced to the power-dependent shifts of the resonance frequencies and the fast oscillations with the repetition rate frequencies being eliminated [12][13][14]. Below, we recapture the main steps of the model derivation and highlight the transition from the four-envelope model accounting for the repetition rates to the two-envelope formulation that eliminates them.
The deeper insight into the independence of the XPM nonlinearity from the phase relations between the interacting modes [12][13][14] has recently led us to the prediction of the soliton-blockade effect [15]. The blockade is achieved when the tuning frequency of one of the counter-rotating fields first disrupts and then restores the soliton transmission in the other field. The refractive index change of the soliton field happens via the XPM effect and is particularly strong in microresonators due to their high finesses providing 3-6 orders of magnitude boost to the circulating powers relative to the input one [16]. This creates a variety of opportunities for efficient nonlinear control of optical signals (see, e.g., [17,18] and references therein), and, in particular, the soliton blockade effect utilizes the high sensitivity of the nonlinear response of one of the fields towards the changes of the driving frequency of the counter-propagating one.
In Section 3, we look into how the power and frequency of the soliton and control pumps can be used to expand the soliton-blockade range. We also demonstrate there that the soliton-crystals, well known in the unidirectionally pumped resonators [19,20], emerge from the edges of the blockade interval. We performed both forward and backwards adiabatic scans of the control field frequency and demonstrate non-reciprocity of the solitonblockade effect (see Section 4). In Section 5, we present equations for perturbations around the solitons, which bare the features specific to the integral nature of the XPM terms and demonstrate the breather states close to the soliton-blockade interval (see, e.g., [21][22][23] for the breather studies in the unidirectional resonators).

Equations, Field Envelopes and Modes
We proceed by postulating that the real electric field of a given transverse mode family in a multi-mode ring microresonator can be expressed as [13] b F(r, z)E (t, θ). (1) The whole of the latter expression has units of volts per meter. θ ∈ (0, 2π] is the angular coordinate along the resonator circumference. t is physical time. b is the normalization constant that re-scales E to the desirable units. E 2 has units of power in what follows. F(r, z) is the transverse mode profile normalized so that its maximum equals one. The spectrum of E (t, θ) is assumed to be sufficiently narrow to allow disregarding the dispersive changes of F(r, z). Separation of the spatial variables in Equation (1) is an approximation that generally works well in a typical microresonator.
E can then be sought as a superposition of the two counter-propagating waves, Here, M is the mode number with the resonance frequency ω 0 and ω + is the frequency of the laser pumping the plus-wave, so that defines the respective detuning. Q ± are the envelope functions, which can be expressed via their mode expansions as Q ± µ are the time dependent mode amplitudes. µ = −N/2 + 1, . . . , 0, . . . , N/2 are the relative mode numbers and ω µ = ω 0 + D 1 µ + 1 2 D 2 µ 2 are the corresponding resonance frequencies. D 1 /2π is the resonator repetition rate and D 2 /2π is dispersion.
The structure of the mode expansion in Equation (4) assumes that the effects of D 1 , D 2 , linewidth and nonlinearity will all be embraced by the equations derived for Q ± µ after Equations (4) are substituted to the Maxwell equations. Using ω + in the exponents in front of the mode expansions for the clockwise and counter-clockwise fields in Equation (4) implies that detuning between the frequencies of the lasers pumping the counter-rotating fields will resurface in the equations for Q ± µ (see Figure 1 for a schematic illustration of the pump arrangements). The presence of the backscattering effects induces coupling between the counterpropagating waves. The authors of [13,14] formally traced what could be expected naturally, namely the Q + µ e iµθ wave couples to the Q − µ e iµθ one. Comparing this coupling structure to the definition of the Q − envelope function in Equation (4) suggests the need to define the envelope functions for the reflected waves, Q (r) While the mode amplitudes entering the Q ± and Q (r) ± envelopes are the same, the signs of the exponential terms are arranged differently (cf. Equations (4) and (6)).
The envelope equations that follow are (see [13] for the first-principle derivation using the present scaling and notations). There are several parameters first used in Equation (7) and hence requiring definitions [13]. κ is the loaded linewidth parameter. γ is the nonlinear coefficients measured, such that the units of γ/2π are Hz/W and of γ|Q ± | 2 /2π are Hz. H ± are the pump parameters, H 2 ± = ηF W ± /π. Here, η < 1 is the coupling coefficient, F = D 1 /κ is the resonator finesse and W ± are the on-chip powers of the two lasers. R is the backscattering coefficient, which is taken equal for all the modes (see [13] for a more general treatment).

Redefining the Envelope Functions and Eliminating the Repetition-Rate Terms
Equations (7) contain the D 1 terms with the opposite signs, therefore the transformation into the rotating reference frame could remove the repetition-rate term, for example, from the plus equations, but then the repetition rate would simply double in the minus equations. At the same time, D 1 is the strongly dominant frequency scale in the problem. Typically, it would be 10 GHz and up to 1 THz, while all the other terms, if taken for the respective range of resonators, would vary from 100 kHz to 100 MHz. This indicates that, if the multiples of D 1 are traced in the nonlinear terms, then these terms would be oscillating with the fastest frequency in the model and could be disregarded [12][13][14].
To reveal the frequency scales associated with D 1 , we replace the Q ± µ mode amplitudes in Equations (4) with so that the electric field in Equation (2) is replaced with Two new sets of the envelope functions play a pivotal role in the theory of the microresonators with the counter-propagating waves [13]. Taking the discrete Fourier transforms of ψ ± allows reconstructing E (cf. Equations (2), (8) and (9)). Unlike the equations for Q ± , the ones for ψ ± do not contain the D 1 terms and replace the system for four envelopes, Q ± , Q (r) ± , with the one for two, The above equations could be supplemented with equations for ψ (r) + , ψ − , but this time they are left as an independent pair. With the fabrication techniques rapidly improving, the surface roughness losses, as well as the associated backscattering, can be reduced to such levels that the frequency scale of R drops below the dispersion and nonlinearity induced effects so that we set R = 0 in what follows. Including R = 0 represents a separate problem [14].
The notable consequence of the elimination of the D 1 oscillations from the dynamics is that the nonlinear coupling, i.e., the cross-phase modulation (XPM) between the counterpropagating waves, has lost its sensitivity to the phases of the individual modes [12][13][14] (see Equation (11c)). Now, the XPM is simply expressed by the nonlinear shifts of the detuning, i.e., by the 2g + and 2g − terms.
In what follows, we remove the exponential term from Equation (11b) by applying an − e iεt , omit the '(r)' superscript, and deal with where

Single-Mode, Single-Soliton and Soliton-Crystal States and Their Role in the Soliton-Blockade Effect
The single-mode, µ = 0, regime of the resonator operation in the plus and minus fields is called here the cw-cw state (where cw stands for continuous wave). The respective modal amplitudes can be expressed as where g ± solve the coupled real algebraic equations Apart from the cw-cw states, there are also the soliton-cw, cw-soliton and solitonsoliton states, where the first and second words in the solution classification characterize the field in the ψ + and ψ − components, respectively [15]. Figure 2a shows how the cw-cw solution varies with δ 0 for a set of the fixed values of the frequency offset parameter, ε = ω + − ω − . Large |ε| separate and make quasiindependent the resonance structures in the two components along the δ 0 axis, i.e., if the g + vs. δ plot has the nonlinearity tilted resonance originating at δ 0 = 0, then g − has the similar resonance starting at δ 0 ≈ −ε. Varying ε from the relatively large negative to the large positive values drags the g − resonance across the effectively immobile g + one (see Figure 2a). Two resonances overlap and interact strongly for |ε|/κ 1. For ε = 0 and H + = H − , the model becomes symmetric and, therefore, its symmetric solution (see the black lines in the ε = 0 panels of Figure 2a) goes through the symmetry breaking bifurcation [5,6,15].
The soliton-cw states are solutions of the coupled differential-algebraic system If |ε|/κ is large, then Equation (16a) becomes quasi-independent from Equation (16b), implying that the system operates in the limit where it is approximately reduced to the Lugiato-Lefever model having the bright soliton solutions around the bistability interval for δ > 0 [1]. However, when |ε| becomes the order of the linewidth, the two equations start interacting strongly, and the low power branch of the cw-state in the plus-field becomes distorted by the growth of the second resonance peak, see the ε < 0 cases in the top row in Figure 2. This interaction also creates the narrow parameter range, where the low power branch of the plus field ceases to exist (see the yellow interval in Figure 2b). Around and substantially beyond this range, the soliton pulses in the plus-field are not able to exist simply because the system does not offer a low-amplitude background state for the bright soliton to nest on (see the grey stripe in Figure 2b and the respective no-soliton range in Figure 2c). This is the essence of the soliton-blockade effect, when the solitons can be switched on and off by tuning the frequency of the counter-propagating field, so that ε is tuned within the soliton-forbidden interval [15]. Figure 3 shows how the soliton-blockade domains are shaped in the (δ 0 , ε) and (H − , ε) parameter spaces, respectively. In particular, the blockade regime is more easily accessible if the minus , i.e., the cw component is pumped harder (see Figure 3b). Figure 4 shows how the XPM coefficients for the plus, g + , and minus, g − , components of the soliton-cw state (see Equation (16c)) vary with ε for the value of δ slightly different from the one in Figure 2c. Figure 4a plots the soliton g + for the one-, four-, and six-soliton states, i.e., for the soliton crystals. The respective g − values corresponding to the cwcomponent are shown in Figure 4b, and the soliton crystal profiles are in Figure 5. While the more detailed understanding of the bifurcations of the soliton crystals [19,20] in the bi-directionally pumped geometry goes beyond our present scope.

Controlling Solitons by Tuning the Cw-Component: Direct and Reverse Scans
We now present the results of the dynamic simulations of Equation (11) by initializing them with the cw-cw states and by scanning the frequency of the control (minus) field adiabatically. The adiabatic scans of ε were implemented in both forward (ε goes from negative to positive) and reverse directions. The outcomes of the two scans are generally different due to the complex and non-symmetric structure of the bistable resonances (see the cw-cw vs. ε plot in Figure 2b).
First, we choose H − /H + = 1, and continue to keep δ 0 on the positive side, so that the plus solitons would always exist in the absence of the interaction with the minus field. The outcomes of the forward and backward ε-scans are shown in the left and right columns of Figure 6, respectively. For the forward scan, the plus solitons keep being generated before ε/2π −4 MHz, then disappear in the soliton-blockade region and then reappear for ε/2π 7 MHz. Throughout the blockade interval, the intraresonator field most typically converges to the turbulent state known for the near-bistability operation in the uni-directional setting. For ε/2π 2 MHz, one can also see that the minus-field starts generating solitons.
In the forward scan, the soliton-blockade interval for the plus field is wider than predicted by the time-independent methods in the previous section. This is because the eight solitons formed in the minus field provide the value of g − sufficient to keep δ 0 − 2g − outside the plus-soliton existence interval. Simultaneously, the value of δ 0 + ε − 2g + is such that the minus field can sustain the solitons. The variations of the net power in either of the fields can lead to different scenarios as ε is increased. The one that is realized in the dataset shown in Figure 6 is that the downwards fluctuation of g + takes the minus field outside the soliton regime and brings it to the low power cw-state. Simultaneously, this transition restores the soliton regime in the plus field itself. The exit from the blockade regime is also sensitive to the chosen value of H + /H − and the scan parameters. Generally, nonlinear dynamics under conditions when a parameter is swept across the bifurcation points is an area of research attracting attention in optics and beyond (see, e.g., [24,25] and references therein). The reverse, i.e., ε-positive to ε-negative, scan shows the much narrower range of the soliton-blockade in the plus field, and no solitons in the minus field (see the right column in Figure 6). This is because the intra-resonator field in the minus component picks and follows the stable low power cw-state until ε comes to near zero. Therefore, in the reverse scan, the dynamically seen soliton-blockade range is practically the same as the one predicted by the time-independent analysis in the previous section.
We now choose H − /H + = 1/2, which makes the soliton-blockade impossible (see Figure 3b). The outcomes of the respective forward scan are shown in Figure 7. Following our expectations, we see only the stable solitons in the plus field. In the minus field, the periodic pattern is induced in the narrow range of δs where the low amplitude cw state is unstable. The uninterrupted solitons seeing in the plus field during the scan are correlated with the continuous existence of the low amplitude cw-state (see the red line in Figure 7). The narrow and low red bistability peak disturbs the pulse but does not switch the plus field into a persistent non-soliton regime (cf. Figure 6). The data in Figures 6 and 7 were generated while adiabatically tuning δ − /2π between −10 and 30 MHz over the 1 ms time interval (the round-trip time ≈ 60 fs, D 1 /2π = 15 GHz).  Figure 3b): (a,b) the plus and minus fields during the scan; and (c) g + for the cw-cw state over the same interval of ε.

Soliton Stability and Breathers
The problem of determining soliton stability is worthy of some attention, in particular, because the g ± terms in the governing equations represent a certain challenge in this regard. The soliton stability is determined by computing the eigenvalue of the Jacobian matrix. In order to work out the latter, we seek solutions of Equation (12) as ψ + = A(θ) + a(θ, t), ψ − = B(θ) + b(θ, t), where a, b are the small perturbations on top of the soliton or any other spatially inhomogeneous state, ψ + = A(θ), ψ − = B(θ).
Truncating the second and higher orders of a and b, we are left with If the integrals are replaced with the sums, then Equation (17) becomes ∂ t ε = −i M ε, where ε = [a, a * , b, b * ] T . The growth rate of the perturbations is introduced by setting ε ∼ e λt . The instability occurs if the real part of λ is positive. We computed the soliton spectrum around the blockade interval and found the oscillatory instabilities that stabilize away from the interval (see Figure 8a). When the instabilities are present, they typically lead to the formation of the soliton breather states with different periods (see Figure 8b,c).

Summary
We further investigated the soliton-blockade effect in the bidirectional microring resonators. In particular, we presented the simulation results on the soliton control by swiping the frequency of the counter-propagating wave in the forward and backward directions and found that the latter provides a much narrower blockade interval and a wider range of the soliton existence. Furthermore, we demonstrated that the blockade interval can be expanded by increasing the power of the control field and that the soliton crystals and breathers exist on the both sides of it.