mathematics Article Integrability Properties of the Slepyan–Palmov Model Arising in the Slepyan–Palmov Medium Muhammad Usman 1 , Akhtar Hussain 2 , F. D. Zaman 2 , Asier Ibeas 3, * and Yahya Almalki 4 1 College of Electrical and Mechanical Engineering (CEME), National University of Sciences and Technology (NUST), H-12, Islamabad 44000, Pakistan; musman.bsh23ceme@student.nust.pk 2 Abdus Salam School of Mathematical Sciences, Government College University, Lahore 54600, Pakistan; akhtarhussain21@sms.edu.pk (A.H.); f.zaman@sms.edu.pk (F.D.Z.) 3 Department of Telecommunications and Systems Engineering, Universitat Autònoma de Barcelona, 08193 Barcelona, Spain 4 Department of Mathematics, College of Science, King Khalid University, Abha 61413, Saudi Arabia; yalmalki@kku.edu.sa * Correspondence: asier.Ibeas@uab.cat Abstract: This study investigates the Slepyan–Palmov (SP) model, which describes plane longitudinal waves propagating within a medium comprising a carrier medium and nonlinear oscillators. The primary objective is to analyze the integrability properties of this model. The research entails two key aspects. Firstly, the study explores the group invariant solution by utilizing reductions in symmetry subalgebras based on the optimal system. Secondly, the conservation laws are studied using the homotopy operator, which offers advantages over the conventional multiplier approach, especially when arbitrary functions are absent from both the equation and characteristics. This method proves advantageous in handling complex multipliers and yields significant outcomes. Keywords: Lie symmetry method; the Slepyan–Palmov model; optimal system; longitudinal wave; conservation laws; homotopy operator MSC: 70S10; 70G65; 35B06; 37K06; 17B81; 35R03 Citation: Usman, M.; Hussain, A.; Zaman, F.D.; Ibeas, A.; Almalki, Y. Integrability Properties of the 1. Introduction Slepyan–Palmov Model Arising in The Cauchy stress principle is a fundamental assumption in classical continuum the Slepyan–Palmov Medium. mechanics which establishes a key equivalence between the action of all internal forces Mathematics 2023, 11, 4545. https:// applied to an elementary area and the action of their resultant force applied at the center of doi.org/10.3390/math11214545 the area. This principle forms the backbone of classical continuum mechanics and plays Academic Editor: Alain Rapaport a crucial role in understanding the behavior of continuous materials under various loads and deformations. However, as engineering and scientific investigations have advanced, Received: 8 October 2023 Revised: 27 October 2023 it has become evident that the simple stress principle embodied by Cauchy’s hypothesis Accepted: 2 November 2023 might not fully capture the complexity of certain materials and deformation scenarios. In Published: 3 November 2023 more general cases, the action of an arbitrary system of forces on a deformable solid cannot be solely represented by stresses alone. In addition to stresses, moment stresses emerge, leading to the formation of asymmetric tensors. These moment stresses significantly influence the mechanical response of materials, especially in situations where classical Copyright: © 2023 by the authors. continuum mechanics fall short of describing real-world phenomena accurately. Licensee MDPI, Basel, Switzerland. To adequately account for these complex behaviors and to broaden the understanding This article is an open access article of material responses, the introduction of additional degrees of freedom in the medium distributed under the terms and becomes necessary. This involves considering a physically infinitesimal volume (over conditions of the Creative Commons which the properties of the medium are averaged) not as a simple material point, but Attribution (CC BY) license (https:// as a more intricate object with new degrees of freedom. This recognition of additional creativecommons.org/licenses/by/ degrees of freedom allows for the inclusion of a microstructure within the material, such as 4.0/). Mathematics 2023, 11, 4545. https://doi.org/10.3390/math11214545 https://www.mdpi.com/journal/mathematics Mathematics 2023, 11, 4545 2 of 13 graininess, fibrousness, or cellular structures present in real-world materials. By embracing this expanded perspective, the range of properties that can be modeled within a continuous medium expands considerably. The presence of internal microstructure and the considera- tion of additional degrees of freedom enable a more comprehensive representation of the mechanical response of materials, particularly in cases where classical continuum mechan- ics might be limited in its accuracy. Both the classical continuum model and the generalized continuum model are extensively utilized in modern deformable solid mechanics [1–3]. In addition to the classical continuum model and the generalized continuum model, there are other models used in deformable solid mechanics. One example is the Cosserat continuum [4], which considers internal rotational degrees of freedom. This theory was originally published in Russian and it has been subject to comments and discussions [5]. Gradient models, like the Leroux continuum [6], Jeremillo model, and Tupin model [7], also exist. These models incorporate gradients of certain fields to handle materials with spatial variations in properties, proving useful in describing complex materials. Furthermore, there are models that account for media with oscillatory degrees of freedom. These models address the presence of oscillations or vibrations within the material. The works of Slepyan [8] in one-dimensional systems and Palmov [9] in three-dimensional systems have significantly contributed to this area of research. These diverse models enrich our understanding of material behaviors and offer valuable insights into various physical phenomena in deformable solids. The SP model combines a linearly elastic carrier medium with non-interacting os- cillators (elastic or viscoelastic) suspended at each point. It postulates that the dynamic behavior of the carrier model is described by the Lame equations and the oscillators fixed within it have continuously distributed eigenfrequencies. The following equations govern the dynamics of the SP model Z ∞ ($ + µ) grad div v + µ∆v − ρv̈ − m(q)ẅq dq + K + Q = 0,   0 (1) ∂  m(q)ẅq + c(q) 1 + R̃ wq − v = Qq . ∂κ where • ρ represents the mass density of the carrier medium. • $ and µ are the Lame elastic moduli characterizing the carrier medium’s elasticity. • v denotes the displacement vector of points within the carrier medium. • K stands for the intensity of the external body force acting on the medium. • wq represents the absolute displacement vector of the oscillator mass with respect to its equilibrium position. • Qq denotes the external force applied to the mass of the oscillator. • The quantity m(q)dq corresponds to the mass of all oscillators with eigenfrequencies lyingRwithin the interval (q, q + dq) multiplied by a unit volume. ∞ • m = 0 m(q)dq is the total mass density of all oscillators fixed to the carrier medium. • c(q) = q2 m(q) represents the static stiffness of the oscillator suspension. The quantity R̃(∂/∂κ ) characterizes the energy dissipation in the oscillator suspension. An interesting characteristic of this model is that even with the low damping of oscillators, the spatial attention of vibrations in the medium is finite [10]. As a result, the model finds effective application in calculating the vibrations of aircraft, rockets, and space technology objects, as well as submarines. It proves to be a valuable tool in analyzing and predicting the behavior of these structures under various vibrational conditions. In the context of an isolated environment with K = 0 and Qq = 0, we focus on a one- dimensional version of the system (1). To account for nonlinearity attributed to the carrier Mathematics 2023, 11, 4545 3 of 13 medium, we introduce nonlinear terms into the system. The resulting one-dimensional system is expressed as "   # ∂2 vy (y, κ ) ∂ ∂vy 2 ∂2 v y Z ∞ ∂2 wqy (y, κ ) ($ + 2µ) + F − ρ − m ( q ) dq = 0, ∂y2 ∂y ∂y ∂κ 2 0 ∂κ 2 (2) ∂2 wqy    ∂  m(q) + c(q) 1 + R̃ wqy − vy = 0, ∂κ 2 ∂κ where F is the coefficient characterizing the nonlinearity of the carrier medium. Under the assumption that the absolute displacement  of the oscillator does not depend on its eigenfrequency, we can express it as wqy = wy Z ∞ ∂2 wqy ∂2 wy (y, κ ) Z ∞ m(q) dq = m(q)dq. 0 ∂κ 2 ∂κ 2 0 The system (2) is simplified to a single equation representing the longitudinal displace- ment of the carrier medium, denoted as vy given by ∂2 v y 2 2 ∂ vy c21 ∂4 vy 1 ∂4 v y ∂3 v y 3 2 ∂ vy − c 1 − + + R̃ − R̃c 1 ∂κ 2 ∂y2 q2 ∂y2 ∂κ 2 α2 q2 ∂κ 4 ∂κ 3 ∂y2 ∂κ " # " # "  # (3) c22 ∂3 ∂vy 2 ∂vy 2 2 ∂vy 2   2 ∂ 2 ∂ − 2 − c2 − R̃c2 = 0, q ∂y∂κ 2 ∂y ∂y ∂y ∂y∂κ ∂y with Z ∞ $ + 2ν F ρ̃ ρ̃ = ρ + m(q)dq, c21 = , c22 = , α2 = · 0 ρ̃ ρ̃ ρ We define the dimensionless parameters as vx αq θ= , x= y, t = αqκ. v0 c1 By introducing v0 as the maximum displacement, within which the deformation of the carrier medium remains elastic, Equation (3) is modified to take the following form ∂2 θ ∂2 θ ∂2 ∂2 θ 2  ∂ ∂2 θ ∂2 θ    2∂ θ − 2+ 2 −α + R̃αq − 2 ∂t2 ∂x ∂t ∂t2 ∂x2 ∂t ∂t2 ∂x "  # "  # "  # 2 2 ∂ ∂θ ∂ 3 ∂θ ∂ 2 ∂θ 2 − Ñ − Ñα2 2 − R̃ Ñαq = 0, ∂x ∂x ∂t ∂x ∂x ∂t∂x ∂x where Ñ = c22 αqv0 /c31 ∂2 θ ∂2 θ ∂2 ∂2 θ 2  ∂ ∂2 θ ∂2 θ    2∂ θ or − + − α + R − ∂t2 ∂x2 ∂t2 ∂t2 ∂x2 ∂t ∂t2 ∂x2 "  # "  # "  # ∂ ∂θ 2 ∂2 ∂θ 2 ∂3 ∂θ 2 − N1 − RN1 −β 2 = 0, ∂x ∂x ∂t∂x ∂x ∂t ∂x ∂x where β/N1 = α2 > 1. Considering the nonlinearity of the medium in the absence of dissipation (R = 0, N1 = 0), Equation (3) can be expressed as follows "  # ∂2 θ ∂2 θ ∂2 ∂2 θ 2  ∂3 ∂θ 2  2∂ θ 2 − 2+ 2 2 −α 2 −β 2 = 0. ∂t ∂x ∂t ∂t ∂x ∂t ∂x ∂x Mathematics 2023, 11, 4545 4 of 13 In this study, we conduct a Lie symmetry analysis of the SP model [11] in the given mathematical form θtt − θ xx + (θtt − α2 θ xx ) xx − β((θ x )2 )ttx = 0. (4) Lie theorists like Ovsyannikov [12], Ibragimov [13], Bluman [14], Olver [15], Hydon [16], Stephani [17], and others [18,19] have played a vital role in connecting group structures to differential equations. This link has been crucial in discovering invariant solutions for nonlinear partial differential equations and gaining deeper insight into their behavior. A significant area of study has been the group analysis of the SP model, where Lie theory techniques have led to the identification of invariant solutions and conservation laws, highlighting essential underlying symmetries and dynamic properties. This study provides group invariant solutions that remain unchanged under a specific Lie symmetry generator, adding novelty to our research as such solutions are not present in the existing literature. Additionally, the conservation laws outlined in this study are reported for the first time, contributing to the enrichment of the SP model’s applications from a divergence perspective. Our study is specifically centered on invariant solutions and nonlocal con- servation laws, which represents a limitation in terms of the methods we have employed. References [20–22] represent valuable contributions to this field, shaping our understanding of the SP model and its applications in various scientific and engineering domains. The article is organized as follows: Section 1 presents a detailed explanation of the SP model. In Section 2, the Lie group method and one-dimensional subalgebras for the SP model are discussed. Section 3 is devoted to exploring the invariant solutions of the SP model. Moving on to Section 4, the article examines the derivation of conservation laws using the homotopy operator. Section 5 focuses on the physical interpretation of the solutions. Lastly, in Section 6, the article concludes by highlighting potential future research directions. 2. Lie Group Method In this section, we investigate the Lie symmetries and optimal system of Equation (4). We consider a one-parameter Lie group of transformations to identify the inherent symme- tries in the equation x̃ → x + ςϕ1 ( x, t, θ ) + O(ς2 ), t̃ → t + ςϕ2 ( x, t, θ ) + O(ς2 ), (5) 2 θ̃ → θ + ς$( x, t, θ ) + O(ς ), where ς is the group parameter. The vector field associated with the above transforma- tions is ∂ ∂ ∂ E = ϕ1 ( x, t, θ ) + ϕ2 ( x, t, θ ) + $( x, t, θ ) · (6) ∂x ∂t ∂θ The coefficient functions ϕ1 , ϕ2 , and $ are to be found, and the operator E fulfills the Lie symmetry condition [23] E [4] (θtt − θ xx + (θtt − α2 θ xx ) xx − β((θ x )2 )ttx )|(4) = 0, (7) where E [4] is the fourth extension of E . This leads to the four-dimensional symmetry algebra for Equation (4) given by (Table 1): ∂ ∂ ∂ ∂ E1 = , E2 = , E3 = , E4 = t · (8) ∂t ∂θ ∂x ∂θ The adjoint representation is given by (Table 2): e2 Ad(exp (eEi ).E j ) = E j − e[Ei , E j ] + [E , [E , E ]] − · · · (9) 2! i i j Mathematics 2023, 11, 4545 5 of 13 Table 1. Commutator table. [Ei , E j ] E1 E2 E3 E4 E1 0 0 0 E2 E2 0 0 0 0 E3 0 0 0 0 E4 −E2 0 0 0 Table 2. Adjoint table. Ad(ee ) E1 E2 E3 E4 E1 E1 E2 E3 E4 − e E2 E2 E1 E2 E3 E4 E3 E1 E2 E3 E4 E4 E1 + e E2 E2 E3 E4 Proposition 1. Let L4 be the Lie algebra of the SP model (4). The optimal system of one-dimensional subalgebras is then generated by the generators listed below J1 = hE4 i, J2 = hE1 + cE4 i, c 6= 0, J3 = hE3 + cE4 i, c 6= 0, J4 = hE1 + cE3 + dE4 i, c, d 6= 0, J5 = hE2 i, (10) J6 = hE1 i, J7 = hE3 i, J8 = hE2 + cE3 i, c 6= 0, J9 = hE1 + cE3 i, c 6= 0. Proof. Take any element E ∈ L4 . We have, E = µ1 E1 + µ2 E2 + µ3 E3 + µ4 E4 (11) µ4 µ4 6 = 0 µ4 = 0 µ3 = 0 µ3 6 = 0 µ3 = 0 µ3 6 = 0 µ1 = 0 µ1 6 = 0 µ1 = 0 µ1 6 = 0 µ1 = 0 µ1 6 = 0 µ1 = 0 µ1 6 = 0 µ2 = 0 µ2 6 = 0 Case 1 Case 2 Case 3 Case 4 Case 5 Case 6 Case 9 Case 7 Case 8 Case 1: µ4 6= 0, µ3 = 0, µ1 = 0. Then we have, E = µ2 E2 + µ4 E4 (12) Mathematics 2023, 11, 4545 6 of 13 Under the adjoint action on E , we have E 0 = Ad(ee E1 )E = µ4 E4 (13) Hence, we obtain J1 = E4 , (14) Case 2: µ4 6= 0, µ3 = 0, µ1 6= 0. Then we have, E = µ1 E1 + µ2 E2 + µ4 E4 (15) E 0 = Ad(ee E4 )E = µ1 E1 + µ4 E4 (16) Hence, we obtain J2 = E1 + cE4 , c 6= 0, (17) Case 3: µ4 6= 0, µ3 6= 0, µ1 = 0. Then we have, E = µ2 E2 + µ3 E3 + µ4 E4 (18) E 0 = Ad(ee E1 )E = µ3 E3 + µ4 E4 (19) Hence, we obtain J3 = E3 + cE4 , c 6= 0, (20) Case 4: µ4 6= 0, µ3 6= 0, µ1 6= 0. Then we have, E = µ1 E1 + µ2 E2 + µ3 E3 + µ4 E4 (21) E 0 = Ad(ee E1 )E = µ1 E1 + µ3 E3 + µ4 E4 (22) So, we obtain J4 = E1 + cE3 + dE4 , c, d 6= 0, (23) Case 5: µ4 = 0, µ3 = 0, µ1 = 0. Then we have, E = µ2 E2 (24) So, we obtain J5 = E2 , (25) Case 6: µ4 = 0, µ3 = 0, µ1 6= 0. Then we have, E = µ1 E1 + µ2 E2 (26) E 0 = Ad(ee E4 )E = µ1 E1 (27) So, we obtain J6 = E1 , (28) Case 7: µ4 = 0, µ3 6= 0, µ1 = 0, µ2 = 0. Then we have, E = µ3 E3 (29) So, we obtain J7 = E3 , (30) Case 8: µ4 = 0, µ3 6= 0, µ1 = 0, µ2 6= 0. Then we have, E = µ2 E2 + µ3 E3 (31) Mathematics 2023, 11, 4545 7 of 13 So, we obtain J8 = E2 + cE3 , c 6= 0, (32) Case 9: µ4 = 0, µ3 6= 0, µ1 6= 0. Then we have, E = µ1 E1 + µ2 E2 + µ3 E3 (33) E 0 = Ad(ee E4 )E = µ1 E1 + µ3 E3 (34) So, we obtain J9 = E1 + cE3 , c 6= 0. (35) 3. Similarity Reductions and Invariant Solutions Vector field J6 = hE1 i. The characteristic equation associated with the vector field E1 = ∂ ∂t is written as dx dt dθ = = , 0 1 0 and provides a transformation θ ( x, t) = k(r ), r = x. With the application of this transfor- mation, we acquire the simplified version of Equation (4) presented as follows, −α2 k(iv) − k00 = 0, (36) which gives, r r k(r ) = c1 + c2 r + c3 sin + c4 cos · α α Hence, the solution of (4) in original variables becomes, x x θ ( x, t) = c1 + c2 x + c3 sin + c4 cos · (37) α α Vector field J7 = hE3 i. The characteristic equation associated with the vector field E3 = ∂ ∂x is written as dx dt dθ = = , 1 0 0 and provides a transformation θ ( x, t) = k(r ), r = x. With the application of this transfor- mation, we acquire the simplified version of Equation (4) presented as follows, k00 = 0, (38) which gives, k (r ) = c1 r + c2 . Hence, the solution of (4) in original variables becomes, θ ( x, t) = c1 t + c2 . (39) Vector field J2 = hE1 + cE4 i. The characteristic equation associated with the vector field E1 + cE4 = ∂ ∂t ∂ + ct ∂θ is written as dx dt dθ = = , 0 1 ct Mathematics 2023, 11, 4545 8 of 13 2 and provides a transformation θ ( x, t) = t2c + k(r ), r = x. With the application of this transformation, we acquire the simplified version of Equation (4) presented as follows, −α2 k(iv) − k00 + c = 0, (40) which gives, r r2 cr k(r ) = −α2 c1 cos − α2 c2 sin + + c3 r + c4 . α α 2 Hence, the solution of (4) in original variables becomes, x x ( x 2 + t2 ) c θ ( x, t) = −α2 c1 cos − α2 c2 sin + + c3 x + c4 . (41) α α 2 Vector field J8 = hE2 + cE3 i. The characteristic equation associated with the vector field E2 + cE3 = ∂ ∂θ ∂ + c ∂x is written as dx dt dθ = = , c 0 1 and provides a transformation θ ( x, t) = xc + k (r ), r = t. With the application of this transformation, we acquire the simplified version of Equation (4) presented as follows, k00 = 0, (42) which gives, k (r ) = c1 r + c2 . Hence, the solution of (4) in original variables becomes, x θ ( x, t) = c1 t + c2 + · (43) c Vector field J3 = hE3 + cE4 i. The characteristic equation associated with the vector field E3 + cE4 = ∂ ∂x ∂ + ct ∂θ is written as dx dt dθ = = , 1 0 ct and provides a transformation θ ( x, t) = cxt + k(r ), r = t. With the application of this transformation, we acquire the simplified version of Equation (4) presented as follows, k00 = 0, (44) which gives, k (r ) = c1 r + c2 . Hence, the solution of (4) in original variables becomes, θ ( x, t) = c1 t + c2 + cxt. (45) Vector field J9 = hE1 + cE3 i. The characteristic equation associated with the vector field E1 + cE3 = ∂ ∂t ∂ + c ∂x is written as dx dt dθ = = , c 1 0 and provides a transformation θ ( x, t) = k(r ), r = t − xc . With the application of this transformation, we acquire the simplified version of Equation (4) presented as follows, (2cβk0 + c2 − α2 )k(iv) + c(c3 + 6βk000 − c)k00 = 0. (46) Mathematics 2023, 11, 4545 9 of 13 We suggest solving the above equation numerically. Vector field J4 = hE1 + cE3 + l E4 i. The characteristic equation associated with the vector field E1 + cE3 + l E4 = ∂ ∂t ∂ + c ∂x + ∂ lt ∂θ is written as dx dt dθ = = , c 1 lt 2 −lx + k (r ), r = t − x . With the application of and provides a transformation θ ( x, t) = 2clxt 2c2 c this transformation, we acquire the simplified version of Equation (4) presented as follows, (−2clβr + 2cβk0 + c2 − α2 )k(iv) + (−6β(l − k00 )k000 + ((c2 − 1)k00 + k)c)c = 0. (47) We suggest solving the above equation numerically. 4. Local Conservation Laws via Homotopy Operator When dealing with complex multipliers and/or equations in the context of inverting divergence operators, an effective approach involves the utilization of homotopy operators derived from differential geometry. These operators help simplify the task of determining fluxes, reducing it to a more manageable problem of integration in single-variable calculus. In this section, we introduce the initial series of such formulas, as outlined in reference [19]. It is crucial to emphasize that when dealing with PDE systems and/or multipliers that are not excessively intricate but involve arbitrary constitutive functions, the direct method of flux computation is typically employed. Let us consider a system of partial differential equations of order k, denoted as 4( x, θ (k) ) = 0, where x = ( x1 , x2 , · · · , x n ) represents the independent variables and θ = (θ 1 , θ 2 , · · · , θ m ) represents the dependent variables. The system can be expressed as a collection of individual equations as follows 4(σ) = 4(σ) ( x, θ, ∂θ, · · · , ∂θ k ) = 0, σ = 1, 2, · · · , N. (48) Now, let us define a local conservation law for the system (48). This law is characterized by a divergence expression Di Ωi = D1 Ω1 + D2 Ω2 + · · · + Dn Ωn , (49) which holds true for all solutions of the system (48). Here, the operators Di represent the total derivatives with respect to the variable xi , while the quantities Di , i = 1, 2, · · · n, correspond to the fluxes associated with the conservation laws. The Euler operator is defined as follows ∂ ∂ ∂ Eθ = − Di + · · · + (−1)s Di1 · · · Dis +···. (50) ∂θ ∂θi ∂θi1 ...is Now, let us define the n-dimensional higher-order Euler operator concerning a function θ ( x1 , x2 , · · · , x n ) ∞ ∞     (s1 ,s2 ··· ,sn ) k1 kn ∂ ∑ ∑ k −s Eθ = ··· ··· D11 1 · · · Dknn −sn (k +k +···k ) , (51) k 1 = s1 k n =sn s1 sn ∂θ 1 2 n ∂(k1 +k2 +···k n ) θ where θ (k1 +k2 +···kn ) = · It is worth noting that when s1 = 0 and sn = 0, we ∂k1 x1 ···∂k n x n (0,0) obtain the original Euler operator = Eθ as defined in (50). Eθ The n-dimensional homotopy operator [24] is introduced for an expression g[θ ] = g( x, θ, ∂θ, · · · ), where θ = (θ 1 ( x ), θ 2 ( x ), · · · , θ m ( x )) and x = ( x1 , x2 , · · · , x n ). This Mathematics 2023, 11, 4545 10 of 13 operator is defined through its n components, each corresponding to a specific independent variable xi , i = 1, 2, · · · , n, and can be expressed as Z 1 m i ( xi ) dΛ H(x ) = 0 j =1 ∑ Ij ( g[θ̂ ]) |θ̂ =Λθ Λ , (52) ( xi ) where Ij ( g[θ̂ ]) is determined using the expression ∞ ∞   ( xi ) 1 + si (s ,··· ,sn ) Ij ( g[θ̂ ]) = ∑ ··· ∑ 1 + s1 + · · · + s n × D1s1 · · · Dsnn (θ̂ j Eθ̂ j 1 ( g[θ̂ ])), (53) s1 =0 s n =0 where j = 1, · · · , m. In essence, the n-dimensional homotopy operator enables us to deal with expressions g[θ ] involving multiple variables θ and x and it provides a systematic way of handling each independent variable xi individually. The main theorem is presented below. Proposition 2. Let g[θ ] be a divergence expression given by [24] g[θ ] = div Ω = D1 Ω1 [θ ] + · · · + · · · , Dn Ωn [θ ], (54) and assume that g[0] = 0. Then the fluxes Ωi can be expressed as i Ωi = H ( x ) ( g[θ ])), i = 1, · · · , n, (55) up to the corresponding fluxes of a trivial conservation law, provided that the integrals (53) converge. By utilizing the direct (multiplier) method to solve Equation (4), we have derived the first-order multipliers as follows x x 1 Ψ( x, t, θ, θ x ) = −α2 (C3 t + C7 ) cos − α2 (C4 t + C8 ) sin + C1 t3 α α 6 (56) 1 1 1 + C2 t + (3C1 x + 6C5 x + 6C6 )t + C2 x2 + C9 x + C10 . 2 2 2 6 2 In the above equation, Ψ( x, t, θ, θ x ) represents the first-order multipliers obtained from the direct (multiplier) method applied to Equation (4). The constants C1 to C10 are coefficients determined during the solution process. Further details and the complete derivation process can be found in the reference [25,26]. • Using the characteristic Ψ1 = 16 t3 + 12 x2 t in Equation (55) and following the integral Formula (52), we obtain  Ω1 = − 31 βθtxx θ x x2 t − βx2 tθtxx + βt2 θ xx θ x + βx2 θ x θ xx − 13 βt2 θtx θ xx  t  − βx2 tθtx θ xx − 21 θx2 + 16 t3 θt + 12 x2 tθt + 12 t2 θ xx + 12 x2 θ xx − 16 t3 θtxx      − 21 x2 tθtxx , (57)   Ω1 = θ xt − α2 tθ x − βtθ x2 − 21 x2 tθ x − 16 α2 t3 θ xxx − tθ x − 16 t3 θ x    x  +α2 xtθ xx − 12 α2 x2 tθ xxx .  • Using the characteristic Ψ2 = 12 t2 + 12 x2 in Equation (55) and following the integral Formula (52), we obtain  Ω2 = − βt θ x θtxx − βx θ x θtxx + 2βtθ x θ xx − βt θtx θ xx − βx θtx θ xx − tθ  t 2 2 2 2  + 21 t2 θt + 12 x2 θt + tθ xx − 12 t2 θtxx − 12 x2 θtxx ,   (58)    Ω2x = α2 xθ xx − 21 α2 t2 θ xxx − 12 α2 x2 θ xxx − 21 t2 θ x − 21 x2 θ x − βθ x2 − xθ − θ − α2 θ .   x x Mathematics 2023, 11, 4545 11 of 13 • Using the characteristic Ψ3 = −α2 t cos ( αx ) in Equation (55) and following the integral Formula (52), we obtain ( Ω3t = α2 cos ( αx )(2βtθtx θ xx + 2βtθ x θtxx − 2βθ x θ xx − tθt + tθtxxx + θ − θ xx ), (59) Ω3x = α3 t(αθ xx cos ( αx ) + θ xx sin ( αx )). • Using the characteristic Ψ4 = −α2 t sin ( αx ) in Equation (55) and following the integral Formula (52), we obtain ( Ω4t = α2 sin ( αx )(2βtθtx θ xx + 2βtθ x θtxx − 2βθ x θ xx − tθt + tθtxx + θ − θ xx ), (60) Ω4x = −α3 t(−αθ xxx sin ( αx ) + θ xx cos ( αx )). • Using the characteristic Ψ5 = xt in Equation (55) and following the integral Formula (52), we obtain ( Ω5t = −2βxtθtx θ xx − 2βxtθtxx θ x + 2βxθ x θ xx + xtθt − xtθtxx − xθ + xθ xx , (61) Ω5x = −α2 xtθ xxx + α2 tθ xx − xtθ x + tθ. • Using the characteristic Ψ6 = t in Equation (55) and following the integral Formula (52), we obtain ( Ω6t = −2βtθtx θ xx − 2βtθtxx θ x + 2βθ x θ xx + tθt − tθtxx − θ + θ xx , (62) Ω6x = −α2 tθ xxx − tθ x . • Using the characteristic Ψ7 = −α2 cos ( αx ) in Equation (55) and following the integral Formula (52), we obtain ( Ω7t = α2 cos ( αx )(2βθtx θ xx + 2βθ x θtxx − θt + θtxx ), (63) Ω7x = α3 (α cos ( αx )θ xxx + θ xx sin ( αx )). • Using the characteristic Ψ8 = −α2 sin ( αx ) in Equation (55) and following the integral Formula (52), we obtain ( Ω8t = α2 sin ( αx )(2βθtx θ xx + 2βθ x θtxx − θt + θtxx ), (64) Ω8x = −α3 (−α sin ( αx )θ xxx + θ xx cos ( αx )). • Using the characteristic Ψ9 = x in Equation (55) and following the integral Formula (52), we obtain ( Ω9t = −2βxθtx θ xx − 2βxθ x θtxx + xθt − xθtxx ), (65) Ω9x = −α2 xθ xxx + α2 θ xx − xθ x + θ. • Using the characteristic Ψ10 = 1 in Equation (55) and following the integral Formula (52), we obtain ( Ω10 t = −2βθtx θ xx − 2βθ x θtxx + θt − θtxx ), (66) Ω10 = −α2 θ xxx − θ x . x 5. Wave Nature of the Obtained Solutions The graphical interpretation of a solution is of paramount importance as it offers a visual representation of complex mathematical relationships. It provides intuitive insights into the behavior, trends, and critical points of the solution. Graphs aid in understanding the sensitivity to varying parameters, validating results, and facilitating communication with a broader audience. They serve as valuable tools for exploration, model selection, and Mathematics 2023, 11, 4545 12 of 13 enhancing our understanding of intricate systems and phenomena. Figures 1 and 2 show the behavior of nonlinear longitudinal waves in an SP medium using the SP model. 40 ◆ ■ ● 30 ◆ ◆ 20 ■ θ ● ■ ● ◆ ● t=0 10 ◆ ◆ ■ t=2 ■ ● ■ ◆ t=4 ■ ● ● 0 -6 -4 -2 0 2 4 6 x (a) 3D plot (b) 2D plot Figure 1. Wave behavior of the SP model (4) with c1 = c2 = c3 = c4 = 0. 30 ◆ 20 ■ ◆ ■ ● 10 ● ◆ ■ ● 0 θ ● ■ ◆ -10 ● ● t=1 ● ■ ◆ ■ t=2 -20 ■ ◆ t=3 ◆ -30 -6 -4 -2 0 2 4 6 x (a) 3D plot (b) 2D plot Figure 2. Wave behavior of the SP model (4) with c1 = c2 = c3 = c4 = 0. 6. Concluding Remarks In this study, we successfully applied the Lie symmetry method to analyze the in- tegrability properties of the SP model. By exploring the group invariant solutions via the reduction in symmetry subalgebras based on the optimal system, we gained valuable insights into the behavior of the system. Additionally, we listed the conservation laws using the homotopy operator, which proved advantageous over the traditional multiplier approach due to the absence of arbitrary functions in both the equation and characteristics. In the existing literature, Erofeev et al. [11] attempted to explore linear and non-linear plane longitudinal waves in the SP medium. However, their work did not present any soliton solutions or invariant ones. Our study focuses on invariant solutions, emphasizing solutions that remain unchanged under specific Lie symmetry transformations. This unique aspect adds novelty to our results. Furthermore, to illustrate the behavior of the solution at specific points in the model, we introduced nonlocal conservation laws reported for the first time. These findings contribute to the enhanced exploration of plane longitudinal waves in the SP medium via the application of the SP model. These findings shed light on the behavior of nonlinear longitudinal waves in the Slepyan–Palmov medium and the impact of various parameters on wave characteristics. These results motivate us to continue using the Lie symmetry method and the homotopy operator in tackling mathematical physics problems to further our understanding and contribute to this field of research. Author Contributions: Conceptualization, A.H., M.U. and F.D.Z.; methodology, M.U., A.H. and F.D.Z.; software, A.H. and M.U.; validation, A.I., F.D.Z. and Y.A.; formal analysis, Y.A. and A.I.; investigation, F.D.Z. and Y.A.; resources, M.U. and A.H.; data curation, A.I.; writing—original draft preparation, M.U. and A.H.; writing—review and editing, M.U., A.H. and F.D.Z.; visualization, A.I. Mathematics 2023, 11, 4545 13 of 13 and Y.A.; supervision, F.D.Z.; project administration, A.I. and Y.A.; funding acquisition, A.I. All authors have read and agreed to the published version of the manuscript. Funding: The authors extend their appreciation to the Deanship of Scientific Research at King Khalid University for funding this work through a large group Research Project under the grant number RGP2/366/44. Data Availability Statement: All data generated or analyzed during this study are included in this published article. Conflicts of Interest: The authors declare no conflict of interest. References 1. Altenbach, H.; Forest, S. (Eds.) Generalized Continua as Models for Classical and Advanced Materials; Springer International Publishing: Cham, Switzerland, 2016. 2. Bagdoev, A.G.; Erofeyev, V.I.; Shekoyan, A.V. Wave Dynamics of Generalized Continua; Springer: Berlin, Germany, 2016. 3. Altenbach, H.; Eremeyev, V.A. Generalized Continua—From the Theory to Engineering Applications; Springer Science & Business Media: Berlin/Heidelberg, Germany, 2012. 4. Cosserat, E.M.; Cosserat, F. Theory of Deformable Bodies; A. Hermann and Sons: New Ulm, MN, USA, 1909. 5. Erofeev, V.I. The Cosserat brothers and generalized continuum mechanics. Vycisl. Meh. Splosh Sred-Comput. Contin. Mech. 2009, 2, 5–10. [CrossRef] 6. Le Roux, J. Research on the geometry of finite deformations. Ann. Ec. Norm. Supér. 1913, 30, 193–245. 7. Tupin, R.A. Theories of elasticity, taking into account moment stresses. Mekh. Sborn. Perevod. 1965, 3, 113–140. 8. Slepyan, L.I. Deformation wave in a rod with damped masses. Inzh. Zh. Mekh. Tverd. Tela. 1967, 5, 34–40. 9. Pal’mov, V. Fundamental equations of the theory of asymmetric elasticity. J. Appl. Math. Mech. 1964, 28, 496–505. [CrossRef] 10. Palmov, V.A. Application of a generalized continuum theory to the problem of space decay in complex mechanical systems. Vychisl. Mekh. Sploshn. Sred. 2009, 2, 105–110. 11. Erofeev, V.I.; Korsakov, M.I.; Leontieva, A.V. Linear and Nonlinear Plane Longitudinal Waves in the Slepyan-Palmov Medium. Mech. Solids 2022, 57, 1385–1395. [CrossRef] 12. Ovsyannikov, L.V. Lectures on the Theory of Group Properties of Differential Equations; World Scientific Publishing Company: Lisbon, Portugal, 2013. 13. Ibragimov, N.H. Optimal system of invariant solutions for the Burgers equation. In Proceedings of the 2nd Conference on Non-Linear Science and Complexity, Session MOGRAN XII, Porto, Portugal, 28–31 July 2008. 14. Bluman, G.; Anco, S. Symmetry and Integration Methods for Differential Equations; Springer Science & Business Media: Berlin/Heidelberg, Germany, 2008. 15. Olver, P.J. Classical Invariant Theory; Cambridge University Press: Cambridge, UK, 1999. 16. Hydon, P.E. Symmetry Methods for Differential Equations; Cambridge University Press: Cambridge, UK, 2000. 17. Stephani, H. Differential Equations: Their Solution Using Symmetries; Cambridge University Press: Cambridge, UK, 1989. 18. Hussain, A.; Kara, A.H.; Zaman, F. Symmetries, associated first integrals and successive reduction of Schrödinger type and other second order difference equations. Optik 2023, 287, 113423. [CrossRef] 19. Darvishi, M.T.; Najafi, M. A modification of extended homoclinic test approach to solve the (3+1)-dimensional potential-YTSF equation. Chin. Phys. Lett. 2011, 28, 040202. [CrossRef] 20. Hussain, A.; Chahlaoui, Y.; Usman, M.; Zaman, F.D.; Park, C. Optimal system and dynamics of optical soliton solutions for the Schamel KdV equation. Sci. Rep. 2023, 13, 15383. [CrossRef] [PubMed] 21. Hussain, A.; Kara, A.H.; Zaman, F.D. New exact solutions of the Thomas equation using symmetry transformations. Int. J. Appl. Comput. Math. 2023, 9, 106. [CrossRef] 22. Usman, M.; Hussain, A.; Zaman, F.D. Invariance analysis of thermophoretic motion equation depicting the wrinkle propagation in substrate-supported Graphene sheets. Phys. Scr. 2023, 98, 095205. [CrossRef] 23. Olver, P.J. Applications of Lie Groups to Differential Equations; Springer Science & Business Media: Berlin/Heidelberg, Germany, 1993. 24. Cheviakov, A.F. Computation of fluxes of conservation laws. J. Eng. Math. 2010, 66, 153–173. [CrossRef] 25. Anco, S.C.; Bluman, G. Direct construction method for conservation laws of partial differential equations Part I: Examples of conservation law classifications. Eur. J. Appl. Math. 2002, 13, 545–566. [CrossRef] 26. Anco, S.C.; Bluman, G. Direct construction method for conservation laws of partial differential equations Part II: General treatment. Eur. J. Appl. Math. 2002, 13, 567–585. [CrossRef] Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.