# Symmetric and asymmetric solitons trapped in H-shaped potentials

###### Abstract

We report results of numerical and analytical studies of the spontaneous symmetry breaking in solitons, both two- and one-dimensional, which are trapped in H-shaped potential profiles, built of two parallel potential troughs linked by a narrow rung in the transverse direction. This system can be implemented in self-attractive Bose-Einstein condensates (BECs), as well as in a nonlinear bulk optical waveguide.We demonstrate that the introduction of the transverse link changes the character of the symmetry-breaking bifurcation (SBB) in the system from subcritical to supercritical (in terms of the corresponding phase transition, it is a change between the first and second kinds). A noteworthy feature of the SBB in this setting is a non-monotonous dependence of the soliton’s norm at the bifurcation point on the strength of the transverse link. In the full 2D system, the results are obtained in a numerical form. An exact analytical solution is found for the bifurcation in the 1D version of the model, with the transverse rung modeled by the local linear coupling between the parallel troughs with the -functional longitudinal profile. Replacing the -function by its finite-width Gaussian counterpart, similar results are obtained by means of the variational approximation (VA). The VA is also applied to the 1D system with a mixed linear and nonlinear transverse localized coupling. Comparison of the results produced by the different varieties of the system clearly reveals basic features of the symmetry-breaking transition in it.

###### pacs:

05.45.Yv, 03.75.Lm, 42.65.Tg## I Introduction

Symmetric double-well potentials is one of fundamental settings studied in quantum mechanics LL () and in the theory of optical guided-wave propagation, which obeys the Schrödinger equation similar to that known in quantum mechanics NLS (). Counterparts of the double-well potential in optics are represented by directional couplers Snyder (); dual-core (), including various dual-core waveguides created in photonic-crystal matrices PCFcoupler (). It is commonly known that the ground state produced by the linear Schrödinger equation keeps the symmetry of the underlying double-well potential LL (). On the other hand, the introduction of the self-attractive nonlinearity, which transforms the linear equation into the Gross-Pitaevskii equation (GPE) for a Bose-Einstein condensate (BEC) of interacting atoms, loaded into the double-well potential BEC (), or the nonlinear Schrödinger equation (NLSE) modeling nonlinear dual-core waveguides in optics NLS (), leads to the ubiquitous effect of the spontaneous symmetry breaking. As a result, the symmetric ground state is replaced, via the symmetry-breaking bifurcation (SBB; in fact, it is a variety of phase transitions), by an asymmetric state providing for a minimum of the system’s energy, when the strength of the nonlinearity exceeds a critical value. This effect was originally discovered in a discrete model of self-trapping Chris (), and later studied in many settings misc (). Manifestations of the spontaneous symmetry breaking were also studied in detail in BEC models misc2 () and demonstrated experimentally in a self-repulsive BEC Markus ().

In nonlinear optics, the SBB was analyzed in Ref. Snyder () for continuous-wave (spatially uniform) states in the model of dual-core waveguides. For self-trapped modes in dual-core optical systems, i.e., solitons, the bifurcation was studied in Refs. dual-core (); Amir (). A specific example is the SBB for gap solitons in dual-core fiber Bragg gratings Mak (). Later, manifestations of the SBB for matter-wave solitons held in dual-trough potential traps (including solitons of the gap type, supported by a periodic optical-lattice potential) were explored too Arik ()-Luca (). The difference of the dual-trough configuration from the usual double well is the presence of an additional free direction, transverse to that in which the double-well potential acts, see Fig. 1 below.

A characteristic property of solitons trapped in the symmetric dual-trough potential is that, under the action of the self-focusing cubic nonlinearity, they exhibit the SBB of the subcritical type. In that case, branches of the asymmetric states emerge as unstable ones, going backward from the SBB point and getting stable after switching their direction forward at turning points bif (). Thus, the system features the bistability before the bifurcation point, and this mode of the spontaneous symmetry change may be understood as the phase transition of the first kind. On the other hand, the addition of a sufficiently strong periodic potential acting along the troughs changes the character of the SBB in the self-attractive medium from subcritical to supercritical. In the latter case, the asymmetric branches emerge as stable ones, immediately going in the forward direction Arik (); Warsaw2 (), which may also be realized as the phase transition of the second kind. The above-mentioned SBB for gap solitons in the dual-core fiber Bragg grating is of the supercritical type too Mak ().

Another realization of effective double-well potentials is provided by settings based on double-peak spatial modulations of the local nonlinearity coefficient, which may be implemented in optics and BEC alike, see recent review Barcelona () of the topic of solitons in nonlinear potentials. The simplest form of this setting has the nonlinearity concentrated at two points, in the form of a symmetric pair of delta-functions or narrow Gaussians Thawatchai (), as well as a two-dimensional (2D) counterpart of the system, based on the set of two parallel stripes Hung (), or two circles Thaw2D () carrying the self-attractive nonlinearity. The SBB of solitons in the double-well nonlinear potentials also features the symmetry breaking of the subcritical type Thawatchai (); Hung (); Barcelona ().

Further, the spontaneous symmetry breaking was analyzed for 1D and 2D solitons in dual-core discrete systems, with the uniform coupling between two parallel chains Herring (), or with the coupling concentrated at a single site Ljupco (). In the former and latter cases, the SBB is subcritical and supercritical, respectively. Another implementation of the SBB in discrete settings was recently reported for a pair of nonlinear sites embedded into or side-coupled to a linear host lattice Almas ().

The objective of the present work is to consider the SBB of self-attractive localized wave fields trapped in H-shaped potential landscapes, i.e., two parallel troughs linked by a transverse rung, as shown in Fig. 1 below. This configuration can be implemented in effectively two-dimensional BEC, using a set of attractive (red-detuned) laser sheets and/or blue-detuned repelling sheet pairs sheet (), or in BEC layers isolated by a strong standing optical wave standing (), that can also be combined with magnetic trapping fields 2Dcombined (). It is also possible to use the techniques allowing one to “paint” complex potential landscapes (in fact, even more complex than the H-shaped ones that we aim to consider) by rapidly moving laser beams writing (), or induce time-averaged adiabatic landscapes created by means of variable magnetic fields TAAP (). Essentially the same effective potentials can be created in nonlinear optics, using properly patterned photonic-crystal media PCFcoupler (), or a transverse trapping structure permanently written in bulk silica Jena (). We aim to consider both the full 2D model with the transverse H-shaped potential (similar to the 2D model with the dual-trough potentials considered in Refs. Warsaw1 (); Warsaw2 ()) and its simplified 1D counterpart (cf. the 1D version of the dual trough introduced in Ref. Arik ()).

The transverse link (rung) added to the dual-trough potential can be used to control the dynamical properties of the system and, eventually, alter the character of the spontaneous symmetry breaking in the system. In particular, the recent results reported for the discrete system Ljupco () suggest a possibility to change the type of the SBB from sub- to supercritical (i.e., the kind of the respective phase transitions from first to second) by gradually increasing the strength of the transverse link.

The paper is structured as follows. The 2D and 1D models are introduced in Section II, and the full 2D system is considered in Section III by means of numerical methods. A set of bifurcation diagrams indeed demonstrates a switch from the sub- to supercritical SBB in the 2D setting. While one may expect that the strengthening of the linear coupling between the parallel troughs should lead to an increase of the critical value of the nonlinearity strength at which the SBB happens, we observe that, with the introduction of the transverse link, the critical nonlinearity strength at first decreases, due to the change of the character of the bifurcation, and only later starts to grow. In Section IV, we deal with the 1D versions of the system. In that case, the transverse link reduces to a localized linear coupling between two one-dimensional GPEs/NLSEs. In that context, we consider different longitudinal profiles of the coupling. First, approximating it by a -function of the longitudinal coordinate, we report exact analytical solutions for the solitons of all the types—symmetric, antisymmetric, and asymmetric. Accordingly, an exact comprehensive solution for the SBB is available. Next, we consider the coupling localized in a finite interval, with the Gaussian profile. In that case, we develop a variational approximation (VA), and verify its predictions by comparison to numerical findings. Finally, we consider the 1D system which combines the linear and nonlinear localized couplings between the one-dimensional GPEs/NLSEs, both with the Gaussian profile. For that purpose, the VA is developed too. In all the versions of the 1D system, the SBB is always found to be of the supercritical type, unlike the full 2D model that reveals the transition to the supercritical type from the subcritical one. The paper is concluded by Section V.

## Ii The model

### ii.1 The two-dimensional setting

We introduce the model in terms of the two-dimensional GPE for the mean-field wave function of the self-attractive BEC trapped in the H-shaped potential. In physical units, the usual form of the equation is

(1) |

where and are the atomic mass and the respective scattering length, is the confinement length in the transverse direction [, if the confinement is provided by the harmonic-oscillator potential, ], and the trapping potential of depth is defined as follows:

(2) |

As shown in Fig. 1, is the width of the two longitudinal troughs, with separation between them, and is the width of the transverse rung. The norm of the wave function gives the total number of atoms in the condensate,

(3) |

Using scaled variables , , , and , we cast Eq. (1) into the dimensionless form:

(4) |

where the rescaled norm is , and the potential takes the form of

(5) |

with rescaled constants

(6) |

Thus, the system is governed by the set of four dimensionless parameters: norm , scaled potential depth the relative separation between the parallel channels, , and the relative width of the rung, . In the next section, we apply a numerical method to search for symmetric and symmetry-broken localized modes supported by the interplay of the H-shaped trapping potential and self-attractive nonlinearity in Eq. (4). We will identify the SBB in this model, and produce the respective bifurcation diagrams. The predicted results can be readily translated back into physical units be undoing the above transformations.

In the application to the optical guided-wave propagation, Eq. (4) plays the role of the NLSE for the evolution of the amplitude of the guided electromagnetic wave, with time replaced by the propagation distance, . In the latter case, the effective potential represents the modulation of the refractive index in the transverse plane, and the cubic term accounts for the Kerr self-focusing.

### ii.2 The one-dimensional system

The 1D limit of the model corresponds to the case of , i.e., , in terms of the relative parameters defined in Eq. (6). Then, approximating the wave function in each relatively narrow trough by 1D wave functions , the tails of the wave functions channeled by the narrow transverse rung in the transverse direction take the form of

(7) |

[which implies ; otherwise (if the rung is very deep), the effective coupling between the parallel troughs will be stronger than given below by expression (9)]. Then, using the general formalism elaborated for the analysis of the interaction between far separated 2D localized modes 2D (), it is easy to calculate the effective Hamiltonian of the interaction between the parallel troughs, mediated by the channeled tails (7), and thus approximate the full 2D model (4), (5) by the system of coupled 1D equations with the local coupling:

(8) | |||||

with the coupling coefficient identified by equating the above-mentioned 2D interaction Hamiltonian to its counterpart corresponding to the 1D system (8):

(9) |

In fact, by an additional rescaling of Eqs. (8) it is possible to fix , which is adopted below, in Section IV. To understand the genericity of the results produced by the 1D system with the -functional coupling profile, we will also consider the 1D system with the -function replaced by the Gaussian profile with a finite width,

(10) |

keeping the coefficient in front of the Gaussian to be .

## Iii Symmetric and asymmetric two-dimensional solitons

Stationary soliton solutions to Eqs. (4), (5) were found by means of the imaginary-time integration method imaginary (). The stability of the so generated solitons was then tested by direct simulations of the perturbed evolution in real time. Typical examples of stable symmetric and asymmetric 2D solitons are shown in Fig. 2.

Varying the set of control dimensionless parameters, , we identified the critical value, , of the scaled norm , at which the SBB occurs and asymmetric states appear. The natural measure for the asymmetry of such states is defined as

(11) |

As mentioned in the Introduction, the 2D system without the transverse rung, which is tantamount to the present 2D model with Warsaw1 (), as well as its 1D counterpart Arik (), based on the system of 1D equations uniformly coupled by linear terms, feature the SBB of the subcritical type. The difference of the present model is that, with the increase of from zero, the character of the bifurcation quickly switches to supercritical (in other words, the kind of the respective phase transition switches from first to second), as shown in detail by means of the set of bifurcation diagrams in Fig. 3, where control parameters are are fixed, while is varied. In terms of the underlying system, this means the gradual increase of the width of the transverse rung in Fig. 1. As verified by direct simulations (not shown here in detail), all the solution branches displayed in Fig. 3 are dynamically stable.

The transition to the supercritical bifurcation is a consequence of the enhancement of the local transverse coupling between the troughs through the rung. A qualitatively similar change of the bifurcation was recently observed in dual discrete chains, with the transition from the uniform transverse linear coupling Herring () to that at a single site Ljupco ().

The bifurcation picture is further characterized, in Figs. 4 and 5, respectively, by dependences of the value of the norm at the bifurcation point, , on the relative width of the transverse rung, , and the relative distance between the troughs, . Because, as said above, the increase of (making the rung wider) implies the strengthening of the linear coupling between the parallel troughs, one should expect the growth of with . This is, generally, observed in Fig. 4, but after an initial decrease. This surprising feature is explained by the above-mentioned change of the character of the bifurcation, from subcritical to supercritical. On the other hand, the monotonous decrease of with the increase of is a natural consequence of the weakening strength of the coupling between the toughs.

The coupling between the troughs is, obviously, sensitive to the relative width of the separating barrier, with respect to the troughs, which is measured by , see Fig. 1 and Eq. (6). The fact that the bifurcation diagrams, and the respective critical values of the norm, are quite similar in panels (a) and (b) in Figs. 3 and 4, which pertain to and , respectively, demonstrates that the type of the symmetry breaking reported in this section comprises a broad parametric area.

Finally, the simulations demonstrate that the 2D solitons suffer the collapse, in the present model, exactly when it is expected, i.e., when the total scaled norm exceeds the well-known threshold value, Berge ().

## Iv The one-dimensional system: exact, variational, and numerical solutions

### iv.1 Exact solutions for the linear coupling with the delta-functional profile

#### iv.1.1 General analysis

A remarkable feature of the 1D system based on Eq. (8) with the transverse coupling accounted for by the -functions is a possibility to find exact stationary solutions for the solitons of all the types, symmetric, antisymmetric, and asymmetric (note that exact solutions are not available in the discrete counterpart of the system considered in Ref. Ljupco ()). Stationary solutions to Eq. (8) with a given (negative) chemical potential, (in the application to optics, is the propagation constant), are looked for as

(12) |

where real functions and obey the following equations:

(13) | |||||

[recall is fixed in Eq. (8) by means of rescaling]. As follows from Eqs. (13), at the solutions must satisfy the following boundary conditions:

(14) |

where stands for the jump of the derivative at . Solutions to Eqs. (13) and (14) are looked for as

(15) |

where and correspond to the symmetric and antisymmetric states, respectively (or their asymmetric deformations), while shifts and are determined by equations following from the substitution of ansatz (15) into Eqs. (14):

(16) | |||||

Remind is treated here as an arbitrary constant parameterizing the family of solutions.

#### iv.1.2 Symmetric and antisymmetric solutions

Symmetric modes correspond to and , in which case Eq. (16) yields

(17) |

As seen from here, the symmetric solution exists for , i.e., , with the total norm which can be readily calculated:

(18) |

Antisymmetric states, with and , are also possible, with given by the same expression (17) as in the symmetric case. However, the antisymmetric states are expected to be completely unstable, as they correspond to a maximum, rather than minimum, of the respective interaction Hamiltonian, therefore they are not be considered below (the instability of the antisymmetric states is corroborated by numerical tests, see, e.g., Fig. 8 below).

#### iv.1.3 Asymmetric solitons

For the most interesting asymmetric solutions with and , Eq. (16) can be transformed into the following system:

(19) | |||||

(20) |

an explicit solution to which is

(21) | |||||

As seen from these expressions, with the increase of , i.e., with the growth of the norm of the symmetric solutions [see Eq. (18)], the asymmetric modes emerge and exist at

(22) |

with the following values of the norm in the two troughs:

(23) | |||||

Accordingly, the total norm of the asymmetric solution is

(24) |

or, inversely, constant may be expressed in terms of the total norm, which is then treated as the intrinsic parameter of the family of the asymmetric solitons:

(25) |

The corresponding asymmetry measure is

(26) |

Substituting here from Eq. (25), we obtain an eventual expression for the asymmetry parameter as a function of the total norm:

(27) |

which provides a full analytical description of the SBB in the 1D system. The norm of the asymmetric solitons assumes values , with corresponding to the symmetry-breaking point (22).

The bifurcation diagram predicted by Eq. (27), i.e., as a function of , is displayed in Fig. 6. Obviously, the bifurcation revealed by the exact solution is supercritical, being quite similar to the supercritical SBB in the 2D model, which is displayed above in Fig. 3. Direct simulations (with the ideal -functions replaced by their regularized versions, see below) confirm that, as expected, all the solution branches shown in Fig. 6 are stable.

### iv.2 The coupling with the Gaussian profile

To check how generic the exact solution found with the -functional coupling profile is, it is natural to compare the results to those generated by the Gaussian profile (10). Because exact solutions are not available in this case, we will tackle the problem by means of the VA (variational approximation), which is an effective tool for the analysis of a broad class of systems similar to the present one VA ().

To apply the VA, we use the expression for the energy of the 1D system (8), (10):

(28) | |||||

and introduce the following ansatz (hereafter we use instead of subscripts 1 and 2):

(29) |

where is the width of the soliton, and is the asymmetry parameter defined as

(30) |

cf. Eq. (11). Substituting ansatz (29) into expression (28) yields the energy as a function of variational parameters, and :

(31) |

The variational equations following from expression (31), , take the following form:

(32) | |||||

These equations were solved in a numerical form, to find parameters and of stationary solutions as functions of the free parameters, and . The stability of the solutions was verified, in the framework of the VA, by checking if they correspond to a local minimum of energy (31).

The so generated bifurcation diagrams are displayed in the left panel of Fig. 7 for several fixed values of width of the Gaussian profile of the coupling. It is clearly seen that the SBB is supercritical in this case too, being quite similar to that produced by the exact solution for -functional profile, cf. Fig. 6, and to the supercritical diagrams found in the 2D model, cf. Fig. 3. The respective dependence of the norm at the bifurcation point on the width of the Gaussian profile, . In fact, the latter plot is a 1D counterpart of the one which, in the 2D model, shows the dependence of on the relative width of the transverse rung, , cf. Fig. 4.

We have also constructed stationary states of the 1D system with the Gaussian-shaped coupling as numerical solutions of Eqs. (8) and (10). Two different algorithms were used for this purpose, the imaginary-time integration and the Newton-Kantorovich iteration method, both yielding identical results. Thus, three types of stationary configurations were found—symmetric, antisymmetric, and asymmetric ones. Their stability was tested by direct simulations of the evolution in real time. It has been confirmed that the respective SBB is indeed supercritical, and the variational results are very close to the numerically exact ones. It was also found that antisymmetric modes (which do not undergo the bifurcation) are unstable (as mentioned above), see an example in Fig. 8. The two-peak profile of the stationary mode in this figure is a consequence of its antisymmetry (the Gaussian profile of the coupling is effectively repulsive, in this case, which also explains the instability of the mode). It is seen that almost the entire norm is spontaneously transferred into a single trough, where the two peaks split into separating single-component solitons.

We have also developed a similar analysis, using both the VA and numerical solutions, for the 1D model with a rectangular profile of the transverse coupling. The results (not shown here) are very similar to those generated by the Gaussian profile.

## V The combined linear-nonlinear coupling with the Gaussian profile

It may also be relevant to take into regard a nonlinear correction to the linear coupling between the two troughs, in the framework of the 1D system (note that the full 2D model automatically takes into regard both linear and nonlinear effects of the coupling). To this end, we adopt the following modification of the 1D system (8), with the Gaussian coupling profile (10):

(33) | |||||

where is the relative strength of the nonlinear coupling. In this case, the energy functional is

(34) |

To apply the VA method, we adopt the same Gaussian ansatz (29) which was used above. The substitution of the ansatz into functional (34) yields

(35) |

the corresponding variational equations, , being

The bifurcation diagram produced by the numerical solution of Eqs. (LABEL:LN1L), and the respective dependence of the total norm at the SBB point on are displayed (for ) in Fig. 9.

The comparison of the dependences , obtained in the three versions of the1D model with the finite width of the coupling profile (Gaussian, rectangular—which is not displayed here in detail—and mixed linear-nonlinear) suggests that the inclusion of the nonlinear coupling into Eq. (33) makes this characteristic essentially more similar to its counterpart obtained in the full 2D model, cf. Fig. 4: While in Fig. 7, and in its counterpart obtained in the model with the rectangular profile, the curves are convex, they are concave in both Figs. 4 and 9. This conclusion is naturally explained by the above-mentioned fact that the full 2D model automatically takes into account both the linear and nonlinear coupling between the parallel troughs.

## Vi Conclusion

The objective of this work is to extend the study of the spontaneous symmetry breaking in 1D and 2D solitons trapped in the systems of two parallel tunnel-coupled potential troughs. Unlike recently studied models, here we introduce the H-shaped potential, and focus on the SBB (symmetry-breaking bifurcation) for solitons in this system. The analysis of both the 2D and 1D versions of the model demonstrates that the transverse link (the rung of the H-shaped potential profile) transforms the bifurcation from subcritical into the supercritical one. In other words, the rung controls the switching of the corresponding symmetry-breaking phase transition from the first into second kind. A nontrivial manifestation of the change of the SBB type is the non-monotonous dependence of the critical value of the soliton’s norm, at the bifurcation point, on the strength of the transverse link: prior to the expected growth, it demonstrates a region of the decrease. In the full 2D model, the results were obtained in the numerical form. On the other hand, the 1D version with the -functional profile of the local linear coupling between the troughs admits the exact analytical solutions for the solitons of all the types—symmetric, antisymmetric, and asymmetric—and, accordingly, the bifurcation diagram was obtained in the exact analytical form, which is a rare possibility. In other variant of the 1D model, with the Gaussian profile of the local transverse linear coupling, the results were obtained by means of the VA (variational approximation). The VA was also used to solve the bifurcation problem in the system combining the linear and nonlinear localized transverse coupling. The set of the results reported in the paper provides for a comprehensive description of the symmetry breaking in the H-shaped system. In particular, adding the nonlinear coupling to the 1D model makes the characteristics of the symmetry breaking closer to those found in the full 2D model, which automatically incorporates linear and nonlinear coupling effects.

The settings considered in this work can be realized in the self-attractive BEC, and in self-focusing bulk optical media, with the appropriate transverse profile of the refractive index. It is relevant to stress that the description of the BEC based on the coupled GPEs is entirely based on the mean-field approximation. While it is well known that this approximation provides for an exceptionally accurate description of virtually all matter-wave patterns observed in experiments, quantum fluctuations and other beyond-mean-field effects being detectable only under special conditions BEC (), one may expect that fluctuations can be amplified in a vicinity of the phase transition (i.e., of the SBB), which suggests a question, whether the GPEs are applicable in this case. A full analysis of this issue should be a subject for a separate extended work; nevertheless, a relevant analogy is suggested by the analysis of fluctuations in the vicinity of the SBB, both in the continuous-wave and soliton settings, which was performed in the framework of the model of the directional coupler in nonlinear optics, based on the coupled 1D NLSEs, in Ref. Amir (). The conclusion was that, on the contrary to the anticipations, the fluctuational effects remain virtually negligible in that case, their amplification in the vicinity of the respective SBBs producing no conspicuous changes against the mean-field description. Based on this analogy, we expect that the mean-field description remains valid in the present case too.

One possibility for further development of the analysis is to perform it in a two-component system. Another interesting extension may be carried out for a configuration with two parallel transverse rungs, in which case one may expect a double spontaneous symmetry breaking: between the parallel potential troughs, and, independently, between vicinities of the two rungs. In particular, a challenging problem would be to find an exact solution in the case when the pair of the rungs is represented by a symmetric set of two -functions.

## References

- (1) D. Landau and E. M. Lifshitz, Quantum Mechanics (Moscow: Nauka Publishers, 1974).
- (2) Y. S. Kivshar and G. P. Agrawal, Optical Solitons: From Fibers to Photonic Crystals (Academic Press: San Diego).
- (3) A. W. Snyder, D. J. Mitchell, L. Poladian, D. R. Rowland, and Y. Chen, J. Opt. Soc. Am. B 8, 2101 (1991).
- (4) S. Trillo, S. Wabnitz, E. M. Wright, and G. I. Stegeman, Opt. Lett. 13, 672 (1988); C. Paré and M. Fłorjańczyk, Phys. Rev. A 41, 6287 (1990); A. I. Maimistov, Kvant. Elektron. 18, 758 [Sov. J. Quantum Electron. 21, 687 (1991)]; M. Romagnoli, S. Trillo, and S. Wabnitz, Opt. Quantum Electron. 24, S1237 (1992); N. Akhmediev and A. Ankiewicz, Phys. Rev. Lett. 70, 239 (1993); P. L. Chu, B. A. Malomed, and G. D. Peng, Opt. Lett. 18, 328 (1993); Yu. S. Kivshar and M. L. Quiroga-Teixeiro, ibid. 18, 980 (1993).
- (5) A. Mostofi, B. A. Malomed, and P. L. Chu, Opt. Commun. 137, 244 (1997); ibid. 145, 274 (1998).
- (6) B. H. Lee, L. B. Eom, J. Kim, D. S. Moon, U. C. Paek, and G. H. Yang, Opt. Lett. 27, 812 (2002); S. Boscolo, M. Midrio, and C. G. Someda, IEEE J. Quant. Electr. 38, 47 (2002); A. Martinez, F. Cuesta, and J. Marti, IEEE Phot. Tech. Lett. 15, 694 (2003); M. Thorhauge, L. H. Frandsen, and P. I. Borel, Opt. Lett. 28, 1525 (2003); K. Saitoh, Y. Sato, and M. Koshiba, Opt. Exp. 11, 3188 (2003); D. Mori and T. Baba, ibid. 13, 9398 (2005). Quiroga-Teixeiro, J. Opt. Soc. Am. B 12, 898 (1995)
- (7) S. Giorgini, L. P. Pitaevskii, and S. Stringari, Rev. Mod. Phys. 80, 1215 (2008); H. T. C. Stoof, K. B. Gubbels, and D. B. M. Dickrsheid, Ultracold Quantum Fields (Springer: Dordrecht, 2009).
- (8) J. C. Eilbeck, P. S. Lomdahl, and A. C. Scott, Physica D 16, 318 (1985).
- (9) E. A. Ostrovskaya, Y. S. Kivshar, M. Lisak, B. Hall, F. Cattani, and D. Anderson, Phys. Rev. A 61, 031601(2000); R. D’Agosta, B. A. Malomed, C. Presilla, Phys. Lett. A 275, 424 (2000); R. K. Jackson and M. I. Weinstein, J. Stat. Phys. 116, 881 (2004); D. Ananikian and T. Bergeman, Phys. Rev. A 73, 013604 (2006); E. W. Kirr, P. G. Kevrekidis, E. Shlizerman, and M. I. Weinstein, SIAM J. Math. Anal. 40, 566 (2008).
- (10) G. J. Milburn, J. Corney, E. M. Wright, and D. F. Walls, Phys. Rev. A 55, 4318 (1997); A. Smerzi, S. Fantoni, S. Giovanazzi, and S. R. Shenoy, Phys. Rev. Lett. 79, 4950 (1997); S. Raghavan, A. Smerzi, S. Fantoni, and S. R. Shenoy, Phys. Rev. A 59, 620 (1999); K. W. Mahmud, H. Perry, and W. P. Reinhardt, Phys. Rev. A 71, 023615 (2005); E. Infeld, P. Ziń, J. Gocalek, and M. Trippenbach, Phys. Rev. E 74, 026610 (2006); G. Theocharis, P. G. Kevrekidis, D. J. Frantzeskakis, and P. Schmelcher, Phys. Rev. E 74, 056608 (2006); G. L. Alfimov and D. A. Zezyulin, Nonlinearity 20, 2075 (2007); C. Wang, P. G. Kevrekidis, N. Whitaker, and B. A. Malomed, Physica D 237, 2922 (2008).
- (11) M. Albiez, R. Gati, J. Fölling, S. Hunsmann, M. Cristiani, and M. K. Oberthaler, Phys. Rev. Lett. 95, 010402 (2005); R. Gati, M. Albiez, J. Fölling, B. Hemmerling, and M. K. Oberthaler, Appl. Phys. B 82, 207 (2006).
- (12) W. C. K. Mak, B. A. Malomed, and P. L. Chu, J. Opt. Soc. Am. B 15, 1685 (1998).
- (13) A. Gubeskys and B. A. Malomed, Phys. Rev. A 75, 063602 (2007).
- (14) M. Matuszewski, B. A. Malomed, and M. Trippenbach, Phys. Rev. A 75, 063621 (2007).
- (15) M. Trippenbach, E. Infeld, J. Gocalek, M. Matuszewski, M. Oberthaler, and B. A. Malomed, Phys. Rev. A 78, 013603 (2008).
- (16) L. Salasnich, B. A. Malomed, and F. Toigo, Phys. Rev. A 81, 045603 (2010).
- (17) G. Iooss and D. D. Joseph, Elementary Stability Bifurcation Theory (Springer-Verlag: New York, 1980).
- (18) Y. V. Kartashov, B. A. Malomed, and L. Torner, Solitons in nonlinear lattices, Rev. Mod. Phys. 83, 247 (2011).
- (19) T. Mayteevarunyoo, B. A. Malomed, and G. Dong, Phys. Rev. A 78, 053601 (2008); N. Dror and B. A. Malomed, ibid. A 83, 033828 (2011).
- (20) L. C. Qian, M. L. Wall, S. Zhang, Z. Zhou, and H. Pu, Phys. Rev. A 77, 013611 (2008).
- (21) N. V. Hung, P. Ziń, M. Trippenbach, and B. A. Malomed, Phys. Rev. E 82, 046602 (2010).
- (22) T. Mayteevarunyoo, B. A. Malomed, and A. Reoksabutr, Spontaneous symmetry breaking of photonic and matter waves in two-dimensional pseudopotentials, J. Mod. Opt., in press.
- (23) G. Herring, P. G. Kevrekidis, B. A. Malomed, R. Carretero-González, and D. J. Frantzeskakis, Phys. Rev. E 76, 066606 (2007).
- (24) Lj. Hadžievski, G. Gligorić, A. Maluckov, and B. A. Malomed, Phys. Rev. A 82, 033806 (2010); M. D. Petrović, G. Gligorić, A. Maluckov, Lj. Hadžievski, and B. A. Malomed, Phys. Rev. E 84, 026602 (2011).
- (25) E. Bulgakov, K. Pichugin, and A. Sadreev, Phys. Rev. B 83, 045109 (2011); V. A. Brazhnyi and B. A. Malomed, Phys. Rev. A 83, 053844 (2011).
- (26) V. Bagnato and D. Kleppner, Phys. Rev. A 44, 7439 (1991).
- (27) R. Dumke, T. Muther, M. Volk, W. Ertmer, G. Birkl, Phys. Rev. Lett. 89, 220402 (2002); D. Rychtarik,B. Engeser, H. C. Nagerl, and R. Grimm, ibid. 92, 173003 (2004); R. Bucker, A, Perrin, S. Manz, T. Betz, C. Koller, T. Plisson, J. Rottmann, T. Schumm, and J. Schmiedmayer, New J. Phys. 11, 103039 (2009).
- (28) M. Greiner, I. Bloch, O. Mandel, T. W. Hansch, and T. Esslinger, Phys. Rev. Lett. 87, 160405 (2001).
- (29) K. Dieckmann, R. J. C. Spreeuw, M. Weidemüller, and J. T. M. Walraven, Phys. Rev. A 58, 3891 (1998); A. Gorlitz, J. M. Vogels, A. E. Leanhardt, C. Raman, T. L. Gustavson, J. R. Abo-Shaeer, A. P. Chikkatur, S. Gupta, S. Inouye, T. Rosenband, and W. Ketterle, Phys. Rev. Lett. 87, 130402 (2001).
- (30) K. Henderson, C. Ryu, C. MacCormick, and M. G. Boshier, New J. Phys. 11, 043030 (2009).
- (31) M. Gildemeister, E. Nugent, B. E. Sherlock, M. Kubasik, B. T. Sheard, and C. J. Foot, Phys. Rev. A 81, 031402(R) (2010).
- (32) A. Szameit, J. Burghoff, T. Pertsch, S. Nolte, A. Tünnermann, and F. Lederer, Opt. Express 14, 6055 (2006).
- (33) B. A. Malomed, Phys. Rev. E 58, 7928 (1998).
- (34) M. L. Chiofalo, S. Succi, and M. P. Tosi, Phys. Rev. E 62, 7438 (2000); W. Z. Bao and Q. Du, SIAM J. Scient. Comput. 25, 1674 (2004).
- (35) L. Bergé, Phys. Rep. 303, 259 (1998).
- (36) B. A. Malomed, Progr. Optics 43, 71 (E. Wolf, editor: North Holland, Amsterdam, 2002); http://www.sciencedirect.com/science/article/pii/S0079663802800269.