Bound States of Light Bullets in Passively-Mode-Locked Semiconductor Lasers
BBound States of Light Bullets in Passively-Mode-Locked SemiconductorLasers
Fabian Dohmen, Julien Javaloyes, and Svetlana V. Gurevich
3, 4, 2, a) Institute for Theoretical Physics, University of Münster, Wilhelm-Klemm-Str. 9, D-48149 Münster,Germany Departament de Física, Universitat de les Illes Balears, & Institute of Applied Computing and Community Code (IAC-3),C/ Valldemossa km 7.5, 07122 Mallorca, Spain Institute for Theoretical Physics, University of Münster, Wilhelm-Klemm-Str. 9, D-48149 Münster,Germany Center for Nonlinear Science (CeNoS), University of Münster, Corrensstrasse 2, D-48149 Münster,Germany
In this paper, we analyze the dynamics and formation mechanisms of bound states (BSs) of light bullets in the output ofa laser coupled to a distant saturable absorber. First we approximate the full three-dimensional set of Haus master equa-tions by a reduced equation governing the dynamics of the transverse profile. This effective theory allows us to performa detailed multiparameter bifurcation study and to identify the different mechanisms of instability of BSs. In addition,our analysis reveals a non-intuitive dependence of the stability region as a function of the linewidth enhancement factorsand the field diffusion. Our results are confirmed by direct numerical simulations of the full system.
Three-dimensional spatio-temporal localized states, so-called light bullets (LBs) were recently predicted theoreti-cally in the output of a laser coupled to a distant saturableabsorber. In this paper, we analyze the stability and therange of existence of bound states consisting of LBs. Inorder to reduce the complexity of the analysis, we first ap-proximate the three-dimensional model by a reduced ef-fective equation governing the dynamics of the transverseprofile. This effective theory provides an intuitive pictureof the BS formation mechanism. Moreover, it allows us toperform a detailed multiparameter bifurcation study andto identify the different mechanisms of instabilities.
I. INTRODUCTION
Localized structures (LSs) appear in a variety of fields rang-ing from chemical, biological and optical systems to plantecology and social sciences . They can interact via expo-nentially decaying tails and form so-called bound states (BSs)or molecules. BSs can emerge e.g., due to the presence of os-cillating tails of the individual interacting LSs, resulting in theformation of structures with stable equilibrium distances andphase differences . If, however, the tails decay monoton-ically or a non-local repulsive interaction between individualLSs is present, they tend to distribute equidistantly in space ortime, leading to periodic pulse trains exhibiting large distancesbetween the consequent LSs . However, even if LSs in anindividual system exhibit strong repulsion, the formation ofBSs can be achieved by arranging several systems in an arraywith nearest-neighbor coupling . Recently it was shown thatBSs can also be created in systems with a pointwise nonlo-cality. There, the resulting molecules are composed by LSswhich are globally bounded but locally independent . a) Electronic mail: [email protected]
In optics, LSs of light have been intensively studied theoret-ically and observed experimentally in both spatial and tempo-ral domains . In particular, temporal LSs were observedin a semiconductor passively mode-locked laser . Passivemode locking (PML) is a well known method for achievingshort optical pulses . It is achieved by combining two el-ements inside an optical cavity, a laser amplifier providinggain and a nonlinear loss element, usually a saturable ab-sorber (SA). For proper parameters, this combination leads tothe emission of temporal pulses much shorter than the cavityround-trip. In it was shown that in the so-called long delaylimit, where the round-trip time is much longer than the semi-conductor gain recovery time, the PML pulses become indi-vidually addressable temporal LSs. Interestingly, this tempo-ral localization regime was recently found to be a basis forthe generation of the long sought three-dimensional spatio-temporal LSs, so-called light bullets (LBs) . These pulses oflight are simultaneously confined in the transverse and alongthe propagation directions. Note that in many cases LBs arefound to be unstable and collapse in three dimensions , orthey were observed in models neglecting the semiconductordynamics . As such, the LBs reported in are stiff mul-tiple timescale objects since the optical pulse is followed bya material trail that can differ in extension by three orders ofmagnitudes. This stiffness, stemming from the temporal do-main, is aggravated in the presence of the transverse dimen-sions making a bifurcation analysis of two and three dimen-sional LBs a challenging problem. The properties of a sin-gle LB were investigated in . There, the dynamics of three-dimensional LBs was approximated by the superposition ofa slowly evolving transverse profile and a short pulse propa-gating inside the cavity. This allowed to reduce the dynamicsto a partial differential equation governing the transverse pro-file. A detailed multi-parameter bifurcation study identifiedseveral mechanisms of instability where the LBs experienceeither homogeneous oscillation or symmetry breaking lateralwaves radiation.In this paper we study the dynamics and formation mecha-nisms of BSs consisting of LBs. We perform our analysis in a r X i v : . [ n li n . PS ] J a n two steps: Starting with the three-dimensional generic Hausmodel, we approximate its solution by the product of a slowlyevolving transverse profile and a short pulse propagating in-side the cavity. The obtained reduced two-dimensional model,governing the dynamics of the transverse profile is then an-alyzed employing the continuation and bifurcation packagepde2path . Starting with the one-dimensional case thatcorresponds to the two-dimensional BS of the full system, weshow that BSs corresponding to different phase differencescan be observed and the stability of these solutions is studied.Our analysis reveals a non-intuitive dependence of the stabil-ity region as a function of the linewidth enhancement factorsand the field diffusion. In the second stage, our predictionsare confirmed by two-dimensional direct numerical simula-tions of the full system. II. MODEL
We consider a passively mode-locked laser where a broadarea gain chip is coupled to a distant saturable absorber (SA)with telescopic optics in self-imaging conditions . We as-sume that the diffraction in the system, resulting from thepropagation within the active sections, is sufficiently small sothat the paraxial approximation can be used . Further, weconsider the uniform field limit i.e., the limit of moderate gain( G ) and saturable absorption ( Q ). In this case, the existenceand the dynamical properties of LBs in passively mode-lockedVCSELs can be theoretically described using the genericHaus partial differential equation (PDE) for the field profile E ( r ⊥ , z , σ ) over the slow time scale σ∂ σ E = (cid:26) √ κ (cid:20) + − i α G ( r ⊥ , z , σ ) − − i β Q ( r ⊥ , z , σ ) (cid:21) − + γ ∂ z + ( d + i ) ∆ ⊥ (cid:27) E ( r ⊥ , z , σ ) , (1)whereas the carrier dynamics is given by ∂ z G = Γ G − G (cid:16) Γ + | E | (cid:17) + D g ∆ ⊥ G , (2) ∂ z Q = Q − Q (cid:16) + s | E | (cid:17) + D q ∆ ⊥ Q . (3)Here, r ⊥ = ( x , y ) are transverse space variables, i.e., ∆ ⊥ = ∂ x + ∂ y is the transverse Laplacian and the longitudinal vari-able ( z ) can be identified as a fast time variable and representsthe evolution of the field within the round-trip. Further, κ isthe fraction of the power remaining in the cavity after eachround-trip, γ is the bandwidth of the spectral filter represent-ing, e.g., the resonance of a VCSEL and α and β denotethe linewidth enhancement factors of the gain and absorbersections, respectively. The parameter d represents the smallamount of field diffusion incurred for instance by the depen-dence of the reflectivity of the VCSEL distributed Bragg re-flectors upon the angle of incidence. In Eqs. (2-3) G de-scribes the pumping rate, Γ = τ − g the gain recovery rate, Q is the value of the unsaturated losses, s the ratio of the satu-ration energy of the gain and of the SA sections and D g , q thescaled diffusion coefficients. However, in it was shown that carrier diffusion plays almost no role in the LBs dynamics, sothat we set both diffusion coefficients to zero.We define G th = √ κ − + Q as the threshold gain valueabove which the off solution ( E , G , Q ) = ( , G , Q ) becomesunstable. Then for the gain values below the threshold G th and for proper system parameters Eqs. (1-3) possess a stablethree-dimensional LB solution . Note that in general, thenon-instantaneous response of the active medium implies alack of parity along ( z ) for the LBs .Using the fact that the LBs are in general composed ofvariables evolving over widely different timescales, an ap-proximate model governing the shape of the transverse pro-file can be derived . In particular, if one considers E ( r ⊥ , z , σ ) = A ( r ⊥ , σ ) p ( z ) , where p ( z ) is a short normal-ized temporal pulse of length τ p that represents the temporallocalized state upon which the LB is built and A ( r ⊥ , σ ) is aslowly evolving amplitude, one can separate the temporal evo-lution into the fast and slow parts corresponding to the pulseemission and the subsequent gain recovery. This leads to thefollowing equation governing the dynamics of A ( r ⊥ , σ ) : ∂ t A = ( d + i ) ∆ A + f (cid:16) | A | (cid:17) A . (4)Here, the nonlinear function f = f ( | A | ) reads f ( I ) = ( − i α ) g ( + q ) h ( I ) − ( − i β ) qh ( sI ) − , (5)where h ( I ) = (cid:0) − e − I (cid:1) / I , I = | A | .Similar to the derivation provided in , the spatial and tem-poral coordinates in Eqs. (4, 5) are scaled as t = (cid:0) − √ κ (cid:1) σ and ( x , y ) : = (cid:112) − √ κ ( x , y ) . Further, the gain normalized tothreshold g and the normalized absorption q are defined as g = G / G th and q = Q / (cid:16) √ κ − (cid:17) , respectively.Note that Eq. (4) governing the dynamics of the transverseprofile is known in the context of static transverse autosolitonsin bistable interferometers . There, the nonlinear function h ( I ) would correspond to a static saturated nonlinearity, i.e. f ( I ) = h ( I ) → / ( + I ) . By analogy with the static case, andalthough our nonlinear function is different, we shall call ourmodel for the transverse profile a Rosanov equation. In thiscase, for small values of the diffusion coefficient d , an inter-action law between two weakly overlapping autosolitons wasderived in . In particular, it was shown that different boundstates corresponding to phase differences ∆ ϕ of 0, π and ± π / ∆ ϕ = { , π } are found to be stable if d exceeds certain threshold, whereasfor d = ∆ ϕ = π corresponding to the minimaldistance between both autosolitons is stable. In addition, allBS solutions with ∆ ϕ = ± π / ∂ z G = ∂ z Q =
0. However, the reactiontime of the gain is known to profoundly affect the stability ofspatio-temporal structures , that is, its adiabatic eliminationapproximation along the propagation direction could be incor-rect for a semiconductor material. Hence, in what follows wekeep the form (5) of the nonlinear function f ( I ) and performa bifurcation analysis of BS solutions of Eq. (4), comparingthe results with whose of the Haus PDE (1)-(3). III. RESULTS
Bifurcation analysis of the Rosanov equation
We start ouranalysis with the one-dimensional case. In order to track thelocalized BS solutions of Eqs. (4–5) in parameter space, wemake use of pde2path , a numerical pseudo-arc-length bi-furcation and continuation package for systems of partial dif-ferential equations. Taking into account the traslational andthe phase-shift symmetries of the system in question, station-ary LSs of Eqs. (4–5) can be found in the form A ( x , t ) = u ( x − vt ) e − i ω t , (6)where u is a complex amplitude with the field intensity lo-calized around some point in space, ω represents the carrierfrequency of the solution and v is a drift velocity which is ingeneral zero for a BS stationary in time. Substituting Eq. (6)into Eqs. (4, 5) we obtain the following equation for unknowns u , ω and v ( d + i ) ∂ x u + v ∂ x u + i ω u + f (cid:16) | u | (cid:17) u = . (7)In the case of Eq. (4) (or (7)), the primary continuation pa-rameter can be, e.g., the gain parameter g , the field diffusion d or the linewidth enhancement factors α , β . However, thespectral parameter ω and the drift velocity v are two addi-tional free parameters that have to be automatically adaptedduring the continuation. In order to determine them, we im-pose additional integral auxiliary conditions (cid:90) ∂ u old ∂ r udx = , where u old denotes the solution obtained in the previous con-tinuation step and r represents the spatial coordinate or thephase of the solution. In addition we imposed periodic bound-ary conditions on the domain of the length L consisting of N x grid points, i.e. u ( ) = u ( L ) . Now one can start a continuationalgorithm at, e.g., a numerically given BS solution, continueit in parameter space, and obtain a solution branch. The re-sult of such a continuation for the in-phase BS solution (i.e.the phase difference ∆ ϕ =
0) is depicted in Fig. 1, where in(a) the integrated intensity P = (cid:82) | u | as a function of the gainparameter g is shown. One can see that on the high powerbranch, the BS is unstable w.r.t. Andronov-Hopf (AH) bi-furcation H for large values of g . However, at the next AHpoint H the BS gains the stability and remains stable for acertain interval of g (black thick line). An example of the sta-ble intensity profile (black) is shown in the panel (d) togetherwith ℜ ( u ) (blue dotted) and ℑ ( u ) (orange dashed) of the field.However, the BS looses its stability at further AH point H close to the fold F . At the lower power branch, the unstableBS solution undergoes the last AH bifurcation at H and re-mains unstable, decreasing its energy with increasing g . (cf.Fig. 1 (b)). However, at some fixed g value BP (green circle) itexperiences another bifurcation. There, a slightly snaking un-stable branch of a moving BS solution emerges (see Fig. 1 (c)for the exemplary profile). This branch itself has the multi-ple AH and fold bifurcation points. Since there are no stable x -0.80.8a)BP b)c)d)H H d)b)H FH c) FIG. 1. (a) A branch of an one-dimensional BS solution of Eq.(4),calculated for ∆ ϕ = g . The in-tegrated intensity P = (cid:82) I is shown. AH bifurcation points H − H are marked by red diamonds, whereas the fold ( F ) is marked by ablue cross. A branch of unstable moving BS bifurcates from thebranching point (BP, green circle). A BS solution is stable be-tween AH points H and H (thick line). (b)-(d) Three intensityprofiles (black solid lines) corresponding to b) steady unstable BSat g = . g = .
663 and d) steady sta-ble BS at g = . ℜ ( u ) and ℑ ( u ) , respectively. Other parameters are( γ , κ , α , β , q , s , d , L , N x ) = ( , . , . , . , . , , . , , ) . g P -0.40.4-50 0 50 x -0.80.8a) b)c)H H H H c)b)H FH H H FIG. 2. (a) A branch of an one-dimensional BS solution of Eq.(4),calculated for ∆ ϕ = π as a function of the gain parameter g . Theintegrated intensity P = (cid:82) I is shown. AH bifurcation points H − H are marked by red diamonds, whereas the fold ( F ) is marked by ablue cross. A BS solution is stable between AH points H and H aswell as between H and H (thick black line). (b), (c) Two intensityprofiles (black solid lines) corresponding to b) steady unstable BSand (c) steady stable BS at g = . ℜ ( u ) and ℑ ( u ) , respectively. Other pa-rameters are as in Fig. 1. solutions along the branch they are not shown here. In thecase of the BS corresponding to ∆ ϕ = π , the solution branchis shown in Fig. 2 (a) for the same diffusion value d = . ∆ ϕ =
0, the BS solution is un-stable on the high power branch at large g w.r.t. three AHbifurcations H − H and gains the stability at AH point H .However, one can see that the region of stability is smallerthan in the case of ∆ ϕ = H − H and H − H . After thefold F and the following AH point H , the BS solution re-mains unstable at the low power branch. Note that the branchwith unstable moving solutions does not appear for ∆ ϕ = π .Figure 2 (b), (c) represent two exemplary intensity profilesat g = . ℜ ( u ) and ℑ ( u ) of the field. Finally,Fig. 3 (b) shows the part of the branch of BS solutions calcu-lated for ∆ ϕ = π ∆ ϕ = { , π } , this branch does not contain sta-ble time-independent BSs. Instead, the BSs are moving (cf.3 (a)) with a fixed velocity v depending on the gain value (seeFig. 3 (c)). FIG. 3. (a) Space-time representation of an one-dimensional movingBS solution of Eq.(4) calculated numerically for ∆ ϕ = π g = . ( g , P ) plane. The whole branch is unstable. (c) The drift velocity ofa moving BS as a function of g . Other parameters are as in Fig. 1. The linear stability of a particular BS solution along thebranch can be obtained directly during the continuation sothat one has access to the critical eigenfunctions of the systemthat inform on the particular shape of the wave form. It re-veals that the eigenfunctions corresponding to the AH pointsfound for both phase differences { , π } have different sym-metry properties w.r.t. the center of each individual LB aswell as w.r.t. the center of the BS. In particular, the AH point H ( H ) corresponding to ∆ ϕ = H , H for ∆ ϕ = π is depicted in Fig. 4 (a)-(b), whereas panels (c)-(d)show a space-time representation of the intensity field evolu-tion obtained by direct numerical simulations of Eqs. (4-5) fortwo different values of g close to the corresponding AH bifur- FIG. 4. (a,b) Real (blue dotted) and imaginary (orange dashed) partsof critial eigenfunctions associated with the AH points ( H , H ) to-gether with normalized intensity profile (yellow solid) caclulated for ∆ ϕ = π at g = .
631 and g = . g cor-responding to (a,b). Other parameters are as in Fig. 1. cation points H , H . As it can be seen in Fig. 4 (a), both realand imaginary parts of the critical eigenfunction are asym-metric towards the center of the BSe w.r.t. the individual LBcenters. Hence, an asymmetrical breathing in space wouldbe expected which is verified by direct numerical simulationsrepresented in the panel (c). In Fig. 4 (b), on the contrary,the real part of the eigenfunction is asymmetric towards thecenter and is nearly symmetrical in the real part and asym-metrical in the imaginary part w.r.t. the peak maxima. Basedon the value of the real and imaginary contributions one wouldassume a breathing oscillation with an asymmetry in the pulseenvelope, cf. Fig. 4 (d). There, a clear pulsation in time isobservable as well as a slight elongation of the pulse towardsthe center. The oscillatory dynamics of the BS can however FIG. 5. Complex spatio-temporal oscillatory dynamics of a BS so-lution observed in the high intensity branch at (a) g = .
668 and (b) g = .
664 between AH bifurcation points H and H for ∆ ϕ = π .Other parameters are as in Fig. 1. become more complex for increasing values of g . Figure 5shows two examples of the temporal evolution of the BS ob-served between H and H bifurcation points. Here, severalfrequencies contribute to the dynamics leading to a complexspatio-temporal behavior. The numerical continuation analy-sis reveals the hysteresis behavior for increasing (decreasing)values of g indicating the subcritical character of the AH bi-furcation H .Besides the gain g , the linewidth enhancement factors of thegain and the absorber sections as well as the field diffusionare important control parameters that determine the stabilityrange of the BS. In order to study the influence of α , β and d on the stability range of the single BS solution we performa two-parameter continuation of the fold and AH bifurcationpoints for the both phase differences ∆ ϕ = , π . Figure 6 (a) H H H H H FIG. 6. Evolution of the fold (black squares) and AH points H − H in the ( g , d ) plane for the phase differences (a) ∆ ϕ = ∆ ϕ = π (blue triangles) (cf. Figs. 1, 2). No stable BSsolutions exist below the threshold d th for the field diffusion with theminimal value of (a) d th (cid:39) .
002 and (b) d th (cid:39) . shows the stability diagram in ( g , d ) plane for ∆ ϕ =
0. Here,black squares indicate the fold F threshold, whereas red tri-angles stand for the boundary of the AH bifurcations H and H . The BS solution is stable in the region bounded by the F , H and H lines. One can see that the AH point H moves inthe direction of the fold F for increasing values of d , so thatabove a certain threshold, the stability region is bounded by F and H lines only. Notice that the stability region is lim-ited from below in d , that is, there is no stable BS solutionsexist below a certain field diffusion threshold d th dependingon the gain value g . The stability diagram for the phase dif-ference ∆ ϕ = π is depicted in Fig. 6 (b). Here, the stabilityregion for the BS solution is limited by the line given by thefold position F (black squares) and the boundaries of the AHbifurcations H , H and H , H (blue triangles). For increas-ing values of d , the AH point H moves in the direction of thefold F , whereas H and H move toward each other until theycollide at some critical value of d . That is, for large valuesof the field diffusion the BS is stable between F and H linesand only one stability region exists. Note that as in the caseof ∆ ϕ =
0, there exist a threshold value of the field diffusionthat bounds the stability region from below. This thresholdvalue also depends on the gain value and is in general higherthen in ∆ ϕ = , where theantiphase BS can be stable for vanishing values of the fielddiffusion. H F H H H H F FIG. 7. One-dimensional stability diagrams showing the evolutionof the threshold of the fold F (black squares) and AH bifurcationpoints (a) H − H and (b) H − H in ( g , α ) plane for (a) ∆ ϕ = ∆ ϕ = π (blue triangles). A BS is stable (a) inthe area bounded by the lines F , H , H and (b) between the lines H , H and below the lines H and H . In the influence of both linewidth enhancement factorson the stability of a single LB was studied. In particular itwas shown that the range of the stability of the LB increasestoward higher α values and decreases for growing β . Fig-ure 7 represents the corresponding bifurcation diagram for theBS in the ( g , α ) plane for the fixed value of β and d for (a) ∆ ϕ = ∆ ϕ = π , respectively. One can see that in thefirst case the fold position F (black squares) remains almostconstant for increasing α and almost coincides with the AH H line (red triangles). Hence, the stability region is boundedby the AH lines H and H and is getting wider for increasing α values. For the phase difference ∆ ϕ = π the stability di-agram is more complicated since more AH points influencethe dynamics of the BS, see Fig. 7 (b). As in the case of ∆ ϕ =
0, the fold position remains almost constant for increas-ing α and almost coincides with the line of the AH point H .Hence, for small values of α , the BS is stable between the H and H lines. However for increasing α two other AH points H and H emerge. The corresponding lines are moving to-wards larger (smaller) gain values so that two stability regionsemerge. However, the H line collides with H at some α , sothat beyond this point, only one region of stability between H and H lines exist. H F F H H H H FIG. 8. One-dimensional stability diagrams showing the evolutionof the fold F (black squares) and AH bifurcation points (a) H − H and (b) H − H in ( g , β ) plane for (a) ∆ ϕ = ∆ ϕ = π (blue triangles). A BS is stable in the area bounded by thelines (a) F , H , H and (b) F , H and H − H . A similar stability diagram in ( g , β ) plane for the fixedvalue of α and d is shown in Fig. 8. The width of the stabilityregion is limited with either the fold F or the correspondingAH lines ( H and H − H , respectively) on the left and theAH line H ( H ) on the right. In both cases, the region ofstability shrinks with increasing β . Bound state dynamics in a 2D Rosanov equation
In thetwo-dimensional case of Eqs. (4)-(5) the dynamics of the BSsolution is similar to the one-dimensional case: A radially-symmetric stationary BS solution exists over the wide rangeof the gain g values and becomes unstable via several AH bi-furcations, which morphology depends on the phase differ-ence, the field diffusion and linewidth enhancement factors.Two examples of different oscillatory behavior calculated for ∆ ϕ = π for two different gain values is shown in Fig. 9, whereFig. 9 (a) represents the time evolution of a cross section of aBS of two LB that both oscillate in time keeping the distancebetween both LBs constant, whereas the panel (b) shows a BSwith an additional oscillation in the distance. FIG. 9. A cross-section of an oscillating BS of the two-dimensionalRosanov equation (4)-(5) calculated on the two-dimensional domain200 ×
200 with 512 ×
512 grid points for ∆ ϕ = π at (a) g = .
725 and(b) g = . Bound states of LB in the 2D Haus model
In order to val-idate our approach, now we compare the results of the bifur-cation analysis of the one-dimensional Rosanov model (4)-(5)with the numerical predictions of the two-dimensional Hausequations (1)-(3); The details of the numerical methods aregiven in . We study Eqs. (1)-(3) in the regimes of LSs, thatis, for values of the gain below threshold. In this regime, sim-ilar to the single LBs, BSs occur as multiple time scale ob-jects as shown in Fig. 10: Here, one notices that the opticalcomponent presented in Fig. 10 (a) is much shorter than thematerial ones in Fig. 10 (b),(c), where the gain G = G ( x , z ) and absorber Q = Q ( x , z ) fields are depicted. One can see thatthe BS presented in Fig. 10 consists in two LBs, which lon-gitudinal components correspond to two temporal LSs. In thetransverse plane the pulses have weakly overlapping oscilla-tory tails, thus forming a BS with a certain distance d andphase difference ∆ ϕ between the pulses. To analyze the dy-namical properties of these two-dimensional BS, a series ofdirect numerical simulations of Eqs. (1)-(3) were conducted.They indicate that for both phase differences stable BSs be-come unstable via AH bifurcations leading to a complex oscil-latory dynamics. The resulting dynamical regimes are similarto those obtained from the one-dimensional Rosanov Eqs. (4)-(5). One example is shown in Fig. 11, where a temporal evolu-tion of the BS’s cross-section is shown for two different valuesof the gain and ∆ ϕ = π . One can see in Fig. 11 (a) that the BS FIG. 10. An exemplary numerical BS solution of Eqs. (1)-(3) showing (a) the intensity profile I (b) gain G and (c)absorber Q fields of the stable BS. The spatial compo-nent is in phase corresponding to ∆ ϕ =
0. Other param-eters are ( γ , κ , α , β , G , Γ , Q , s , d , L x , L z , N x , N z )=(40, 0.81.5, 0.5, 0.736, 0.04, 0.3, 30, 0.03, 300, 3, 512, 512).FIG. 11. Spatio-temporal dynamics of a cross-section of a BS con-sisting of two LBs obtained from a direct numerical simulations of atwo-dimensional Haus Eqs. (1)-(3) for ∆ ϕ = π and (a) g = .
738 and(b) g = . can, e.g., experience oscillation w.r.t. the center of each LBs,whereas the distance between them remains fixed in time (cf.Fig. 4 (d)) or can show complex quasi-periodic explosion-likeoscillations as in Fig. 11 (b) similar to the dynamics observedfor the one-dimensional Rosanov model in Fig. 5 (b). IV. CONCLUSION
In conclusion, in this paper we discussed the dynamics andformation mechanisms of bound states consisting of light bul-lets. We have shown that the dynamics of a three-dimensionalBS can be successfully approximated by a simplified modelgoverning the dynamics of the transverse profile of the BS.The bifurcation analysis of this effective Rosanov equationallowed us to obtain guidelines regarding the existence andstability of the BS. Starting with the case of one spatial di-mension, we have found that BSs corresponding to differ-ent phase differences can exist. While BSs correspondingto ∆ ϕ = { , π } can be stable in a range of system param-eters, the stationary BSs with ∆ ϕ = π where the dynamics of transverse autosolitons in bistableinterferometers corresponding to the static saturated nonlin-earity in Eqs. (4)-(5) was studied. However, in our case forboth ∆ ϕ = { , π } the existence of a threshold value of thefield diffusion that bounds the stability region from below wasdemonstrated, whereas in the antiphase BS can be stable forvanishing values of the field diffusion. This threshold valuealso depends on the gain value and is in general higher thenin ∆ ϕ = are expectedto influence the dynamics of BSs. The analysis and charac-terization of the full two-dimensional dynamics of the BSs ofthe Haus model can effectively be done employing the path-continuation methods. However, the multiscale nature of in-dividual LBs and the additional transverse spatial dynamicsmakes the implementation technically involved and will be atopic of further studies. ACKNOWLEDGMENTS
S.G. thanks PRIME programme of the German AcademicExchange Service (DAAD) with funds from the German Fed-eral Ministry of Education and Research (BMBF). J.J. ac-knowledge the financial support of the MINECO ProjectMOVELIGHT (PGC2018-099637-B-100 AEI/FEDER UE). N. Akhmediev and A. Ankiewicz.
Dissipative Solitons: From Optics toBiology and Medicine Series , volume 751 of
Lecture Notes in Physics .Springer Berlin Heidelberg, 2008. M. Tlidi, P. Mandel, and R. Lefever. Localized structures and localizedpatterns in optical bistability.
Phys. Rev. Lett. , 73:640–643, Aug 1994. W. J. Firth and A. J. Scroggie. Optical bullet holes: Robust controllablelocalized states of a nonlinear cavity.
Phys. Rev. Lett. , 76:1623–1626, Mar1996. V. K. Vanag and I. R. Epstein. Localized patterns in reaction-diffusionsystems.
Chaos: An Interdisciplinary Journal of Nonlinear Science ,17:037110, 2007. H. G. Purwins, H. U. Bödeker, and S. Amiranashvili. Dissipative solitons.
Adv. Phys. , 59:1193, 2010. A. W. Liehr.
Dissipative solitons in reaction-diffusion systems . Springer,2013. Ehud Meron.
Nonlinear physics of ecosystems . CRC Press, 2015. David J. B. Lloyd, Christian Gollwitzer, Ingo Rehberg, and ReinhardRichter. Homoclinic snaking near the surface instability of a polarisablefluid.
Journal of Fluid Mechanics , 783:283–305, 2015. David J.B. Lloyd and Hayley O’Farrell. On localised hotspots of an urbancrime model.
Physica D , 253:23 – 39, 2013. N. N. Rosanov and G. V. Khodova. Diffractive autosolitons in nonlinearinterferometers.
J. Opt. Soc. Am. B , 7:1057–1065, 1990. Boris A. Malomed. Bound solitons in the nonlinear schrödinger–ginzburg-landau equation.
Phys. Rev. A , 44:6954–6957, Nov 1991. N. N. Akhmediev, A. Ankiewicz, and J. M. Soto-Crespo. Multisoliton solu-tions of the complex ginzburg-landau equation.
Phys. Rev. Lett. , 79:4047–4051, Nov 1997. D. Y. Tang, W. S. Man, H. Y. Tam, and P. D. Drummond. Observation ofbound states of solitons in a passively mode-locked fiber laser.
Phys. Rev.A , 64:033814, Aug 2001. A. S. Moskalenko, A. W. Liehr, and H. G. Purwins. Rotational bifurca-tion of localized dissipative structures.
EPL (Europhysics Letters) , 63:361,2003. Dmitry Turaev, Andrei G. Vladimirov, and Sergey Zelik. Chaotic boundstate of localized structures in the complex ginzburg-landau equation.
Phys.Rev. E , 75:045601, Apr 2007. P. Grelu and J.M. Soto-Crespo. Temporal soliton "molecules" in mode-locked lasers: Collisions, pulsations, and vibrations. In
Dissipative Soli-tons: From Optics to Biology and Medicine , volume 751 of
Lecture Notesin Physics , pages 1–37. Springer Berlin Heidelberg, 2008. Hidetsugu Sakaguchi, Dmitry V. Skryabin, and Boris A. Malomed. Sta-tionary and oscillatory bound states of dissipative solitons created by third-order dispersion.
Opt. Lett. , 43(11):2688–2691, Jun 2018. C. Elphick, E. Meron, J. Rinzel, and E.A. Spiegel. Impulse patterning andrelaxation propagation in excitable media.
Journal of Theoretical Biology ,146(2):249 – 268, 1990. Michel Nizette, Dmitrii Rachinskii, Andrei Vladimirov, and Matthias Wol-frum. Pulse interaction via gain and loss dynamics in passive mode locking.
Physica D: Nonlinear Phenomena , 218(1):95 – 104, 2006. P. Camelin, J. Javaloyes, M. Marconi, and M. Giudici. Electrical addressingand temporal tweezing of localized pulses in passively-mode-locked semi-conductor lasers.
Phys. Rev. A , 94:063854, Dec 2016. D. Puzyrev, A. G. Vladimirov, A. Pimenov, S. V. Gurevich, and S. Yanchuk.Bound pulse trains in arrays of coupled spatially extended dynamical sys-tems.
Phys. Rev. Lett. , 119:163901, Oct 2017. J. Javaloyes, M. Marconi, and M. Giudici. Nonlocality induces chains ofnested dissipative solitons.
Phys. Rev. Lett. , 119:033904, Jul 2017. M. Brambilla, L. A. Lugiato, F. Prati, L. Spinelli, and W. J. Firth. Spatialsoliton pixels in semiconductor devices.
Phys. Rev. Lett. , 79:2042–2045,1997. S. Barland, J. R. Tredicce, M. Brambilla, L. A. Lugiato, S. Balle, M. Giu-dici, T. Maggipinto, L. Spinelli, G. Tissoni, T. Knödl, M. Miller, andR. Jäger. Cavity solitons as pixels in semiconductor microcavities.
Nature ,419(6908):699–702, Oct 2002. F. Leo, S. Coen, P. Kockaert, S.P. Gorza, P. Emplit, and M. Haelterman.Temporal cavity solitons in one-dimensional kerr media as bits in an all-optical buffer.
Nat Photon , 4(7):471–476, Jul 2010. T. Herr, V. Brasch, J. D. Jost, C. Y. Wang, N. M. Kondratiev, M. L. Gorodet-sky, and T. J. Kippenberg. Temporal solitons in optical microresonators.
Nature Photonics , 8(2):145–152, 2014. M. Marconi, J. Javaloyes, S. Balle, and M. Giudici. How lasing localizedstructures evolve out of passive mode locking.
Phys. Rev. Lett. , 112:223901,Jun 2014. H. A. Haus. Mode-locking of lasers.
IEEE J. Selected Topics QuantumElectron. , 6:1173–1185, 2000. J. Javaloyes. Cavity light bullets in passively mode-locked semiconductorlasers.
Phys. Rev. Lett. , 116:043901, Jan 2016. Y. Silberberg. Collapse of optical pulses.
Opt. Lett. , 15:1282–1284, 1990. N.A. Veretenov, A.G. Vladimirov, N.A. Kaliteevskii, N.N. Rozanov, S.V.Fedorov, and A.N. Shatsev. Conditions for the existence of laser bullets.
Optics and Spectroscopy , 89(3):380–383, 2000. M. Brambilla, T. Maggipinto, G. Patera, and L. Columbo. Cavity lightbullets: Three-dimensional localized structures in a nonlinear optical res-onator.
Phys. Rev. Lett. , 93:203901, 2004. S. V. Gurevich and J. Javaloyes. Spatial instabilities of light bullets inpassively-mode-locked lasers.
Phys. Rev. A , 96:023821, Aug 2017. T. Dohnal, J. Rademacher, H. Uecker, and D. Wetzel. pde2path-version 2.0:faster fem, multi-parameter continuation, nonlinear boundary conditions,and periodic domains-a short manual. arXiv:1409.3119 , 2014. H. Uecker, D. Wetzel, and J. Rademacher. pde2path - a matlab package forcontinuation and bifurcation in 2d elliptic systems.
Numerical Mathemat-ics: Theory, Methods and Applications , 7:58–106, 2014. P. Genevet, S. Barland, M. Giudici, and J. R. Tredicce. Cavity solitonlaser based on mutually coupled semiconductor microresonators.
Phys. Rev.Lett. , 101:123905, Sep 2008. M. Marconi, J. Javaloyes, S. Balle, and M. Giudici. Passive mode-lockingand tilted waves in broad-area vertical-cavity surface-emitting lasers.
Se- lected Topics in Quantum Electronics, IEEE Journal of , 21(1):85–93, Jan2015. C. Schelte, J. Javaloyes, and S. V. Gurevich. Dynamics of temporally lo-calized states in passively mode-locked semiconductor lasers.
Phys. Rev. A ,97:053820, May 2018. G. New. Pulse evolution in mode-locked quasi-continuous lasers.
QuantumElectronics, IEEE Journal of , 10(2):115 – 124, February 1974. N. N. Rosanov and G. V. Khodova. Autosolitons in nonlinear interferome-ters.
Opt. Spectrosc. , 65:449–450, 1988. A. G. Vladimirov, S. V. Fedorov, N. A. Kaliteevskii, G. V. Khodova, andN. N. Rosanov. Numerical investigation of laser localized structures.
Jour-nal of Optics B: Quantum and Semiclassical Optics , 1(1):101, 1999. A. G. Vladimirov, G. V. Khodova, and N. N. Rosanov. Stable bound statesof one-dimensional autosolitons in a bistable laser.
Phys. Rev. E , 63:056607,Apr 2001. L. Columbo, I. M. Perrini, T. Maggipinto, and M. Brambilla. 3d self-organized patterns in the field profile of a semiconductor resonator.