Perspective Chapter: Lattice Solitons in a Nonlocal Nonlinear Medium with Self-Focusing and Self-Defocusing Quintic Nonlinearity

Abstract

The fundamental lattice solitons are explored in a nonlocal nonlinear medium with self-focusing and self-defocusing quintic nonlinearity. The band-gap boundaries, soliton profiles, and stability domains of fundamental solitons are investigated comprehensively by the linear stability spectra and nonlinear evolution of the solitons. It is demonstrated that fundamental lattice solitons can stay stable for a wide range of parameters with the weak self-focusing and self-defocusing quintic nonlinearity, while strong self-focusing and self-defocusing quintic nonlinearities are shortened the propagation distance of evolved solitons. Furthermore, it is observed that when the instability emerges from strong quintic nonlinearity, increasing anisotropy of the medium and modification of lattice depth can be considered as a collapse arrest mechanism.

1. Introduction

Optical solitons arise as localized optical fields that preserve their shapes during propagation when the medium’s diffraction and self-phase modulation are balanced [1]. Since they were first observed experimentally in 2003 [2], spatial solitons [3] in nonlinear optical systems with additional optical lattices have received a significant interest. Currently, fundamental solitons and vortices on a variety of lattices including real and complex structures in the media with cubic (Kerr) [4, 5], saturable [6] and with competing nonlinearities [7] have been studied. Solitons may also exist in aperiodic or quasicrystal lattices [8, 9, 10, 11, 12] and in lattices with defects [13] and dislocations [14, 15].

In the studies above, the governing equations are nonlinear cubic Schrödinger (NLS) equation and their variants, e.g., NLS with additional terms and/or with external lattices. On the other hand, numerous optical materials, such as potassium niobate (KNbO₃) [16] or lithium niobate (LiNbO₃) [17], acquire quadratic nonlinear responses as well [18, 19, 20, 21, 22, 23, 24]. NLS equation with coupling to a mean term (NLSM) system governs the pulse dynamics in quadratically polarized media. The general NLSM system is given by [18, 25, 26].

$${\mathit{iu}}_z+\mathrm{\Delta }u+{\left|u\right|}^{2}u-\mathit{\rho u\varphi }=0,{\varphi }_{\mathit{xx}}+{\mathit{\nu \varphi }}_{\mathit{yy}}={\left({\left|u\right|}^2\right)}_{\mathit{xx}}$$E1

Here $u\left(x,y\right)$ is the normalized amplitude of the envelope of the static electric field propagating in the $z$ direction. $\mathrm{\Delta }u\equiv u_{\mathit{xx}}+u_{\mathit{yy}}$ corresponds to diffraction, and the nonlinear cubic term comes from the nonlinear Kerr effect. The parameter $\rho$ is a coupling constant that emerges from the combined optical rectification and electro-optic effects formed by the $\varphi \left(x,y\right)$ field, and $\nu$ is the coefficient that originates from the anisotropy of the material [26]. These equations arise due to the interaction between the fundamental and dc fields when second-harmonic-generation is not phase-matched. In such a situation, the second harmonic component can be solved explicitly and generates an additional self-phase modulation contribution due to cascaded nonlinearity. Consequently, the NLSM system is obtained as a result of the nonlocal nonlinear coupling between the first field (with the cascaded effect from the second harmonic) and a static field that arises from the zeroth harmonic (mean term). We would like to note that solitons of a media with solely quadratic nonlinear response can also be governed by the NLSM type Equations [25, 26].

NLSM equations were derived by Benney and Roskes [27] for a finite water depth, neglecting the surface tension. In 1974, Davey and Stewartson [28] investigated the evolution of a three-dimensional wave packet for a finite water depth and obtained an equivalent form of these equations. The integrability of the NLSM equations that were derived for the shallow water limit was studied by Ablowitz and Haberman [29]. Djordevic and Reddekopp [30] improved the study of Benney and Roskes by including the surface tension in 1977. Later, Ablowitz et al. [18, 25, 26] derived from first principles that NLSM type equations model the evolution of the electromagnetic field in a quadratically polarized media.

One major drawback of the dynamics of solitons under the NLS equation in 2D is that they exhibit collapse. As a matter of fact, in [31], it was revealed that the collapse dynamics of NLSM solitons are similar to the collapse of NLS solitons. Possible collapse arrest mechanisms have been studied extensively in nonlinear optics, e.g., nonlinear saturation [32, 33]. Merle and Raphael [34] investigated the collapse of the NLS solitons and its variants in depth. Furthermore, Gaeta and coworkers [35] carried out detailed experiments in order to expose the nature of the singularity formation in cubic optical media and demonstrated experimentally that collapse occurs with a self-similar profile.

Recently, wave collapse in the NLSM system was studied in [36] and it is shown that in both water waves and optics, collapse occurs with a quasi self-similar profile. By evolving the NLSM solutions that are computed by the Spectral Renormalization method which is essentially a fixed point iteration scheme. In this study, it is also revealed that the NLSM solitons have astigmatic profiles and their collapse can be arrested by adding nonlinear saturation into the system. Wave self-rectification and beam ellipticity as a collapse arrest mechanism for NLSM solitons was put forward by in [16]. In the aforementioned work, for simplicity, the authors considered Gaussian profiles with various input powers and astigmatism and then evolve these profiles for long distances. More recently, wave collapse in the NLSM system was arrested by adding a periodic external lattice in [23] and it has been shown that deeper lattices may serve as a collapse arrest mechanism.

In this chapter, the NLSM model (1) is extended by adding a quintic term and an external lattice. Using this new model, the dynamics of the fundamental lattice solitons in a nonlocal nonlinear medium with quintic nonlinearity are explored. The study is concentrated on effects formed by the imposed lattice depth and the strength of quintic nonlinearity. The stability of the obtained lattice solitons is examined by the nonlinear evolution and linear spectra. In the light of the conducted numerical analysis, it is confirmed that fundamental lattice solitons can exist in nonlocal nonlinear media with both self-focusing and self-defocusing quintic nonlinearities and stability of these solitons is achieved for a broad range of system parameters. Since many optical materials such as chalcogenide glasses are engineered to reveal fifth and seventh-order effects in addition to cubic nonlinear effects [37], and high-order nonlinearities can arise even with pure Kerr materials [38, 39, 40], it is crucial to consider the soliton dynamics consists of higher-order nonlinearities.

The outline of the chapter is as follows: In Section 2, we present the model equations and compute the periodic lattice solitons numerically. In Section 3, we explore the linear spectra and nonlinear evolution of the solitons. Results are outlined in Section 4.

2. The model and its soliton solutions

It is known that, the steady-state solutions of the NLSM system (1) collapse [36], and the collapse of solitons can be arrested by adding external lattices [23]. The NLSM model (1) is extended as follows to describe the nonlocal qubic-quintic nonlinear medium with an external lattice,

$${\mathit{iu}}_z+\frac{1}{2}\mathrm{\Delta }u+\beta {\left|u\right|}^{2}u-\mathit{\rho u\varphi }+\gamma {\left|u\right|}^{4}u-V\left(x,y\right)u=0,{\varphi }_{\mathit{xx}}+{\mathit{\nu \varphi }}_{\mathit{yy}}={\left({\left|u\right|}^2\right)}_{\mathit{xx}}$$E2

where $\gamma$ denotes the strength of quintic nonlinearity and $V(x,y$) shows the optical lattice. We consider lattices that are created by the sum of $N$ phase-modulated plane waves [8].

$$V\left(x,y\right)=\frac{V_0}{N^2}{\left|\sum _{n=0}^{N-1}{e}^{i\left(k_x^{n}x+k_y^{n}y\right)}\right|}^2$$E3

where $V_0>0$ is the depth of lattice and the wave vector $\left(k_x^n,k_y^n\right)=\left[\mathit{\text{Kcos}}\left(2\mathit{\pi n}/N\right),\mathit{\text{Ksin}}\left(2\mathit{\pi n}/N\right)\right]$. The lattices for $N=\mathrm{2,3,4,6}$ correspond to crystal (periodic) structures, while $N=5,7$ correspond to aperiodic (Penrose) quasi-crystals. Contour plot and diagonal cross-section of the lattice $V\left(x,y\right)$ is displayed in Figure 1 for $V_0=12.5,N=4$ and $k_x=k_y=2\pi$. It can be seen that the lattice is periodic, and the center ($\left(x,y\right)=\left(0,0\right)$) of lattice is a local maximum.

Figure 1.

The top view and diagonal cross-section of the lattice $V\left(x,y\right)$ when $V_0=12.5,N=4$ and $\left(x,y\right)\in \left[-15,15\right]$.

2.1 Numerical solution for the fundamental solitons

The squared operator (SOM) method (that was proposed by Yang and Lakoba [41]) is utilized to calculate the soliton solutions of the (2 + 1)D NLSM model (2). The SOM method is based on integrating squared-operators of evolution equations. Derivation of these operators and the scheme of the method are explained below.

A real-valued nonlinear evolution equation can be presented in the following form

$${\mathrm{L}}_0\mathrm{u}\left(\mathrm{x}\right)=0$$E4

where $\mathbf{x}$ is the multidimensional spatial variable, $\mathbf{u}\left(\mathbf{x}\right)$ is a real valued function and the operator ${\mathbf{L}}_0$ includes the solitary wave’s propagation constant. Let the operator ${\mathbf{L}}_1$ denotes linearization of Eq. (4) around the solution $\mathbf{u}$, and given by

$${\mathrm{L}}_1\left(\mathrm{u}+\stackrel{\sim }{\mathrm{u}}\right)={\mathrm{L}}_1\stackrel{\sim }{\mathrm{u}}+O\left({\stackrel{\sim }{\mathrm{u}}}^2\right)$$E5

where $\stackrel{\sim }{u}\ll 1$. In order to get a solitary wave solution, the following time-dependent squared operator evolution equation will be integrated

$${\mathrm{u}}_t=-{\mathrm{M}}^{-\mathbf{1}}{{\mathrm{L}}_{\mathbf{1}}}^{†}{\mathrm{M}}^{-\mathbf{1}}{\mathrm{L}}_{\mathbf{0}}\mathrm{u}$$E6

here ${{\mathbf{L}}_{\mathbf{1}}}^{†}$ denotes the Hermitian of the operator and M is a real valued positive definite Hermitian operator that is introduced to accelerate the convergence. Using the forward Euler method, the solitary wave solution $\mathbf{u}$ will be calculated by the following iteration.

$${\mathrm{u}}_{\mathrm{n}+\mathbf{1}}={\mathrm{u}}_{\mathrm{n}}-{\left[{\mathrm{M}}^{-\mathbf{1}}{{\mathrm{L}}_{\mathbf{1}}}^{†}{\mathrm{M}}^{-\mathbf{1}}{\mathrm{L}}_{\mathbf{0}}\mathrm{u}\right]}_{\mathrm{u}\mathbf{=}{\mathrm{u}}_{\mathrm{n}}}\mathrm{\Delta }t.$$E7

It was shown that the SOM method converges to a solitary wave solution for a broad-range of nonlinear evolution equations when a convenient initial condition is given and time step $\mathrm{\Delta }t$ is small enough [41, 42].

To calculate a solitary wave solution of the model Eq. (2) by the SOM, the following scheme is constructed. Inserting the solution suggestion $u=U\left(x,y\right)exp\left(\mathit{i\mu z}\right)$ into the model (2), the following sub-operators are obtained

$$\begin{array}{l}{F}_0=-\mu +\beta {\left|U\right|}^2-\mathit{\rho \varphi }+\gamma {\left|U\right|}^4-V\left(x,y\right),{\varphi }_{\mathit{xx}}+{\mathit{\nu \varphi }}_{\mathit{yy}}={\left({\left|U\right|}^2\right)}_{\mathit{xx}}\\ F_1=-\mu +3\beta U^2-\mathit{\rho \varphi }+5\gamma U^4-V\left(x,y\right)\end{array}$$E8

where $\mu$ is propagation constant. Using $F_0$ and $F_1$, we get ${\mathbf{L}}_0$ and ${\mathbf{L}}_1$ operators as follows:

$$\begin{array}{l}{\mathbf{L}}_{0}U=\frac{1}{2}\mathrm{\Delta }U+F_{0}U,M_0={\mathcal{F}}^{-1}\left(\frac{\mathcal{F}\left({\mathbf{L}}_{0}U\right)}{K^2+c}\right),\\ {\mathbf{L}}_{1}U=\frac{1}{2}\mathrm{\Delta }{M}_0+F_1{M}_0,M_1={\mathcal{F}}^{-1}\left(\frac{\mathcal{F}\left({\mathbf{L}}_{1}U\right)}{K^2+c}\right)\end{array}$$E9

where

$$\mathcal{F}\left\{f\left(x,y\right)\right\}=\widehat{f}\left(k_x,k_y\right)={\int }_{-\infty }^{\infty }f\left(x,y\right)e^{i\left(k_{x}x+k_{y}y\right)}\mathit{\text{dxdy}}$$E10

$${\mathcal{F}}^{-1}\left\{\widehat{f}\left(k_x,k_y\right)\right\}=f\left(x,y\right)=\frac{1}{2\pi }{\int }_{-\infty }^{\infty }\widehat{f}\left(k_x,k_y\right)e^{-i\left(k_{x}x+k_{y}y\right)}{\mathit{dk}}_x{\mathit{dk}}_y,$$E11

$K^2=k_x^2+k_y^2$ and the parameter $c$ that has a considerable effect on convergence of the SOM method is utilized for parametrizing the numerical procedure.

After obtaining $M_1$, the iteration scheme is executed as follows,

$$U_{n+1}=U_n-M_1\mathrm{\Delta }t.$$E12

while the mean term $\varphi \left(x,y\right)$ is calculated by

$${\varphi }_n={\mathcal{F}}^{-1}\left(\frac{k_x^2}{k_x^2+\nu k_y^2}\mathcal{F}\left({\left|U_n\right|}^2\right)\right).$$E13

In order to avoid division by zero error, the first element of $K^2$ is set to be 1. Starting from an initial guess, this numerical scheme is iterated until the error $E=\sqrt{\parallel U_{n+1}-U_n{\parallel }^2}<{10}^{-8}$. It is noted that $c$ and $\mathrm{\Delta }t$ are chosen heuristically as positive real numbers to obtain a convergent solution, and we take $c=1.3$ and $\mathrm{\Delta }t=0.2$ in this chapter. The initial condition of the SOM method is chosen as a Gaussian $u\left(x,y\right)=exp\left(-\left[{\left(x-x_0\right)}^2+{\left(y-y_0\right)}^2)\right]\right)$. The location of the initial condition is determined by $x_0$ and $y_0$. In previous studies [8, 15, 23], it was shown that the solitons centered at the maxima of the lattices cannot stay stable. Thus, in this study, the initial condition is located on a local minimum of the periodic lattice where $\left(x_0,y_0\right)=\left(\pi ,0\right)$ (that is shown in Figure 1). Unless otherwise stated, the parameters in the model (2) are fixed to the following values

$$\left(\mu ,\rho ,\nu ,\beta ,\gamma ,V_0\right)=\left(-\mathrm{0.1,0.5,1.5,2},\pm \mathrm{0.2,12.5}\right).$$E14

It should be noted that $\rho =0.5$ and $\nu =1.5$ are specifically chosen to characterize the electro-optical effects in potassium niobate (KNbO₃) [16], and $\gamma =0.2$ and $\gamma =-0.2$ cases are chosen to investigate soliton dynamics with self-focusing and self-defocusing quintic nonlinearity, respectively.

The fundamental solitons of the NLSM model (2) are obtained by the SOM method for the considered parameters (14) and are shown in Figure 2.

Figure 2.

3D profiles (first column), top views (second column) of fundamental solitons that are obtained (a) for $\gamma =-0.2$ (self-defocusing) and (b) for $\gamma =0.2$ (self-focusing) quintic nonlinearity. All other parameters are taken as in Eq. (14).

As shown in previous studies [23, 24, 36], fundamental solitons of the NLSM system (2) (with or without lattice) are not radially symmetric due to the anisotropy in the medium. To investigate the effect of the quintic nonlinearity on the level of astigmatism in the solitons, the following formulation is defined as a measure of astigmatism

$$e=\frac{\text{radius}\text{along}y-\text{axis}}{\text{radius}\text{along}x-\text{axis}}.$$E15

Here, $e=1$ when the soliton is radially symmetric, and if $e<1$ and $e>1$ the obtained solitons are comparatively wider along the $x$ and $y$ axes, respectively. Thus, the solitons are elliptical if $e\ne 1$.

Contour plots of solitons are depicted in Figure 3 for $\rho =0$, $\rho =0.5$ and $\rho =1$ when the quintic nonlinearity coefficient $\gamma =-0.5$, $\gamma =0$ and $\gamma =0.5$. Figure 3 shows that astigmatism along the $x$-axis is emphasized when the fundamental solitons are obtained with a larger coupling parameter $\rho$ (see the first column). It is also observed that the fundamental soliton becomes less astigmatic along $x$-axis with a larger quintic nonlinearity parameter $\gamma$ (for a fixed $\rho$), and when $\rho =0$ and the quintic nonlinearity is self-focusing (i.e., when $\gamma >0$) the solitons become comparatively wider along $y$-axis (see the first row).

Figure 3.

Contour image of fundamental solitons that are obtained for $\rho =0$ (the first row), $\rho =0.5$ (the second row) and $\rho =1$ (the third row) while $\gamma =-0.5$ (the first column), $\gamma =0$ (the second column) and $\gamma =0.5$ (the third column). All other parameters are taken as in Eq. (14) and $e$ shows the degree of astigmatism for each soliton.

To see the impact of the anisotropy parameter $\nu$ on the soliton profile, the on-axes amplitudes and contour images of fundamental solitons are displayed in Figure 4. Figure 4 shows that the amplitude of fundamental solitons are decreasing with increasing $\nu$ (from 0 to 3), and solitons become more extended along $x$-axis with larger $\nu$ values.

Figure 4.

The on-axes amplitudes of fundamental solitons along x-axis (first column), along y-axis (second column) and contour image of solitons (third row) (a) for $\nu =0$; (b) for $\nu =1.5$ and (c) $\nu =3$. All other parameters are taken as in Eq. (14) and $e$ shows the degree of astigmatism for each soliton.

Comparing the contour images in Figures 3 and 4, it is seen that the astigmatism of the fundamental solitons changes significantly with $\rho$, whereas it depends weakly on $\nu$.

3. Stability analysis

The stability properties of fundamental solitons (that are obtained by the SOM method) are examined by the linear eigenvalue spectra and the direct simulation of the nonlinear model (2).

The linear stability of solitons are studied by the linear spectra of fundamental solitons. To obtain the linear spectrum, the fundamental soliton is perturbed as follows

$$U=e^{\mathit{i\mu z}}\left[u_0\left(x,y\right)+g\left(x,y\right)e^{\mathit{\lambda z}}+h^{\ast }\left(x,y\right)e^{{\lambda }^{\ast }z}\right]$$E16

where $u_0\left(x,y\right)$ is the fundamental soliton and $g,h\ll 1$ are infinitesimal perturbations. Substituting the solution $U$ into (2) and linearizing, it is found that these normal modes satisfy the following linear eigenvalue problem

$$\mathcal{L}\mathbf{V}=\lambda \mathbf{V}$$E17

where

$$\mathcal{L}=i\left(\begin{array}{cc}0& {\mathcal{L}}_{+}\\ {\mathcal{L}}_{-}& 0\end{array}\right),\mathbf{V}=\left(\begin{array}{l}g\\ h\end{array}\right)$$

and the matrix coefficients of $\mathcal{L}$ are

$$\begin{array}{l}{\mathcal{L}}_{+}=-\mu +\frac{1}{2}\mathrm{\Delta }-\mathit{\rho \varphi }+3\beta u^2+5\gamma u^4-V,\\ {\mathcal{L}}_{-}=-\mu +\frac{1}{2}\mathrm{\Delta }-\mathit{\rho \varphi }+\beta u^2+\gamma u^4-V.\end{array}$$E18

We can solve this eigenvalue problem by the Fourier collocation method [42]. Any eigenvalue with a positive real part in the spectrum indicates the linear instability of the fundamental soliton considered.

The power (or energy is some other contexts) of solitons can be computed by P=${\iint }_{-\infty }{\left|u\right|}^2\mathit{\text{dxdy}}$. In [43], Vakhitov and Kolokolov (VK) demonstrated that there is a critical power value above which the solitons cannot be linearly stable, and a soliton can be linearly stable only if its power increases as the propagation constant $\mu$ is increased, i.e.,

$$\mathit{dP}/\mathit{d\mu }>0.$$E19

Moreover, in [44, 45], it was shown that the solitons can be stable nonlinearly only if the slope condition (19) is satisfied.

In the light of these results, prior to the nonlinear stability analysis, the power and linear stability of solitons are studied in Figure 5. These power-eigenvalue diagrams are calculated for the parameters that are given in Eq. (14), and the linear (in)stability intervals of gap solitons are determined via the computation of eigenvalue spectra of the solitons at each point of the power curves (see Figure 5(a1)–(b1)). Furthermore, the maximum real part in the eigenvalue spectra (spectral radius) is examined in Figure 5(a2)–(b2) for the gap solitons. The solitons are found to be linearly stable for self-defocusing quintic nonlinearity ($\gamma =-0.2$) when the power $P\in \left[\mathrm{1.61,2.39}\right]$ and propagation constant $\mu \in \left[-\mathrm{0.35,0}\right]$ (see Figure 5(a1)–(a2)), and the solitons are linearly stable when $P\in \left[\mathrm{1.76,2.17}\right]$ and $\mu \in \left[-\mathrm{0.15,0.3}\right]$ for self-focusing ($\gamma =0.2$) quintic nonlinearity (see Figure 5(b1)–(b2)). These results are consistent with VK stability criterion since the slope of both power-eigenvalue diagrams is positive.

Figure 5.

The power-eigenvalue ($P-\mu$) diagram of lattice solitons (a1) for $\gamma =-0.2$ (self-defocusing) and (b1) for $\gamma =0.2$ (self-focusing). Maximum real part (spectral radius) of the eigenvalue spectra of solitons (a2) for $\gamma =-0.2$ and (b2) for $\gamma =0.2$. Solid blue and red dotted lines show stable and unstable regions for the gap solitons, respectively.

In addition to the power-eigenvalue diagram, the power and linear stability of solitons are studied in Figure 6 for varied values of $\gamma ,\rho$ and $\beta$. The $x$-axis of each panel shows the existence domain for the varied parameter when other parameters are fixed to the values in Eq. (14). For instance, when $\mu =-0.1,\rho =0.5,\nu =1.5,\beta =2$ and $V_0=12.5$, fundamental solitons can be generated for $\gamma \in \left[-0.7,25.1\right]$ (see the left panel in Figure 6(a)). In Figure 6(a), stable (solid blue line) and unstable (red dotted line) regions are determined in the left panel, and in the right panel zoom-in view of the black rectangular region in the left panel is depicted.

Figure 6.

(a) Soliton power for varying quintic nonlinearity coefficient $\gamma$, left panel shows the domain of existence, and right panel shows zoom-in view of the black rectangular region. Soliton powers for $\gamma =-0.2$ (left panel) and $\gamma =0.2$ (right panel) (b) with varying $\rho$ and (c) with varying $\beta$ parameters. The linear stability and instability regions are shown by solid blue and red dotted lines, respectively.

Figure 6 shows that the power of solitons is decreased as $\gamma$ and $\beta$ are increased. The right panel in Figure 6(a) shows that solitons are linearly stable when $\gamma$ is between $-0.24$ and $0.26$, and comparing the left and right panels in Figure 6(b) and (c), it is seen that stability regions of $\rho$ and $\beta$ parameters for self-focusing ($\gamma =0.2$) case of quintic nonlinearity is narrow than that of self-defocusing ($\gamma =-0.2$) case.

To confirm obtained linear stability results, we study nonlinear evolution of solitons by means of direct simulation of the governing Eq. (2). To this end, using the fourth-order Runge–Kutta (RK4) scheme fundamental soliton is advanced in $z$ direction, and the spatial domain $\left(x,y\right)$ is discretized by the finite-difference method. It should be noted that the RK4 method is widely used for numerical analysis of nonlinear evolution equations, in this regard we apply RK4 to compare the results of this chapter with previous studies. The starting point of nonlinear evolution is chosen as the fundamental soliton, and $1\%$ random noise is inserted to amplitude of the soliton to test the nonlinear stability under perturbations.

In Figure 7, linear stability spectrum (first column), nonlinear evolution of the peak amplitude (second column), the evolution of soliton power (third column), and 3D view of the evolved soliton (fourth column) are plotted for the fundamental solitons that correspond to “a”, “b”, “c”, “d” and “e” points in the right panel of Figure 6(a). Here, the impact of quintic nonlinearity is examined by increasing $\gamma$ from $-0.3$ to $0.3$.

Figure 7.

Linear spectra (first column), nonlinear evolution of peak amplitudes (second column), the evolution of soliton power (third column), and 3D view of the soliton after evolution (fourth column) for (a) $\gamma =-0.3$; (b) $\gamma =-0.2$; (c) $\gamma =0$; (d) $\gamma =0.2$; and (e) $\gamma =0.3$. All other parameters are fixed to the values in Eq. (14).

In Figure 7(b)–(d), it is shown that the linear spectrum of the fundamental soliton is purely-imaginary (none of their eigenvalues have a real part), and the peak amplitudes oscillate during the propagation, thus stable evolution can be achieved for the considered solitons in the quadratic-cubic-quintic nonlinear medium when $\gamma =-0.2$, $\gamma =0$ and $\gamma =0.2$. On the other hand, when the coefficient of quintic nonlinearity is increased, there are eigenvalues with positive real part in linear stability spectra, and the peak amplitudes blow-up after a short propagation distance for $\gamma =-0.3$ and $\gamma =0.3$ as shown in Figure 7(a) and (e), respectively. This fact reveals that both strong self-focusing and self-defocusing quintic nonlinearities have an adverse effect on the stability of fundamental solitons.

From previous studies, it is known that increasing potential depth $V_0$ improves the stability of solitons in quadratic [23] and cubic (Kerr) medium [4]. In Figure 8, to examine the impact of the deeper lattice on the soliton stability in the quadratic-cubic-quintic nonlinear medium when $\gamma =-0.3$ (see panel (a)) and $\gamma =0.3$ (see panel (b)), the evolution of the peak amplitudes are displayed for varied potential depths. To compare the propagation distance, the amplitude of solitons are divided by their peak amplitudes. Figure 8 shows that larger lattice depth $V_0$ supports the stability of solitons for $\gamma =0.3$ (self-focusing), whereas shallow potential extends and eventually stabilizes the soliton for $\gamma =-0.3$ (self-defocusing). Thus, when the instability emerges from strong self-focusing nonlinearity, a deeper lattice can be utilized to arrest (or delay) the collapse. This result is consistent with the previous studies, and it is meaningful because the potential in our model (2) is defocusing and it must be balanced with a self-focusing term to obtain stable modes.

Figure 8.

Peak amplitude of the evolved solitons for varying lattice depths when the fundamental soliton is obtained for (a) $\gamma =-0.3$; and (b) $\gamma =0.3$.

Furthermore, it is known that increasing anisotropy parameter $\nu$ assists the stability of solitons when quadratic electro-optic effects are strong [23, 46, 47]. In order to investigate the impact of anisotropy on soliton dynamics, we plot the evolution of peak amplitude when $\nu$ is varying between 0 and 1000 in Figure 9. As shown in Figure 9, increasing the anisotropy coefficient $\nu$ (from 0.001 to 100) extends the propagation distance of the evolved solitons, and it stabilizes the soliton in the medium with a strong self-defocusing nonlinearity (when $\gamma =-0.3$) eventually (see black dotted line in panel (a)). Similarly, increased $\nu$ supports the stability of the solitons for strong self-focusing nonlinearity (when $\gamma =0.3$), and although it cannot stabilize the soliton robustly, it provides a significant extension of propagation distance when $\nu =1000$ (see black dotted line in panel (b)).

Figure 9.

Peak amplitude of the evolved solitons for varying anisotropy coefficients ($\nu$) between 0.001 and 1000 when the fundamental soliton is obtained for (a) $\gamma =-0.3$; and (b) $\gamma =0.3$.

It must be noted that even though soliton solutions of the NLSM system (2) can be obtained when $\nu \in \left(0,\infty \right)$ and stability of solitons are improved by higher values of $\nu$, they cannot be utilized in practical optical applications, due to the fact that $\rho$ and $\nu$ are predetermined constants associated with the type of optical materials.

4. Conclusions

Fundamental lattice solitons were obtained in a nonlocal nonlinear medium with self-focusing and self-defocusing quintic nonlinearity. The steady-state solutions are calculated by the SOM algorithm and the stability properties of solitons have been explored by linear spectra and nonlinear evolution of the amplitude by direct simulation of the nonlinear model. The band-gap boundaries and stability intervals of solitons were found and the power of solitons was investigated for varied $\rho ,\beta$ and $\gamma$ parameters. It has been seen that the power of fundamental solitons is decreased by increasing cubic ($\beta$) and quintic ($\gamma$) nonlinearity, while the power is increasing as the quadratic nonlinear response ($\rho$) is increased, and stability regions of $\rho$ and $\beta$ parameters for self-focusing ($\gamma =0.2$) case of quintic nonlinearity is narrow than that of self-defocusing ($\gamma =-0.2$) case.

The stability analysis showed that fundamental lattice solitons can stay stable for the weak self-focusing and self-defocusing quintic nonlinearity. Strong self-focusing ($\gamma =0.3$) and self-defocusing ($\gamma =-0.3$) quintic nonlinearities shortened the propagation distance of evolved solitons.

Furthermore, it has been observed that when the quintic nonlinearity is strong in the medium, solitons collapse, and the collapse of solitons can be arrested by a deeper lattice for self-focusing quintic nonlinearity ($\gamma =0.3$), whereas unstable solitons can be stabilized by shallow lattice for strong self-defocusing quintic nonlinearity ($\gamma =-0.3$).

In conclusion, the numerical existence and stability of fundamental lattice solitons have been presented in a nonlocal nonlinear medium with the quintic nonlinear response, and it has been demonstrated that stability of unstable solitons can be improved by modification of potential depth and strong anisotropy coefficient.

Conflict of interest

The authors declare no conflict of interest.