13383 FACTA UNIVERSITATIS Series: Electronics and Energetics Vol. 38, No 3, September 2025, pp. 553 - 562 https://doi.org/10.2298/FUEE2503553P © 2025 by University of Niš, Serbia | Creative Commons License: CC BY-NC-ND Original scientific paper SOLITARY AND TRAVELING WAVE SOLUTIONS TO NEMATIC LIQUID CRYSTAL EQUATIONS WITH CUBIC-QUINTIC NONLINEARITY USING THE JACOBI ELLIPTIC FUNCTION EXPANSION METHOD Nikola Z. Petrović Institute of Physics, Belgrade, Serbia Texas A&M University at Qatar, Doha, Qatar ORCID iD: Nikola Z. Petrović https://orcid.org/0000-0002-1297-3163 Abstract. In this paper, the Jacobi elliptic function (JEF) expansion method is applied to the system of equations governing nematic liquid crystals with a cubic-quintic nonlinearity. Solutions that are first order polynomials of the JEFs for the wave function and second order for the angle function are obtained. The solutions impose constraints on only two parameters and include a wide range of functions. Both solitary and traveling wave solutions are possible, as well as solutions both with and without chirp. Key words: Jacobi, nematic, liquid, crystal, photonics 1. INTRODUCTION Nonlinearities in optics are studied perhaps more than most other nonlinear system as there is a pressing need to support application in optical communications [1]. In particular, nonlinear behavior may be well controlled and defined by different kinds of optical materials such as nematic liquid crystals (NLCs) that have been recently produced and studied [1, 2]. Nematic liquid crystals are extremely versatile materials with a large range of practical uses in modern photonics [1]. They are an important system in nonlinear optics as they allow the study of many nonlinear phenomena at low power due to a very large nonlinear response via the light-induced reorientation of the NLC molecules [2], in particular the study of spatial solitons, which when propagating through NLCs are also known as nematicons [3]. The study and modeling of the behavior of nematicons, in particular finding the exact solutions describing their form, has numerous potential practical applications, such as optical information processing [4], molding of optical waveguides [5], beaming and control of the so-called random lasers [6] and many others [7,8]. Received January 12, 2025; revised January 28, 2025; accepted January 31, 2025 Corresponding author: Nikola Z. Petrović Institute of Physics, Belgrade, Serbia and Texas A&M University at Qatar, Doha, Qatar E-mail: nzpetr@ipb.ac.rs https://orcid.org/0000-0002-1297-3163 554 N. Z. PETROVIĆ NLCs are generally described by a pair of interconnected nonlinear differential equations describing the time evolution of the wave function of light and the angular function which describes the tilt of the molecules of the crystal [9]. There are several forms of nonlinearity which can occur in the second equation determining the angular function. The most common form of nonlinearity studied is the third-order nonlinearity, also known as the Kerr nonlinearity [9]. Several papers have produced solutions for the NLC system of equations (NLCSOE) with Kerr nonlinearity and basic solitary wave solutions have been obtained [10-12]. This paper will focus on the NLCSOE with the so- called cubic-quintic nonlinearity. Cubic quintic nonlinearity is a form of nonlinearity where the third and fifth order nonlinearities compete against each other [13]. It has emerged as an important topic of study in nonlinear optics due to the possibility of stabilizing solitary wave solutions with multiple transverse dimensions due to the competing signs of nonlinearities [14]. Several papers have used various techniques, such as the trial equation method [15], the sinh- Gordon expansion method [16] and others [9, 17-19] to find solutions for the NLCSOE with a cubic-quintic nonlinearity, often referred to in the papers as the parabolic law [17]. Recently, there has been a lot of progress in applying the JEF expansion method to find solutions to the Nonlinear Schrodinger equation with various forms of nonlinearity [20-22], as well as the Gross-Pitaevskii equation [23-24]. The method has also successfully been applied to two-component systems such as the Davey-Stewartson equation [25] and the two-component NLSE [26]. The first application of the JEF expansion method on NLCs was made in [27] where solutions were found for the NLC system of equations with a third-order nonlinearity. In this work, we generalize the Jacobi elliptic function (JEF) expansion method that was developed in [22] and [27] to find exact solutions to the NLC system of equations (NLCSOE) for the cubic-quintic (CQ) nonlinearity. As in [27], we apply the principle of harmonic balance to both the wave function and the angular tilt and apply matching conditions to obtain the polynomial degrees of these two functions in terms of the JEF. These degrees will depend on the degree of the nonlinearity inside the liquid crystal and it turns out will differ from the degrees obtained in [27]. 2. METHOD The NLCSOE for the CQ nonlinearity has the general form as follows [9]: 0, 2 t xxiu u pu  + + = (1) 2 4 1 2 0xxcp lp u u + + + = (2) where u is the wave function, p is the angle function determined by the orientation of NLCs, β is the diffraction parameter, χ is the coupling parameter, c and l are parameters describing the strength of the non-local response of the NLCs and α1 and α2 are parameters determining the strength of the nonlinear response to the propagating light. In the special case where the parameter c is equal to 0, the system of equations reduces to the standard cubic-quintic NLSE. Solitary and Traveling Wave Solutions to Nematic Liquid Crystal Equations with Cubic-Quintic... 555 Following [22], the function u is split into the real and imaginary parts: iBu Ae= (3) where A is the amplitude and B is the phase of the solution. Plugging in the equations and splitting the real and imaginary parts we obtain: (2 ) 0 2 t x x xxA A B AB  + + = (4) 2(2 ) 0 2 t xx xAB A AB AP  − + + + = (5) 2 4 1 2 0xxcp lp A A + + + = (6) We now assume the following forms for A and B: 1 1 1( ) ( ) ( ) ( )A f t F f t F − −= + (7) ( ) ( )k t x t = + (8) 2( ) ( ) ( ).B a t x b t x e t= + + (9) where F is a Jacobi elliptic function satisfying the following differential equations: 2 2 4 3 0 2 4 2 42 and 2 dF d F c c F c F c F c F d d  = + + = + (10) where c0, c2 and c4 are coefficients that depend on the choice of the Jacobi elliptic function and the so-called JEF parameter M. For F=dn we have c0=M-1, c2=2-M and c4=-1, while for F=sn we have c0=1, c2=-(1+M), c4=M. The remaining parameters f1, f-1, k, ω, a, b and e are functions of time to be determined. We note that the phase contains the quadratic term a with respect to the transverse variable that is known as the chirp [20]. We now apply the matching principle to find the needed degree of F in p. Since in Eq. (7), the highest degree of F is 3 in the term , the matching conditions indicate that the degree of AP should also be 3 and, therefore, the angle function p should be a second order function of F: 2 2 2 0 2( ) ( ) ( )p g t F g t g t F − −= + + (11) The terms of odd degree are omitted because they add too many new equations without any benefit. It is worth noting that for the ordinary Kerr nonlinearity the matching conditions imposed second degree functions in F for both A and p [27]. We now plug Eqs. (7-9) and Eq. (11) into Eqs. (4-6) to obtain a polynomial function of F. Taking care to equate each coefficient of the polynomial to 0, we obtain a series of algebraic and ordinary differential equations: 0, 1, 1it if a f i+ = = − (12) 22 0,ta a + = (13) 2 0,tb ab+ = (14) 556 N. Z. PETROVIĆ 2 0,tk ak+ = (15) 0.t bk + = (16) For the parameter e, we obtain a pair of equations that will have to be equivalent, i.e. matched, for the solution to be valid: 2 2 0 2 2 1 1 e 0, 1, 1. 2 2 t i i i i i if b f f g f g c f k i   −− − + + + = = − (17) Finally, we obtain several additional constraits between parameters which can be thought of as integrability conditions: 2 2 2 0, 1, 1,i i i if g c f k i  ++ = = − (18) 4 2 2 2 26 0, 1, 1,i i if cc g k i ++ = = − (19) 2 3 2 1 2 2 2 24 4 0 1, 1,i i i i if f f cc g k g l i  −+ + + = = − (20) 2 2 2 2 1 1 1 2 1 1 0 2 4 2 02 6 2 2 0.f f f f cc g k cc g k g l − − −+ + + + = (21) We now proceed to solve Eqs (12-21). Solutions to Eqs. (12-16) are obtained using standard techniques and are as follows: 1 2 0 ,i if f = (22) 0 ,a a = (23) 0 ,b b = (24) 0 ,k k = (25) 0 0 0 0 , t b k dt   = −  (26) where 0 0 1 1 2 t a dt   = +  is the so-called chirp function [17]. In the absence of chirp, i.e. a0 = 0, we have η = 1. Without loss of generality, we can now assume f1 ≠ 0. For f-1 = 0, we obtain g-2 = 0 from Eq. (17) for i = -1. From solving Eqs. (18)-(21), we obtain: 2 4 2 , c k g   − = (27) 4 4 2 4 1 6 , cc k f    = (28) 2 24 1 22 1 (4 ), c k cc k l f    = + (29) 4 0 4 0 2 . l cc c k g   = (30) Solitary and Traveling Wave Solutions to Nematic Liquid Crystal Equations with Cubic-Quintic... 557 For f-1 ≠ 0, from matching the two equations for e in Eq. (17), we obtain: 1 2 1 2 1 1 f g f g f f − − − = (31) and therefore since: 2 4 2 , c k g   − = (32) 2 0 2 , c k g   − − = (33) we obtain: 01 1 4 , 1. cf f c − = =  (34) The formula for α2 is the same as in Eq. (28). The remaining formulas are: ( ) 2 24 1 2 0 42 1 4 ( 6 ) , c k ck c c c l f    = − + (35) 2 2 0 2 0 2 4 4 4 0 0 2 (1 7 ) 4 . l cc c c k g c c l c k c      = + − −        (36) As can be seen, the solutions impose constraints on only two parameters, and , while the remaining parameters β, χ, c and l, are completely arbitrary. This allows for a wide range of flexibility in constructing our solutions. Finally, the formulas for e in both cases will be complicated and dependent on the form of β, χ, c and l chosen. 3. RESULTS We now present the solutions we obtained with this method. We will first select F=dn for our Jacobi elliptic function. This function is convenient because the reciprocal function F=nd doesn’t contain singularities, thus allowing us to obtain novel nonsingular solutions for non-zero ϵ. In Fig. 1a we see a standard bright solitary wave solution. The position can be altered by changing ω0 and the extent of oscillations can be controlled by changing b0. We see in Fig. 1c that the NLC acts as a wave guide for the signal. In Figs 1b and d we see the effects of the chirp function. The chirp function will deform the solution in the transverse direction and introduce oscillations in the amplitude. Solitary waves with such oscillations in amplitude are often called breathers. Since c0=0 for M=1, from Eq. (34) we have f-1=0 and therefore we do not have any solutions for M=1 which combine F=dn and F=nd. 558 N. Z. PETROVIĆ Fig. 1 Bright solitary wave solutions to the NLCSOE as a function of k0x and t for F=dn, β(t) = β0cos(Ωt) and M=1. Graphs (a), (b) depict the square of the angle function |u|2 and graphs (c), (d) depict the angle function p. The values of the parameters are: β0 = Ω = k0 = b0 = f10 = l = c = χ = 1, ω0 = e0 = ϵ = 0 and (a),(c): a0 = 0, (b),(d): a0 = 0.2. Fig. 2 Traveling wave solutions to the NLCSOE as a function of k0x and t for M=0.9, a0 = 0. Graphs (a), (b), (c) depict the square of the angle function |u|2 and graphs (d), (e), (f) depict the angle function p. We have for (a),(d): ϵ = 0, , (b),(e): ϵ=1, (c),(f): ϵ=-1. All the other parameters are the same as in Fig. 1. Solitary and Traveling Wave Solutions to Nematic Liquid Crystal Equations with Cubic-Quintic... 559 In Fig. 2, we see the periodic, so-called traveling wave, solutions to the NLCSOE. For M<1, the JEF no longer produces a solitary wave but a periodic wave structure. We see that both the wave and the angle functions (Figs 2a and 2d) become periodic in the transverse direction. In Figs 2b, c, e and f we see the effects of a non-zero value of ϵ. We see that the overall effect of combining F=dn and F=nd is to double the periodicity of the solutions. The variation of the angle functions also becomes more prominent in the longitudinal direction. Solutions in Fig 2b and Fig 2c are qualitatively alike except for the shift in the overall background amplitude due to the sign of ϵ. The forms of the angle function p are, however, far more complicated and the two solutions in Fig 2e and Fig 2f are quite different from each other. Fig. 3 Traveling wave solutions to the NLCSOE with chirp as a function of k0x and t for M=0.9, a0 = 0.2. All the other parameters are the same as in Fig. 2. In Fig. 3 we see the effects of chirp on the traveling wave solutions. The wave fronts in Figs 3a, b and c are stretched out in the transverse direction and no longer periodic. The more one deviates from an equilibrium point which is near the axis, the more extreme the stretching of the wave front. We also note that the solutions in Fig 3b and Fig 3c are no longer qualitatively alike due to the disruption of symmetry caused by the chirp function. In Fig 3d we see the effect of chirp on the angle function. One can clearly see the orientation of the wave crests change with the change in the transverse variable. This is less noticeable in Figs 3e and f where there is a strong background component coming from g0. Finally, in Fig 4, we see the dark soliton solutions to the NLCSOE where we have used the JEF F=sn. The standard dark solitary wave solution is shown in Fig 4a. We note the periodic structure of the angle function in Fig 4b in the presence of the background with almost a small deviation from it that produces the dark solitary wave. In Fig 4b we see that the background of the solution is affected by the chirp. There is a similar effect on the angle function in Fig 4d of pushing the wave to one side as in Figs 3d, e and f. Lastly, we see an example of the traveling wave solution for the dark solitary wave in 560 N. Z. PETROVIĆ Figs 4c and f. We see, unlike F=dn in Fig 2a, that F=sn reaches 0 and therefore including non-zero ϵ would produce singularities. We see the structure in Fig 2d that produced a dark solitary wave now repeated periodically in Fig 2f. Fig. 4 Dark solitary wave solutions to the NLCSOE as a function of k0x and t for F=sn. All the other parameters for (a), (b), (d), (e) are respectively the same as in Fig. 1(a), (b), (c) and (d). All the other parameters for (c) and (f) are respectively the same as in Fig 2(a) and (d). 4. DISCUSSION We now compare the presently obtained results with the previous results obtained in other papers. In [9], the authors apply a simple ansatz to obtain the basic A=sech solution for the amplitude corresponding to the bright solitary wave for M=1. In [16], the sinh- Gordon expansion method is applied and various forms of hyperbolic trigonometric functions are obtained for the wave function, although it has to be mentioned that the form for the angle function is not necessarily related to the wave function. In [17], the exp(-φ) method is applied and a couple of solutions are obtained, usually in the form of a fraction with a complicated denominator. In [18], the so-called simple equation method is applied. Solutions for the amplitude of the form tanh and coth are obtained as well as various ratios of exponential functions. In [19], two solutions are obtained based on the F=tan and the F=tanh function. In [28], the Lie point symmetry method is applied to the dual-power law nonlinearity which reduces to the cubic-quintic nonlinearity for n=1 and several solutions related to the csch, sec and cos functions are obtained. In [29], the W-shaped solutions, which occur in the case of Kerr nonlinearity [27], are studied using the general exponential rational function method and some solutions for the case of the parabolic nonlinearity are obtained. The exponential rational functions that are used to construct both the wave and angle functions are ratios of either exponential or trigonometric functions. In [30], the Solitary and Traveling Wave Solutions to Nematic Liquid Crystal Equations with Cubic-Quintic... 561 Kudryashov’s approach and the tanh–coth technique are used to construct solutions for various nonlinearities, including parabolic nonlinearity. Finally, in [15] large classes of solutions are obtained although typically they also involve complicated expressions in the denominator, usually some form of a function plus constant terms. Several solutions involving hyperbolic trigonometric functions are obtained and there is even one mention of the Jacobi sn function squared, albeit in one such denominator. There are several advantages to the methods in this paper. First, the method is conceptually simple and doesn’t require many complicated parameters like in the other approaches. Second, the method works for arbitrary functions of β, χ, c, l, α1 and α2 with respect to time, whereas the previous papers treat these functions as constant parameters. This allows greater flexibility in the study of various NLC systems, especially those that employ dispersion management. We note that our method is completely applicable to the case where β, χ, c, l, α1 and α2 are arbitrary constants in the case of no chirp. The only difference in this case is that constraints are then imposed on k and f1 via Eqs (28)-(29). Third, the JEFs are extremely flexible functions, containing both solitary and periodic waves. By varying the choice of the JEF and the parameter M, many different qualitative forms of solutions can be obtained. Finally, no previous paper covers solutions with chirp, which is an important phenomenon in understanding pulse propagation. Here we will briefly discuss the limitations of the method. First, the method is limited by the forms of functions that satisfy Eq. (10). The functions satisfying (10) either have a solitary wave or are periodic. Modeling multiple, but finite, number of waves is difficult with this method. Second, as mentioned before, the method is not applicable for a completely arbitrary set of parameters, but does have two constraints. Finally, the form of the free parameter e in the phase is complicated and needs to be calculated for each set of functions individually because the form of the differential equations in (17) will greatly differ based on the forms of β, χ, c and l. Nevertheless, for most practical applications only the amplitude of u is needed. 5. CONCLUSION We applied the Jacobi elliptic function expansion method to the NLC system of equations with a cubic quintic nonlinearity and obtained abundant classes of solutions to the system. Both solitary and traveling wave solutions were obtained, as well as solutions that contain chirp. In particular, the second derivative of the angle function allowed the wave function to have a better match with the degree of nonlinearity in the system. In addition, the fact that there are only two constraints in the system, allow systems of NLCs to be flexibly tuned to allow the propagation of the wave function through them. This could potentially have many applications in the fields of photonics and nonlinear optics. There are many potential systems to which this method can further be applied, including two-component NLC systems and NLC systems with different forms of nonlinearity, especially the so-called septic (seventh order) nonlinearity. Acknowledgement: Work at the Institute of Physics is supported by the project OI 171006 of the Serbian Ministry of Education and Science. Work at the Texas A&M University in Qatar (now the Hamad Bin Khalifa University College of Science and Engineering) was supported by the project NPRP13S-0121-200126 with the Qatar National Research Fund (a member of the Qatar Foundation.) 562 N. Z. PETROVIĆ REFERENCES [1] A. d’Alessandro and R. Asquini, "Light Propagation in Confined Nematic Liquid Crystals and Device Applications", Appl. Sci., vol. 11, no. 18, p. 8713, 2021. [2] G. Assanto, A. Fratalocchi and M. Peccianti, "Spatial Solitons in Nematic Liquid Crystals: From Bulk to Discrete", Opt. Express, vol. 15, no. 8, pp. 5248-5259, 2007. [3] G. Assanto, M. Peccianti and C. Conti, "Nematicons: Optical Spatial Solitons in Nematic Liquid Crystals," Opt. Photonics News, vol. 14, no. 2, pp. 44-48, 2003. [4] G. Assanto, "Nematicons: Reorientational Solitons from Optics to Photonics", Liq. Cryst. Rev. vol. 6, no. 2, pp. 170-194, 2018. [5] N. Karim, et al., "Molding Optical Waveguides with Nematicons", Adv. Opt. Mater., vol. 5, no. 14, p. 1700199, 2017. [6] S. Perumbilavil, et al., "Beaming Random Lasers with Soliton Control", Nat. Commun., vol. 9, no. 1, p. 3863, 2018. [7] S. Zhang, et al., "Deflecting and Routing Nematicons via Orientation Programmable Liquid Crystal Array", Opt. Express, vol. 31, no. 13, pp. 20955-20964, 2023. [8] U. A. Laudyn, M. Kwasny and M. A. Karpierz, "Nematicons in Chiral Nematic Liquid Crystals", Appl. Phys. Lett., vol. 94, no. 9, p. 091110, 2009. [9] M. Savescu, et al., "Nematicons in Liquid Crystals", J. Comput. Theor. Nanosci., vol. 12, no. 11, pp. 4667-4673, 2015. [10] L. Kavitha, et al. "Optically Induced Switching of Nematic Deformations", Physica Scripta, vol. 88, no. 6, p. 065015, 2013. [11] S. Duran and B. Karabulut, "Nematicons in Liquid Crystals with Kerr Law by Sub-Equation Method", Alexandria Eng. J., vol. 61, no. 2, pp. 1695-1700, 2022. [12] J. M. L. MacNeil, N. F. Smyth and G. Assanto, "Exact and Approximate Solutions for Optical Solitary Waves in Nematic Liquid Crystals", Physica D, vol. 284, pp. 1-15, 2014. [13] A. Desyatnikov, A. Maimistov and B. Malomed, "Three-Dimensional Spinning Solitons in Dispersive Media with The Cubic-Quintic Nonlinearity," Phys. Rev. E, vol. 61, no. 3, p. 3107, 2000. [14] F. Ndzana II, A. Mohamadou and T. C. Kofané, "Modulational Instability in The Cubic–Quintic Nonlinear Schrödinger Equation Through the Variational Approach", Opt. Commun., vol. 275, no. 2, pp. 421-428, 2007. [15] M. Ekici, et. al., "Nematicons in Liquid Crystals by Extended Trial Equation Method", J. Nonlinear Opt. Phys. Mater., vol. 26, no. 1, p. 1750005, 2017. [16] D. Kumar, et. al., "Investigation of Dynamics of Nematicons in Liquid Crystals by Extended Sinh‑Gordon Equation Expansion Method", Opt. Quantum Electron., vol. 51, p. 212, 2019. [17] N. Raza, et al., "Optical Solitons in Nematic Liquid Crystals with Kerr And Parabolic Law Nonlinearities", Opt. Quantum Electron., vol. 51, p. 107, 2019. [18] A. H. Arnous, et. al., "Nematicons in Liquid Crystals by Modified Simple Equation Method", Nonlinear Dyn., vol. 88, p. 2863-2872, 2017. [19] O. A. Ilhan, et al., "New Exact Solutions for Nematicons in Liquid Crystals by The Tan (Φ/2) -Expansion Method Arising in Fluid Mechanics", Eur. Phys. J. Plus, vol. 135, no. 3, pp. 1-19, 2020. [20] V. I. Kruglov, A. C. Peacock and J. D. Harvey, "Exact Self-Similar Solutions of the Generalized Nonlinear Schrödinger Equation with Distributed Coefficients", Phys. Rev. Lett., vol. 90, p. 113902, 2003. [21] W. Zhong, et al., "Exact Spatial Soliton Solutions of the Two-Dimensional Generalized Nonlinear Schrödinger Equation with Distributed Coefficients", Phys. Rev. A, vol. 78, no. 2, p. 023821, 2008. [22] M. Belić, et al., "Analytical Light Bullet Solutions to the Generalized (3+ 1)-Dimensional Nonlinear Schrödinger Equation," Phys. Rev. Lett., vol. 101, no. 12, p.123904, 2008. [23] N. Z. Petrović, M. Belić, and W. Zhong, "Spatiotemporal Wave and Soliton Solutions to the Generalized (3+ 1)- dimensional Gross-Pitaevskii Equation", Phys. Rev. E, vol. 81, no. 1, p. 016610, 2010. [24] N. Z. Petrović, et al., "Analytical Traveling-Wave and Solitary Solutions to the Generalized Gross-Pitaevskii Equation with Sinusoidal Time-Varying Diffraction and Potential", Phys. Rev. E, vol. 83, no. 3, p. 036609, 2011. [25] N. Z. Petrović, "Solitary and Traveling Wave Solutions for the Davey–Stewartson Equation Using the Jacobi Elliptic Function Expansion Method", Opt. Quantum Electron., vol. 52, no. 6, p. 319, 2020. [26] N. Petrović and H. Zahreddine, "Exact Traveling Wave Solutions to Coupled Generalized Nonlinear Schrödinger Equations", Physica Scripta, vol. 2012, no. T149, p. 014039, 2012. [27] N. Petrović, M. Belić and W. Krolikowski, "Solitary and Traveling Wave Solutions to Nematic Liquid Crystal Equations Using Jacobi Elliptic Functions", Chaos, Solitons & Fractals: X, vol. 13, p. 100121, 2024. [28] R. Ramaswamy, A. H. Abdel Kader and A. Elsonbaty, "New Exact Nematicon Solutions of Liquid Crystal Model with Different Types of Nonlinearities", IEEE Access, vol. 9, pp. 107909-107916, 2021. [29] H. F. Ismael, H. Bulut and H. M. Baskonus, "W-shaped Surfaces to the Nematic Liquid Crystals with Three Nonlinearity Laws", Soft. Comput., vol. 25, no. 6, pp. 4513-4524, 2021. [30] N. Raza, et al., "Investigation of New Solitons in Nematic Liquid Crystals with Kerr and Non-Kerr Law Nonlinearities", J. Nonlinear Opt. Phys. Mater., vol. 32, no. 02. p. 2350020, 2023.