3 tional to the pulse energy The conservation of E is used to eliminate A from the remaining variational equations where the equations for β ′ and T ′ were used to obtain the single second order equation for T $^{′′}$. E ≡ A $^{2}$WT = 1 π ∫ ∫ | u | $^{2}$dxdt. (16) W ′ = λbW, (17) b ′ = λ W 4 − λb 2 − E 2 W $^{3}$T , (18) T ′′ = 1 T 3 − E 2 WT $_{2}$, (19) The systematic results obtained using the variational model Eq. (18) are shown in Fig. 1(c). The dotted area in Fig. 1(c) corresponds to the values of parameters E and B for which stable light bullet solutions exist which oscillate with the period of the waveguide modulation. The solutions of the variational equations were obtained using Newton relaxation method. The two boundaries of the stability region can be found analytically, they are indicated in Fig. 1(c) by the two dashed lines . The first vertical line on the right corresponds to the value of waveguide bending B for which the average diffraction 〈 λ 〉 is equal to zero. The bright soliton solutions ex- ist only for positive average diffraction since we assumed that the nonlinearity σ$_{γ}$ is positive. When approaching this limiting value, the spatial extent of solutions grows towards infinity. The second limiting curve in Fig. 1(c) indicates the values corresponding to the critical energy E$_{cr}$ = 2 √ 〈 λ 〉 of the Townes fundamental soliton [31–33] for the average diffraction 〈 λ 〉 . Stable solutions exist only for values of E > E$_{cr}$ . FIG. 1. (a) Schematic of a periodically curved one- dimensional waveguide array. (b) Profile of a spatiotempo- ral light bullet with Gaussian ansatz. (c) Stability diagram in the two-dimensional case. The dashed lines indicate the value of the critical energy E cr and the boundary between positive and negative average diffraction regions 〈 λ 〉 . (d) The dependence of E vs. Δ φ (phase increment after one period) for B = 2 . 1. Solid line corresponds to stable solutions, dashed line to unstable solutions. The dash-dotted line indicates the diffraction map strength S . C. Variational approach We now look for approximate analytical solutions of the Eq. (10) using the variational method [30]. We con- sider the Gaussian Ansatz for the light bullet spatiotem- poral profile [see Fig. 1(b)], For a fixed value of the structural parameter B , there exists a family of soliton solutions. These solutions can be characterized by the pules energy E , the phase accu- mulated over one period Δ φ (this parameter has a similar meaning to the propagation constant for solitons in ho- mogeneous structures), and the diffraction map strength S . We define the values of S analogous to dispersion map strength in structures with varying temporal dis- persion [34, 35]: where A and φ are the amplitude and phase of the soliton, T and W are its temporal and transverse spatial widths, and β and b are the temporal and spatial chirps. The Eq. (10) can be derived from the Lagrangian where W$_{FWHM}$ is the full width at half maximum of the wave packet at the point where it is the narrowest along the x direction. Typical dependencies of E and S on Δ φ are shown in Fig. 1(d). The value of E at Δ φ → 0 is equal to the critical energy of the Townes soliton E$_{cr}$ . We observe that the energy is a monotonously growing func- tion of Δ φ , which indicates that the Vakhitov-Kolokolov criterion [36] for soliton stability is fulfilled. However this criterion is not sufficient to predict stability for period- ically modulated structures, and we find that the solu- tions become unstable for high values of the energy E and diffraction map S . We substitute the ansatz Eq. (13) into Eq. (14) and in- tegrate over spatial coordinates to obtain the reduced Lagrangian which variation with respect to φ yields the energy con- servation dE/dz = 0 with the constant E being propor- S = ∫ L 0 | λ ( z ) | dz W 2 FWHM , (20) u = A ( z ) exp { i φ ( z ) − 1 2 [ x 2 W $^{2}$( z ) + t 2 T $^{2}$( z ) ] + + i 2 [ b ( z ) x 2 + β ( z ) t $^{2}$] } , (13) L = 1 2 ∫ ∫ d x d t [ i ( u$_{z}$u ∗ − u ∗ $_{z}$u ) − λ | u$_{x}$ | 2 − | u$_{t}$ | 2 + | u | $^{4}$] . (14) (4 /π ) L eff = A $^{2}$WT [ 4 φ ′ − b $^{′}$W 2 − β $^{′}$T 2 − W − 2 − DT − 2 − λ ( z ) b $^{2}$W 2 − β $^{2}$T 2 + A 2 $^{$_{/}$√$_{2}$}$] . (15)