Electronic Journal of Differential Equations, Vol. 2024 (2024), No. 68, pp. 1–22. ISSN: 1072-6691. URL: https://ejde.math.txstate.edu, https://ejde.math.unt.edu DOI: 10.58997/ejde.2024.68 INSTABILITY OF ENERGY SOLUTIONS, TRAVELLING WAVES, AND SCALING INVARIANCE FOR A FOURTH-ORDER P-LAPLACIAN OPERATOR WITH SUPERLINEAR REACTION JOSÉ LUIS DÍAZ PALENCIA Abstract. This analysis explores the oscillatory behavior of traveling wave solutions for a higher-order p-Laplacian operator with a superlinear reaction term. The study employs an energy-based approach, incorporating general- ized Sobolev spaces to examine relevant properties of the solutions, including oscillations, diffusive mollification, and compact support. Based on this en- ergy framework, the regularity of the involved operator is established. The problem is then reformulated using a traveling wave approach, revealing the oscillatory nature of solutions near the null solution. Numerical simulations are conducted for each wave speed to validate the analytical results, yield- ing the corresponding traveling profiles. Notably, one of the most significant findings is the attraction towards the null critical point, which helps prevent blow-up formation. Finally, the study delves into the equation’s scale-invariant properties, leading to the derivation of self-similar solutions. 1. Problem description and objectives Reaction-diffusion problems have been studied using various forms of diffusive operators, including the classical Gaussian second-order operator, p-Laplacian, higher-order spatial derivatives, porous medium, and thin film operators, which are among the most notable. The choice of a specific diffusion model for a given problem necessitates a deep understanding of the underlying mechanisms involved in the exchange of molecules, temperature, or energy. In many studies, diffusion has been modeled using the concept of random walks of particles (for a comprehensive discussion, refer to [37] and the references therein). The random walk framework allows for a microscopic-level description of particle interactions and this fact may result in further accurate formulations of diffusion at larger scales. Other significant approaches to describing diffusion phenomena can be found in [14, 15], where the concept of free energy, first introduced by Landau and Ginzburg, is employed. This formulation generalizes the classical second-order diffusion by postulating a motion energy that results in a higher-order operator. For example, the analysis in [14] provided a non-homogeneous diffusion expression derived from the free energy concept. In this case, the authors proposed a free energy function 2020 Mathematics Subject Classification. 35K92, 35K91, 35K55. Key words and phrases. Higher order p-Laplacian operator; travelling waves; homotopy; superlinear reaction. ©2024. This work is licensed under a CC BY 4.0 license. Submitted June 18, 2024. Published November 7, 2024. 1 2 J. L. DÍAZ PALENCIA EJDE-2024/68 dependent on the gradient of the substance concentration 1 2k(∇v)2, being v the particles concentration. By transforming the energy into motion gradients through the use of chemical potential, the authors arrived at a fourth-order operator. The properties of this non-homogeneous operator suggest a loss of regularity, partic- ularly posing challenges in establishing a maximum principle. These aspects of regularity loss have been thoroughly discussed in several studies, as referenced in [19, 25, 38]. There are other noteworthy applications where non-regular diffusion principles have been considered. For example, the Keller-Segel equation, which is signifi- cant in biology for modeling cell motion through chemotaxis, is one such case [30]. Furthermore, various regularity analyses have been conducted to obtain smooth solutions to the complex dynamics of chemotaxis, which involve diffusion, reaction, and absorption processes (see [2, 13, 44, 45]). In other fields, non-homogeneous diffusion, such as that of the porous medium type, has been used to model coagulation phenomena in complex vessel geometries [9] and to simulate porosity in peristaltic transportation within Jeffrey-type fluids [20]. These studies underscore the importance of examining the specific diffusion char- acteristics of a phenomenon when developing a model. It is worth noting that diffusion is typically formulated using the regular Gaussian operator derived from the classical Fick’s law. However, the aforementioned studies highlight the value of exploring alternative diffusion mechanisms and gaining a deeper understanding of their mathematical properties. The p-Laplacian operator has been widely applied in physics, chemistry, and engineering. As a representative example, in [10], the authors provide numerical and analytical findings to model in fluid mechanics with a p-Laplacian operator. It is of interest to mention the p-Laplacian formulation in Emden-Fowler equation (see [12]) and the study of p-Laplacian with heterogeneous reaction (see [17]). In the presented analysis, a p-Laplacian diffusion of higher order is considered. The single p-Laplacian operator exhibits the property of finite propagation in com- pact supports. In addition, it is a monotone operator (see [29] and references therein). The formulated problem provides a superlinear reaction for which blow-up pat- terns have been shown to exist ([23]) and is given by: ut = −∆(|∆u|m∆u) + |u|p−1 u, u0(x) ∈ Hn 0 (RN ), n ≥ 1, N > 1. (1.1) Note that the defined operator is referred as the fourth order p-Laplacian (see [23]), where m = p1 − 2, indeed ∆p1,2u = −∆(|∆u|p1−2∆u) . (1.2) As described in [42], the p-Laplacian higher order operator preserves the finite propagation feature from the single p-Laplacian. In addition, we mention that a space of functions Hn 0 is introduced in (1.1) to provide a mollifying space (this will be further defined afterward) to account for potential higher order instabilities close the null critical point. The analysis of problem (1.1) begins with the definition of energy solutions, as proposed for general diffusion in [26]. The problem is subsequently explored EJDE-2024/68 INSTABILITY OF ENERGY SOLUTIONS 3 using Traveling Wave (TW) solutions, an approach introduced in the 1930s by Fisher [22] and Kolmogorov, Petrovskii, and Piskunov [32]. Both studies were grounded in second-order diffusion and a bi-stable nonlinear reaction term of the form f(u) = u(1− u). The primary question they addressed was the existence of a TW velocity that produces a monotone front profile, free from oscillations. Since its inception, the Fisher-KPP model has garnered significant attention, with applications extending to diverse fields such as ecology, biology (see [4, 5, 6]), and non-Newtonian fluids [33]. Furthermore, various bi-stable equations have been analyzed similarly to the Fisher-KPP model (see [41] and references therein). In some cases, these equations have been extended to incorporate non-homogeneous diffusion characterized by higher-order diffusive phenomena, as seen in the Ex- tended Fisher-Kolmogorov equation in bi-stable systems [11, 16, 39]. Additionally, a p-Laplacian operator has been introduced to model a Fisher-KPP type equation in [7]. TW solutions have become increasingly important in applied sciences (see [40], [36], and [18]), and their analysis has been extended to the study of higher-order operators (see [28, 34, 24, 25]). This study aims to analyze the oscillatory profiles of TW solutions using both analytical and numerical approaches, the latter serving to validate the analytical results. The analysis is based on the energy formulation introduced in [26], and it employs functional spaces to represent the properties of solutions. In particular, generalized Sobolev spaces are utilized to account for compact support functions, mollifiers, and oscillatory patterns. Problem (1.1) is then reformulated in the TW framework, where the oscillatory properties of solutions—often referred to as so- lution instabilities—are examined. The spectrum of the higher-order p-Laplacian operator is explored using a homotopy representation near the null solution. Fi- nally, a numerical approach confirms the attractive properties of the null solution and the stabilizing effect of this null solution as the TW speed increases. It is important to note that the mathematical treatment in the TW domain begins with step-like initial data. Although this step function is not compactly supported, as required by problem (1.1), it is assumed that such an initial function provides a suitable starting condition to study the evolution of a positive mass alongside a null state. The methodology employed in this work combines an analytical approach to demonstrate the instabilities of the traveling wave (TW) profiles with a numerical validation using the bvp4c function in Matlab to support the analysis. It is impor- tant to note that the analysis introduced by Galaktionov in [23] discusses blow-up profiles, yielding exact patterns. In the present analysis, we provide evidence of instabilities in the TW solutions prior to the onset of blow-up. Additionally, we characterize the null solution as an attractor, which significantly influences the single-point blow-up behavior. Specifically, while solutions exhibit oscillations as they approach the null solution, this state acts as an attractor, thereby preventing the formation of blow-up patterns. 2. Preliminaries Firstly, we consider the definition of a generalized energy solution. 4 J. L. DÍAZ PALENCIA EJDE-2024/68 Definition 2.1. u(x, t) is said to be an energy solution to problem (1.1) in (0, t) if for any 0 < τ < t, the following holds (see [26])∫ τ 0 ∫ RN utξ, dx dt+ ∫ τ 0 ∫ RN −|∆u|m∆u∆ξ dx dt+ ∫ τ 0 ∫ RN |u|p−1uξ dx dt = 0, (2.1) for any arbitrary function ξ ∈ L1(0, τ ;H 2 0 (RN )). Consider the following proposition based on results in [1] and [8]. Proposition 2.2. Let F,G,H ∈ C n 0 (RN ), with (n ≥ 1), the following anisotropic Sobolev inequality holds,∫ RN |FGH| dx ≤ K∥F∥ α−1 α Lq ∥∇F∥1/αLs ∥G∥ α−2 α L2 ∥∇G∥1/αL2 ∥H∥L2 , (2.2) where K > 0, α−1 q + 1 s = 1, α > 2, and 1 ≤ q, s < ∞. According to [26], the asymptotic behaviour of solutions exhibiting blow-up de- pends on the asymptotic behaviour of a rescaled kernel for a linear fourth-order parabolic operator of the form ut = −∆2u. This rescaled kernel exhibits an expo- nential behaviour of the form ∼ e−x4/3 . In [26], it is shown that the asymptotic behaviour of any possible blow-up profile for the higher-order p-Laplacian operator exhibits a similar exponential profile. The asymptotic estimates lead to the pro- posal of bundles of blow-up profiles following S and HS regimes. As our intention is to characterize oscillatory exponential bundles of solutions close to the critical points, a norm is introduced with an appropriate weight to scale out the oscillatory profiles close to the null condition. Definition 2.3. Let ∥h∥2Θ = ∫ RN Θ(z) 4∑ j=0 |Djh(z)|2dz, (2.3) where D = d dz , h ∈ HΘ(RN ) ⊂ L2 Θ(RN ) ⊂ L2(RN ). The weight Θ(z) is defined in [35, 26, 25] as Θ(z) = ec0|z| 4/3 , (2.4) where c0 > 0. Definition 2.4. The following fundamental problem is defined as ht = ∆m,2h, (2.5) where ∆m,2 = −∆(|∆ · |m∆). The next definition provides the mollifying exponential kernel. Definition 2.5. Consider the weighted Sobolev norm defined as ∥h∥2Hm ρ = ∫ ∞ −∞ emξ2 |θ(ξ, t)|2dξ, (2.6) that satisfies the Ap-condition of mollifying kernels for p = 1 (refer to [27]). Finally, if the solutions are sufficiently far from the critical null condition and are sufficiently smooth, eliminating the need for the previously defined mollifier, we introduce the classical Sobolev norm. EJDE-2024/68 INSTABILITY OF ENERGY SOLUTIONS 5 Definition 2.6. The usual Sobolev order n functional space is defined as n(RN ) = {h ∈ L2(RN ) : ∇n(h) ∈ L2(RN )} with the norm ∥h∥n = ∥h∥L2 + ∥∇nh∥L2 . (2.7) It should be mentioned that the above definition applies as well to compactly supported functions in accordance with the similar definition of the space Hn 0 (RN ). Consider, now, a sequence of open bounded domains B(0, β) ⊂ RN , β ∈ N, β = 1, 2, 3 . . . Proposition 2.7. Given the Sobolev space Wn,p(B(0, β)), Define l = int{n− N p }. The following inclusion is continuous (see [31, p. 79]), Wn,p(B(0, β)) ↪→ Cl(B(0, β)). (2.8) Given the particular problem in (1.1), any solution is, at least, weakly differen- tiable up to order four, then n = 4 along with p = 2, leading to l = int{4− d 2}. 3. Boundedness of solutions The following lemma provides a-priori bounds. Lemma 3.1. Lett h be a solution to the fundamental equation ht = ∆m,2h. Given h0 ∈ L2(RN ), and assuming that m ∈ 2N, then the following bounds hold: ∥h∥L2 ≤ ∥h0∥L2 , ∥h∥Hn 0 ≤ ∥h0∥Hn 0 , ∥h∥Hm ρ ≤ A0∥h0∥L2 , A2 0 = e mξ2M−tξ4M γ(m+1) π(i ξM )m+1 , ξM = ( 2mπim+1 tγ(m+ 1)(3−m) ) 1 1−m , ∥h∥Hm ρ ≤ ∥h0∥Hm ρ , ∥h∥Hm ρ ≤ A0∥h0∥Hn 0 , ∥h∥Θ ≤ σ∥h0∥Hm ρ , ∥h∥Θ ≤ ( sup z∈Bβ {ec0|z| 4/3 })1/2σ∥h0∥Hn 0 , where σ2 = 25 sup x∈Bβ {h,D1h,D2h,Dh, D4h} and n can take values 1, 2, 3, 4. Note that Bβ = B(0, β) refers to the ball used in the Proposition 2 that is considered for β ≫ 1 and numerically ordered. These last results state that the scaling norm HΘ is bounded by the mollifier norm Hm ρ and the compacting norm Hn 0 . Proof. Departing from the fundamental problem ht = ∆m,2h, a general solution is expressed as h(x, t) = et∆m,2h0(x). The Fourier transformation (in the ξ variable) admits the following inequality obtained by a direct convolution of the Fourier transformation for each term in (|∆ · |m ∆) and weighted by the −∆ term, ĥ(ξ, t) ≤ e −tξ4 γ(m+1) 2π(i ξ)m+1 ĥ0(ξ), (3.1) where γ(m+ 1) refers to the Gamma function and i to the imaginary unit. Firstly, the following bound is shown for m ∈ 2N, ∥h∥2L2 ≤ ∫ |e−tξ4 γ(m+1) π(i ξ)m+1 |∥ĥ0(ξ)∥2dξ ≤ sup ∀ξ∈RN {|e−tξ4 γ(m+1) π(i ξ)m+1 |} ∫ ∥ĥ0(ξ)∥2dξ = ∥h0∥2L2 . 6 J. L. DÍAZ PALENCIA EJDE-2024/68 Then ∥h∥L2 ≤ ∥h0∥L2 . As a direct consequence, and considering the expression in (2.7), we have ∥h∥Hn 0 ≤ ∥h0∥Hn 0 . (3.2) Now, assume that h0 ∈ L2(RN ). Then ∥h∥2Hm ρ = ∫ ∞ −∞ emξ2 |ĥ(ξ, t)|2dξ ≤ sup ∀ξ∈RN {emξ2e −tξ4 γ(m+1) π(i ξ)m+1 } ∫ ∥ĥ0(ξ)∥2dξ. (3.3) Making standard operations, the following holds ∥h∥2Hm ρ ≤ e mξ2M−tξ4M γ(m+1) π(i ξM )m+1 ∥h0∥2L2 , (3.4) where ξM = ( 2mπim+1 tγ(m+ 1)(3−m) ) 1 1−m . (3.5) Since h0 ∈ L2(RN ), the rapid decay of e−tξ4·c(m,iξ), where c(m, iξ) is a constant that depend on m and the imaginary unit, ensures the integral ∫ emξ2 |ĥ0(ξ, t)|2dξ is finite; thus, h0 ∈ Hm ρ . We show then that any solution, h, to the fundamental equation in Hm ρ satisfies ∥h∥2Hm ρ = ∫ emξ2 |ĥ(ξ, t)|2dξ ≤ sup ∀ξ∈RN {|e−tξ4 γ(m+1) π(iξ)m+1 |} ∫ emξ2∥ĥ0(ξ)∥2dξ ≤ ∥h0∥2Hm ρ . (3.6) The following bound is also applicable: ∥h∥2Hm ρ = ∫ emξ2 |ĥ(ξ, t)|2dξ ≤ sup ∀ξ∈RN {|emξ2−tξ4 γ(m+1) π(i ξ)m+1 |} ∫ ∥ĥ0(ξ)∥2dξ ≤ sup ∀ξ∈RN {|emξ2−tξ4 γ(m+1) π(i ξ)m+1 |} ∫ (∥ĥ0(ξ)∥2 + |∇nh0|2)dξ = e mξ2M−tξ4M γ(m+1) π(iξM )m+1 ∥h0∥2Hn 0 , (3.7) where ξM is given in (3.5). Now, the intention is to show a bound for the norm defined in (2.3). The bounding term is given by the mollifier norm introduced in (2.6): ∥h∥2Θ = ∫ Θ(z) 4∑ j=0 |Djh(z)|2dz ≤ ∫ emz2 4∑ j=0 |Djh(z)|2dz ≤ σ2 ∫ emz2 |h(z)|2dz = σ2∥h∥2Hm ρ ≤ σ2∥h0∥2Hm ρ , (3.8) EJDE-2024/68 INSTABILITY OF ENERGY SOLUTIONS 7 being σ2 = 25 supx∈Bβ {h,D1h,D2h,Dh, D4h}, for β ≫ 1 and ordered. The conti- nuity inclusion in the Proposition 2.7 provides the conditions to ensure the existence of the derivatives of a function h ∈ W 4,2(B(0, β)). In addition, the following holds: ∥h∥2Θ = ∫ Θ(z) 4∑ j=0 |Djh(z)|2dz ≤ sup ∀z∈Bβ {ec0|z| 4/3 }σ2 ∫ (|h(z)|2 + |∇nh|2)dz ≤ sup ∀z∈Bβ {ec0|z| 4/3 }σ2∥h∥2Hn 0 ≤ sup ∀z∈Bβ {ec0|z| 4/3 }σ2∥h0∥2Hn 0 , (3.9) where n can take a value 1, 2, 3, 4 and β ≫ 1 and ordered to expand up to the whole RN. Based on the proposed arguments, the Lemma postulations are proved. □ The next objective is to show the local bound properties of compactly supported solutions. To this end, we assume that the following conditions hold: η ∈ L1(0, τ ;H 2 0 (RN )) ∩ C2(0, τ ;H2 0 ∩ C n 0 (RN )), u0(x) ∈ Hn 0 (RN ) ∩ Cn 0 (RN ), (3.10) with n ≥ 1. Lemma 3.2. Each energy solution satisfying (2.1) is bounded in H4 0 (RN ) (see norm (2.7)), i.e. the compact support is locally preserved. Proof. Firstly, the expression (2.1) is rewritten as∫ τ 0 ∫ RN utη dx dt+ ∫ τ 0 ∫ RN |u|p−1uη dx dt = ∫ τ 0 ∫ RN |∆u|m∆u∆η dx dt. (3.11) Note that under the conditions stated in (3.10), Proposition 2.2 can be used to further develop the left-hand side integral. To this end, admit q = 1, s = 2 and α = 2 in (2.2), then∫ RN ∆η∆u|∆u|m dx ≤ K∥∆η∥1/2L2 ∥∇∆η∥1/2L2 ∥∇∆u∥1/2L2 ∥∆u∥m+1 L2 . (3.12) Furthermore, based on Definition 2.6 and Lemma 3.1, the following bounds apply: K∥∆η∥1/2L2 ∥∇∆η∥1/2L2 ∥∇∆u∥1/2L2 ∥∆u∥m+1 L2 ≤ K∥η∥1/2 H2 0 ∥η∥1/2 H3 0 ∥u∥1/2 H3 0 ∥u∥m+1 H2 0 ≤ K∥η∥1/2 H2 0 ∥η∥1/2 H3 0 ∥u0∥1/2H3 0 ∥u0∥m+1 H2 0 (3.13) In addition, ∫ RN |u|p−1uη dx ≤ k∥η∥1/2L2 ∥∇η∥1/2L2 ∥u∥1/2L2 ∥∇u∥1/2L2 ∥u∥p−1 L2 ≤ k∥η∥1/2L2 ∥η∥1/2 H1 0 ∥u∥1/2 H1 0 ∥u∥p− 1 2 L2 ≤ k∥η∥1/2L2 ∥η∥1/2 H1 0 ∥u0∥1/2H1 0 ∥u0∥ p− 1 2 L2 (3.14) Now, we considering the natural Sobolev embedding, k∥η∥1/2L2 ∥η∥1/2 H1 0 ∥u0∥1/2H1 0 ∥u0∥ p− 1 2 L2 ≤ k∥η∥1/2L2 ∥η∥1/2 H1 0 ∥u0∥pH1 0 . (3.15) 8 J. L. DÍAZ PALENCIA EJDE-2024/68 Next we return to the first integral in (3.11), and by using the Gronwall inequality, we hav ∫ τ 0 utη dt ≤ ∫ τ 0 g(t)ηu dt, (3.16) where g(t) is a continuous function coming from the application of the conditions required by the Gronwall inequality. Note that the auxiliary function η is C2(0, τ) as expressed in (3.10). Then, η can be selected such that ηg(t) ≤ 1. Consequently,∫ τ 0 utη dt ≤ ∫ τ 0 u dt. (3.17) Assuming that the constant K in (3.13) is sufficiently large while the constant k is considered sufficiently small, by (3.11),∫ τ 0 ∥u∥H4dt+ k ∫ τ 0 ∥η∥1/2L2 ∥η∥1/2 H1 0 ∥u0∥pH1 0 dt ≤ K ∫ τ 0 ∥η∥1/2 H2 0 ∥η∥1/2 H3 0 ∥u0∥1/2H3 0 ∥u0∥m+1 H2 0 dt. (3.18) Note that the functions η and u0 satisfy the conditions expressed in (3.10), in particular the compact support. Consequently and locally in time (i.e. locally in the proximity of a sufficiently small τ to keep the support): ∥u∥H4 0 + k∥η∥1/2L2 ∥η∥1/2 H1 0 ∥u0∥pH1 0 ≤ K∥η∥1/2 H2 0 ∥η∥1/2 H3 0 ∥u0∥1/2H3 0 ∥u0∥m+1 H2 0 , (3.19) which permits to account for the bound properties of any solution in H4. To this end, it suffices to write ∥u∥H4 0 ≤ K∥η∥1/2 H2 0 ∥η∥1/2 H3 0 ∥u0∥1/2H3 0 ∥u0∥m+1 H2 0 + k∥η∥1/2L2 ∥η∥1/2 H1 0 ∥u0∥pH1 0 , (3.20) which allows us showing the lemma postulations. □ The next lemma aims at exploring the bound of oscillating solutions by the defined mollifier in (2.6). To this end, the support is kept free and oscillating. As a consequence, the following is required previously, η ∈ L1(0, τ ;H m ρ (RN )) ∩ C2(0, τ ;Hm ρ ∩ C n(RN )), u0(x) ∈ Hm ρ (RN ) ∩ C n(RN ), (3.21) Lemma 3.3. Each oscillating energy solution satisfying (2.1) is globally bounded by the mollification introduced in the norm (2.6). Proof. To show the proposed lemma, we use (3.11) and operate with the inequality (2.2) (with q = 1, s = 2, α = 2) along with the norm (2.3). Then∫ RN ∆η∆u |∆u|m dx ≤ K∥∆η∥1/2L2 ∥∇∆η∥1/2L2 ∥∇∆u∥1/2L2 ∥∆u∥m+1 L2 ≤ K∥η∥1/2 Θ ∥η∥1/2 Θ ∥u∥1/2 Θ ∥u∥m+1 Θ . (3.22) Now, based on Lemma 3.1, the following bound in Hm ρ is obtained. K∥η∥1/2 Θ ∥η∥1/2 Θ ∥u∥1/2 Θ ∥u∥m+1 Θ ≤ Kσ∥η∥Hm ρ ∥u0∥ m+ 3 2 Hm ρ . (3.23) EJDE-2024/68 INSTABILITY OF ENERGY SOLUTIONS 9 Operating similarly,∫ RN |u|p−1uη dx ≤ k∥η∥1/2L2 ∥∇η∥1/2L2 ∥u∥1/2L2 ∥∇u∥1/2L2 ∥u∥p−1 L2 ≤ k∥η∥1/2Θ ∥η∥1/2Θ ∥u∥1/2Θ ∥u∥p− 1 2 Θ ≤ kσ∥η∥Hm ρ ∥u0∥pHm ρ (3.24) Now, considering (3.11), applying the Gronwall inequality similarly as in (3.17), the following inequality holds for K sufficiently large and k small,∫ τ 0 ∥u∥Θdt+ kσ ∫ τ 0 ∥η∥Hm ρ ∥u0∥pHm ρ dt ≤ Kσ ∫ τ 0 ∥η∥Hm ρ ∥u0∥ m+ 3 2 Hm ρ dt. (3.25) Note that the functions η and u0 satisfy (3.21). Then any oscillating solution, under the norm (2.3), in (0, τ) is mollified by the norm (2.6). Based on this, the following bound holds, ∥u∥Θ ≤ kσ∥η∥Hm ρ ∥u0∥pHm ρ +Kσ∥η∥Hm ρ ∥u0∥ m+ 3 2 Hm ρ . (3.26) □ 4. Travelling waves The travelling waves (TW) solutions are given by the change u(x, t) = γ(ξ) where ξ = x · nd − λt. Note that ξ ∈ R and the vector nd ∈ RN represents the TW propagating direction. In addition, note that λ represents the TW velocity and γ : R → (0,∞) is the TW profile that complies with the norm (2.3) i.e. γ ∈ HΘ(R) ⊂ L2 Θ(R) ⊂ L2(R). It is to be noted that two TWs are equivalent under discrete symmetry (ξ → −ξ) and translation (ξ → ξ + ξ0). We assume that the TW direction of motion is nd = (1, 0, 0, . . . , 0), such that ξ = x − λt and u(x, t) = γ(ξ) ∈ R. Using the described transformation into the TW-domain, the problem in (1.1) is then reformulated as −λγ′ = −(|γ′′|mγ′′)′′ + |γ|p−1γ. (4.1) The next lemma aims at characterizing the TW propagating direction. This step is relevant to properly describe the wave dynamics when performing the numerical assessments. Lemma 4.1. The TW velocity λ is positive, equivalently, the TW motion departs from ξ → −∞ and ends in ξ → ∞. Proof. Multiply (4.1) by γ′:, −λ(γ′)2 = −(|γ′′|mγ′′)′′γ′ + |γ|p−1γγ′, (4.2) and consider the integration from −∞ to ∞. We start the evaluation of the integrals by the diffusive term:∫ (|γ′′|mγ′′)′′γ′ = γ′(|γ′′|mγ′′)′ − ∫ (|γ′′|mγ′′)′γ(2) = γ′(|γ′′|mγ′′)′ − ( γ′(|γ′′|mγ′′)γ(2) − ∫ γ′(|γ′′|mγ′′)γ(3) ) . (4.3) 10 J. L. DÍAZ PALENCIA EJDE-2024/68 Note that the given integrals are determined in the limit −∞ to +∞, such that the following asymptotic conditions are Assumed to hold: γ′(−∞) = γ(2)(−∞) = γ(3)(−∞) = 0 , γ′(∞) = γ(2)(∞) = γ(3)(∞) = 0. (4.4) Consequently, ∫ (|γ′′|mγ′′)′′γ′ = 0. (4.5) Now, for the term |γ|p−1γγ′, we shall considered that there exists a jump of magnitude c in the TW profile. This fact will be further developed in the numerical exercise, but as a first description, it shall be noted that such jump represents a transition from a finite mass initial condition to a null condition (at infinity) where the instabilities are characterized. Then∫ |γ|p−1γγ′ = γpγ − p ∫ γp = γpγ − p ∫ γp−2γ2 = γpγ − p p− 1 γp+1 + 2p p− 1 ∫ γpγ′. (4.6) Consequently, based on the consideration of the mentioned jump of magnitude c, we have ∫ |γ|p−1γγ′ = γpγ − p p−1γ p+1 1− 2p p−1 = γp+1 p+ 1 = 1 p+ 1 (γp+1(∞)− γp+1(−∞)) = 1 p+ 1 (0− cp+1). (4.7) Finally and upon recovery of the expression (4.2), the following holds, −λ ∫ (γ′)2 = 0− 1 p+ 1 cp+1, (4.8) which leads to λ = 1 p+ 1 cp+1∫ (γ′)2 . (4.9) This last expression permits us to conclude that λ > 0 as claimed. □ 4.1. Travelling wave instabilities. The TWs formulation in expression (4.1) has one critical solution at γ = 0. The aim of this subsection is to analyze the oscillating properties of any solution in the proximity of the null solution. The oscillations induced by the higher order non-linear diffusion are studied with the introduced norm in (2.3). Afterward, this chapter aims at showing that the null solution acts as an attractor of solutions hindering the possible nucleation of blow-up profiles. The study of oscillations close to the null solution follows from a theorem intro- duced to study the Kuramoto-Sivashinsky equation (see [43] and references therein) along with other equations, particularly the Cahn-Hilliard equation (see [28]) and a sixth order diffusion equation (see [34]). Nonetheless, for our present case, the higher order p-Laplacian operator induces a set of changes in the mentioned theo- rem for the cited equations. To this end, the theorem is divided into four Lemmas, so that a instability statement is proved. EJDE-2024/68 INSTABILITY OF ENERGY SOLUTIONS 11 The fist lemma introduces the principle of oscillations (or instabilities) for any energy solution. This is shown based on the already introduced norms in (2.3) and (2.7). Lemma 4.2. Each oscillating energy solution u(x, t) ∈ L2(RN ) is bounded by the Sobolev norms (2.3) and (2.7) i.e. ∥u∥L2 ≤ K1∥u∥H4 , ∥u∥L2 ≤ K2∥u∥Θ. (4.10) Proof. The first condition is trivially shown in virtue of the norm H4 defined in (2.7), ∥u∥H4 = ∥u∥L2 + ∥∇4u∥L2 ≥ ∥u∥L2 . (4.11) Given the positivity of any norm, this last expression concludes on ∥u∥L2 ≤ ∥u∥H4 , i.e. K1 = 1. The next inequality is shown considering the expression (2.3): ∥u∥2L2 ≤ ∫ RN 4∑ j=0 |Dju(z)|2dz ≤ ∫ RN Θ(z) 4∑ j=0 |Dju(z)|2dz = ∥u∥2Θ, (4.12) a.e. in RN . Then, it suffices to consider K2 = 1. □ Now, to introduce the TW convergence analysis, we introduce the function w(x, t) = u(x, t)− φ(x, t), where φ(x, t) represents a perturbation randomly small so as to ensure the TW profiles convergence. Particularly and close to the null solution, φ(x, t) is requested to satisfy (see Lemma 3.1 along with inequality (3.8)): ∥φ∥L2 ≤ ∥φ∥Θ ≤ σ∥φ∥Hm ρ . (4.13) Convergence requires σ → 0, i.e. a mollification on the function φ. Now, the problem (1.1) is hence formulated in terms of w(x, t) and φ(x, t) as wt + φt = −∆( ∞∑ k=0 ( m k ) (∆w)m−k(∆φ)k ∆(w + φ)) + ∞∑ j=0 ( p j ) wp−jφj . (4.14) Now, assume that for any stationary perturbation it holds that 0 < ∥φ∥Θ ≤ C. In addition, wt = F (w), (4.15) where F (w) = −∆( ∑∞ k=0 ( m k ) (∆w)m−k(∆φ)k∆(w + φ)) + ∑∞ j=0 ( p j ) wp−jφj . Lemma 4.3. The mapping F : Hm ρ → L2 is continuously bounded. In addition, there exist α0 > 0, K3 > 0 and α0 > 1 such that ∥F (w)∥L2 ≤ K3∥w∥α0 Hm ρ , provided 0 < ∥w∥Hm ρ < α0. Proof. We have ∥F (w)∥L2 ≤ ∥F (w)∥Θ ≤ ∞∑ k=0 ( m k ) ∥w∥m−k Θ ∥φ∥kΘ(∥w∥Θ + ∥φ∥Θ) + ∞∑ j=0 ( p j ) ∥w∥p−j Θ ∥φ∥jΘ . Considering that 0 < ∥φ∥Θ ≤ C and inequality (3.8), ∥F (w)∥L2 ≤ ∥F (w)∥Θ ≤ ∞∑ k=0 ( m k ) 2∥w∥m−k+1 Hm ρ Ck+1 + ∞∑ j=0 ( p j ) ∥w∥p−j Hm ρ Cj . (4.16) 12 J. L. DÍAZ PALENCIA EJDE-2024/68 For m− k + 1 > p− j, the following holds ∥F (w)∥L2 ≤ ∥F (w)∥Θ ≤ ∞∑ k=0 ( m k ) 3∥w∥m−k+1 Hm ρ Ck+1 = K3∥w∥α0 Hm ρ (4.17) Then, it suffices to consider K3 = ∑∞ k=0 ( m k ) 3Ck+1 and α0 = max{m− k+1} > 1. Considering that p− j > m− k + 1, we have ∥F (w)∥L2 ≤ ∥F (w)∥Θ ≤ ∞∑ j=0 ( p j ) 3∥w∥p−j Hm ρ Cj+1 = K3∥w∥α0 Hm ρ . (4.18) In this case, it suffices to admit K3 = ∑∞ j=0 ( p j ) 3Cj+1 and α0 = max{p− j} > 1. The continuity can be shown as a consequence of the proved inequalities by considering a pair of sequences sufficiently close. This can be done by standard assessments. □ Consider the problem wt = −∆ ( ∞∑ k=0 ( m k ) (∆w)m−k(∆φ)k∆(w + φ) ) + ∞∑ j=0 ( p j ) wp−jφj = Lw +G(w), (4.19) such that Lw = −∆( ∞∑ k=0 ( m k ) (∆w)m−k(∆φ)k∆(w + φ)) and G(w) = ∑∞ j=0 ( p j ) wp−jφj . Based on this, we consider the abstract evolution w(x, t) = et L w0(x). Lemma 4.4. L is the infinitesimal representation of a strongly continuous semi- group given by etL that satisfies∫ 1 0 ∥etL∥L2→Hm ρ = K5 < ∞, ∫ 1 0 ∥etL∥L2→HΘ = K6 < ∞ (4.20) Proof. The proof of this lemma is based on Lemma 3.1. Then ∥w∥Hm ρ ≤ ∥w0∥Hm ρ ≤ A0∥w0∥L2 . (4.21) Based on the abstract evolution w(x, t) = etLw0(x), we have ∥w∥Hm ρ ≤ ∥etL∥L2→Hm ρ ∥w0∥L2 . (4.22) Then ∫ 1 0 ∥etL∥L2→Hm ρ = ∫ 1 0 A0 = K5 < ∞. (4.23) It is easy to check that the value of A0 provides a finite value of K5 upon integration in t ∈ (0, 1]. Operating similarly, it is possible to conclude on the bound of the abstract evolution in HΘ. To this end, it suffices to consider the Lemma 3.1 along with the inequality (3.8), ∥w∥Θ ≤ σ ∥w∥Hm ρ ≤ σ∥w0∥Hm ρ ≤ σA0∥w0∥L2 . (4.24) Again, based on the abstract evolution, ∥w∥HΘ ≤ ∥etL∥L2→HΘ ∥w0∥L2 , (4.25) EJDE-2024/68 INSTABILITY OF ENERGY SOLUTIONS 13 so that ∫ 1 0 ∥etL∥L2→HΘ = ∫ 1 0 A0σ = K6 < ∞, (4.26) which is finite upon integration in t ∈ (0, 1]. □ The next step is to analyze the spectrum of the operator L (as defined in (4.19)) using the norm HΘ. Lemma 4.5. Given the condition (4.13) to the perturbation terms and that 0 < ∥φ∥Θ ≤ C, the spectrum of L (see (4.19)) in HΘ, close to the null solution, has at least an eigenvalue (ϕ) such that Re(ϕ) > 0. Proof. This proposed lemma can be shown through the theory of Evans functions. Indeed, this theory permits determining the location of positive eigenvalues in the proximity of the null solution. The roots of the Evans functions coincide with the characteristic polynomial roots of a linearized operator close to the null solution (see [3] and references therein). In the presented analysis, the eigenvalues are obtained by using the characteristic polynomial near the equilibrium γ = 0. Additionally, the particular dynamics, affected by the TW speed λ, are assessed with a homotopy representation. To this end, a computational exercise is introduced. A first integral can be derived in the TW problem (4.1) close to the null solution, so that −λγ′ = −(|γ′′|mγ′′)′′ → −λγ = −(|γ′′|mγ′′)′ + c1, (4.27) where the constant c1 = 0 for the sake of simplicity and without impacting the lemma results. Now, the problem is converted into the matrix formulationγ0 γ2 γ3 ′ =  0 1 0 0 0 1 λ mγm−1 3 +γm 3 0 0 γ0 γ2 γ3  (4.28) Note that the characteristic polynomial (in the assumption that γ3 is a free param- eter) for the matrix is Q(ϕ) = −ϕ3 + λ mγm−1 3 + γm 3 = 0. (4.29) Note that γ3 = γ′′, as per the standard change of variables to build the matrix rep- resentation. For our purposes (i.e. to show that there exists at least one eigenvalue with positive real part) we consider that in the asymptotic approximation to the null state |γ3| ≪ 1 then the following variable Υ ≫ 1 is introduced, −ϕ3 +Υλ = 0. (4.30) By a standard resolution of the last characteristic polynomial, it is easy to conclude on the existence of at least one eigenvalue with positive real part. In addition, it is necessary to determine the effect of the TW speed. To this end, the different homotopy graphs containing the eigenvalues are given for different values in the TW speed (see Figures 1, 2 for positive TW speeds and Figures 3, 4 for negative TW speeds). □ 14 J. L. DÍAZ PALENCIA EJDE-2024/68 Figure 1. Eigenvalues representations in the complex plane for Q(ϕ) roots. The value of Υ has been taken arbitrary big. Note that Υ affects only on the scale while keeping the structure of the eigenvalues (one with positive real part). λ = 1 (left) and λ = 10 (right). Figure 2. Eigenvalues representations in the complex plane for Q(ϕ) roots. The value of Υ has been taken arbitrary big. Note that Υ affects only on the scale while keeping the structure of the eigenvalues (one with positive real part). λ = 100 (left) and λ = 1000 (right). Through the provided series of lemmas, it was shown that any oscillating energy solution is bounded by Sobolev norms, which prevents the formation of blow-up profiles. Additionally, the study of the spectrum of the operator L in the HΘ norm indicated the presence of at least one eigenvalue with a positive real part, suggesting an instability in the system. This instability is influenced by the speed of the traveling wave, as demonstrated by the analysis of eigenvalues for different TW speeds. 4.2. Exact travelling wave profiles and characteristic propagation speed. In this section, we provide the TWs profiles to validate the results obtained in the previous section 4.1. To this end, a numerical approach has been followed using the solver bvp4c in Matlab. This function consists on a Runge-Kutta implicit algorithm supported by interpolant extensions [21]. To build the solution, the bvp4c requires to solve a collocation method for which the conditions at ξ → −∞ and ξ → ∞ shall be specified. In the presented analysis, the condition at −∞ is admitted to be positive (to this end it suffices to consider γ(−∞) = 1) while the condition at EJDE-2024/68 INSTABILITY OF ENERGY SOLUTIONS 15 Figure 3. Eigenvalues representations in the complex plane for Q(ϕ) roots. The value of Υ has been taken arbitrary big. Note that Υ affects only on the scale while keeping the structure of the eigenvalues (one with positive real part). λ = −1 (left) and λ = −10 (right). Figure 4. Eigenvalues representations in the complex plane for Q(ϕ) roots. The value of Υ has been taken arbitrary big. Note that Υ affects only on the scale while keeping the structure of the eigenvalues (one with positive real part). λ = −100 (left) and λ = −1000 (right). ∞ is kept free. This last condition is particularly relevant to characterize the null solution as an attractor. The numerical exploration has been done over a large interval in ξ ∈ [−100, 1000] so that the problem is not governed by the collocation method values at −∞ and ∞. In addition, and to make the problem tractable, dedicated values have been introduced for the parameters involved in Problem (1.1) without loss of generality. Particularly, it has been considered m = 3 and p = 2. It should be noted that the value of m has been chosen as an odd number to complement the hypothesis in Lemma 3.1, where m is assumed to be even. This choice is made to provide additional information and to explore the numerical behavior of solutions with an odd value of m. Then, for these values, solutions are represented for a wide interval of TW-speeds. It is possible to check that in all cases the null solution acts as an attractor (remind that the collocation required by the bvp4c solver is kept free 16 J. L. DÍAZ PALENCIA EJDE-2024/68 at ∞ where the null state occurs). In addition, the overall instabilities magnitude increases for decreasing values in the TW-speed. As a consequence of the exposed analysis, it is concluded that it is not possible to find a suitable TW-speed for which a positive inner region can be shown (see [25] for a complete discussion) impeding the possibility of formulating a maximal kernel with purely monotone behaviour. Even further, it is not possible to conclude on a TW-speed for which the TW profile is positive in the whole space as in the classical order two KPP-problem [32]. Other values of m have also been considered (even and odd) in the analysis, but the conclusions remain the same. Regardless of whether m is chosen to be odd or even, the null solution consistently acts as an attractor, and the overall magnitude of instabilities increases as the TW-speed decreases. Hence, this particular behaviour in the solutions is not dependent on the specific choice of m. Note that these additional results are not included in this work to maintain conciseness. Figure 5. Solution profiles for low values of TW-speed. Note that for TW-speed values close to zero, the null solution does not behave as an attractor. The increasing of the TW-speed (up to a moderate value of 0.1) stabilizes the attractor behaviour of the null critical point. Figure 6. For increasing values in the TW-speed, it is possible to conclude on the attractor properties of the null solution. EJDE-2024/68 INSTABILITY OF ENERGY SOLUTIONS 17 5. Scaling invariance and symmetry We continue our computation of solutions by considering symmetries in the equa- tion (1.1). Lemma 5.1. The equation (1.1) possesses self-similar solutions under the scaling transformation u(x, t) → λu(λbx, λct) with the exponents b = m− (p− 1) 2(m+ 1) , c = 1− p. Under this transformation, the original equation remains invariant, ensuring the existence of self-similar solutions. Proof. We start by considering the scaling transformation: u(x, t) → λu(λbx, λct). Applying this transformation to each term in the equation, we obtain ut → λ ∂ ∂t u(λbx, λct) = λλ−c ∂ ∂(λct) u(λbx, λct) = λ1−cut, ∆u → λ1−2b∆u, −∆(|∆u|m∆u) → λ1+m(1−2b)−2b∆(|∆u|m∆u), |u|p−1u → λp|u|p−1u. To ensure the equation is invariant under the scaling transformation, the exponents of λ must be equal for each term. This leads to the system of equations 1− c = p, 1− c = 1 +m(1− 2b)− 2b. From the first equation, c = 1− p. Then substituting into the second equation, we have 1− (1− p) = 1 +m(1− 2b)− 2b, and p = 1 +m(1− 2b)− 2b, Simplifying p− 1 = m(1− 2b)− 2b, p− 1 = m− b(2m+ 2), b(2m+ 2) = m− (p− 1), b = m− (p− 1) 2(m+ 1) . Thus, we have b = m− (p− 1) 2(m+ 1) , c = 1− p. This completes the proof that the scaling transformation leaves the original equation invariant, hence ensuring the existence of self-similar solutions. □ 18 J. L. DÍAZ PALENCIA EJDE-2024/68 Next, we construct the self-similar solutions. We assume a self-similar form u(x, t) = (T − t)−αf( x (T − t)β ), where T is the blow-up time, and α and β are constants to be determined. By comparing this form with the scaling transformation, we set α = 1 c , β = b c . From c = 1− p, we have α = 1 1− p , β = m−(p−1) 2(m+1) 1− p = m− (p− 1) 2(m+ 1)(1− p) . We now compute the derivatives. ut = α(T − t)−α−1f(ξ) + βξ(T − t)−α−1f ′(ξ), uxx = (T − t)−α−2βf ′′(ξ). For the nonlinear term, we have |uxx|muxx = (T − t)(−α−2β)(m+1)|f ′′(ξ)|mf ′′(ξ), ∆(|uxx|muxx) = (T − t)(−α−2β)(m+1)[|f ′′(ξ)|mf ′′(ξ)]′′. Substituting these into the PDE and equating the powers of (T − t), so that we derive a proper equation in terms of the single variable ξ, αf(ξ) + βξf ′(ξ) = [|f ′′(ξ)|mf ′′(ξ)]′′ + |f(ξ)|p−1f(ξ). Now, we introduce a resolution for this equation in the proximity of f → 0+, hence let us consider the expression βξf ′(ξ) = [|f ′′(ξ)|m f ′′(ξ)]′′. This equation was solved based on a numerical procedure using finite difference methods and for m = 3 and p = 2. The domain [0, ξmax] was discretized into N = 100 points, and central finite differences were employed to approximate the first, second, and fourth derivatives of the function f(ξ). An initial guess, of a Gaussian-like function (f0(ξ) = exp(−ξ2)), was used to start the iterative process. The system of nonlinear equations was then solved using the fsolve function from the scipy.optimize library in Python. The numerical solution ob- tained was subsequently compared with an assumed analytical solution of the form fanalytical = Ce−aξ, and the parameters C and a were optimized to minimize the mean absolute error between the numerical and analytical solutions, leading to the following fitted parameters, C ≈ 1.208, a ≈ 1.258. To achieve this, we minimize the difference between the numerical solution and the analytical form fanalytical = Ce−aξ by solving the optimization problem min C,a N∑ i=1 (fnumerical(ξi)− Ce−aξi)2. EJDE-2024/68 INSTABILITY OF ENERGY SOLUTIONS 19 The least squares method involves finding the parameters C and a that minimize the sum of the squared differences between the numerical solution fnumerical(ξi) and the analytical solution Ce−aξi over all discretized points ξi. To implement this, we use the curve fit function from the scipy.optimize library in Python, which em- ploys non-linear least squares to fit the exponential model to the numerical data. The curve fit function returns the optimal values of C and a that minimize the squared error. For this, we started with an initial guess for the parameters C = 1 and a = 1. Hence, we computed the error function that was defined as the sum of the squared differences between the numerical solution and the analytical form. This error function is what the optimization algorithm seeks to minimize. Par- ticularly, the curve fit function applies an optimization algorithm to adjust the parameters C and a iteratively. In each iteration, the algorithm evaluated the error function and updated the parameters to reduce the error. The process continued until the algorithm converges to a solution where the parameters C and a yield the minimum error. The resulting fitted parameters provide an exponential func- tion that closely approximates the numerical solution, as illustrated in the Figure 5. The comparison demonstrates that the fitted analytical solution matches the numerical solution with a mean absolute error of approximately 0.0238, indicating a high degree of accuracy. Figure 7. Comparison of Numerical Solution and Fitted Ana- lytical Solution for m = 3 and p = 2. The mean absolute error between for solutions is of approximately 0.0238 showing a high degree of accuracy specially in the asymptotic with ξ ≫ 1. Once the self-similar profile f(ξ) is obtained numerically, we can reconstruct the original solution u(x, t) using the self-similar form. If we consider the cases ofm = 3 and p = 2, the solution is not expected to blow-up in finite time as α = −1. To illustrate this process, we assume for simplicity that T = 1. Using the numerical solution for f(ξ), shown in Figure 5, we can reconstruct the original solution u(x, t) as provided in Figure 5. 20 J. L. DÍAZ PALENCIA EJDE-2024/68 Figure 8. Reconstructed Solution u(x, t) for Different Times 6. Conclusions The analysis followed in the presented study has permitted to introduce results on instabilities of TWs for a kind of problem with a super-linear reaction and higher order p-Laplacian operator. The analysis started by the introduction of appropriate functional spaces to account for the particular oscillating behaviour of solutions along with compact support properties. The regularity of the involved operator has been shown in the defined functional spaces concluding that any oscillating solution (in the normHΘ) can be bounded by mollified solutions (through the normHm ρ ) and compact support solutions (norm Hm 0 ). Afterward, the boundedness of solutions was provided making use of the defined norms and considering energy solutions. Once the regularity results were presented, the problem (1.1) was analyzed in the TW domain making use of different lemmas to show the permanent instabilities of solutions. In the TWs analysis, a special emphasis was set in the null critical poin. A numerical assessment was introduced to provide the homotopy graphs for a wide interval of TW-speeds, along with the solutions exact profiles. This approach permits to validate our analytical assessments and our first postulated intuition: This is that the null solution acted as an attractor to any oscillating flow, leading potentially to hinder blow-up formation. Eventually, we introduced the scaling invariant properties of the equation and obtained self-similar solutions. References [1] Adams, R. A.; Anisotropic Sobolev inequalities. Casopis pro pestovan’i matematiky. 113.3 (1988): 267-279. http://eudml.org/doc/19616 [2] Ahn, J.; Yoon, C.; Global well-posedness and stability of constant equilibria in para- bolic–elliptic chemotaxis system without gradient sensing. Nonlinearity, 32 (2019), 1327-1351. [3] Alexander, J.; Gardner, R.; Jones, C.; A topological invariant arising in the sta- bility analysis of travelling waves. J. Reine Angew. Math. 410 (1990), 167–212. DOI 10.1515/crll.1990.410.167 [4] Aronson, D.; Density-dependent interaction-diffusion systems. Proc. Adv. Seminar on Dy- namics and Modeling of Reactive System, Academic Press, New York, 1980. EJDE-2024/68 INSTABILITY OF ENERGY SOLUTIONS 21 [5] Aronson, D.; Weinberger, H.; Nonlinear diffusion in population genetics, combustion and nerve propagation. Partial Differential Equations and Related Topics. (1975) Pub., New York, 5–49. [6] Aronson, D.; Weinberger, H.; Multidimensional nonlinear diffusion arising in population ge- netics. Adv. in Math. 30 (1978), 33–76. [7] Audrito, A.; Vázquez, J. L.; The Fisher–KPP problem with doubly nonlinear “fast” diffusion, Nonlinear Analysis, 157 (2017), 212-248. [8] Benedek, A.; Panzone, R.; The spaces Lp with mixed norm. Duke Math. J. 28 (1961), 301-324. [9] Bhatti, M.; Zeeshan, A.; Ellahi,R.; Anwar B’eg, O.; Kadir, A.; Effects of coagulation on the two-phase peristaltic pumping of magnetized prandtl biofluid through an endoscopic annular geometry containing a porous medium, Chin. J. Phys. 58 (2019), 222-23. DOI 10.1016/j.cjph.2019.02.004. [10] Bognar, G.; Numerical and Analytic Investigation of Some Nonlinear Problems in Fluid Mechanics. Comp. and sim. in modern sci. Vol.II (2008), pp.172-179. [11] Bonheure, D.; Sánchez, L.; Heteroclinics Orbits for some classes of second and fourth order differential equations. Handbook of differential equations. 3(06) (2006), 103-202. [12] Carelman, T.; Problemes mathematiques dans la theorie cinetique de gas, AlmquistWiksells, Uppsala, 1957. [13] Cho, E.; Kim, Y. J.; Starvation driven diffusion as a survival strategy of biological organisms Bull. Math. Biol., 75 (2013), 845-870. [14] Cohen, D. S.; Murray, J. D.; A generalized diffusion model for growth and dispersal in a population. J. Math. Biology 12 (1981), 237–249. DOI 10.1007/BF00276132 [15] Coutsias, E. A.; Some effects of spatial nonuniformities in chemically reacting systems. Cal- ifornia Institute of Technology, 1980. [16] Dee, G. T.; Van Sarloos, W.; Bistable systems with propagating fronts leading to pattern formation. Physical Review Letter, Volume 60, 1998. [17] Dı́az, J. L.; (2022) Non-Lipschitz heterogeneous reaction with a p-Laplacian operator. AIMS Mathematics, 7(3) (2022): 3395-3417. DOI 10.3934/math.2022189 [18] Durham, A. C.; Ridgway, E. B.; Control of chemotaxis in physarum polycephalum. J. Cell. Biol. 69 (1976), 218–223. DOI 10.1083/jcb.69.1.218. [19] Egorov, Y.; Galaktionov, V.; Kondratiev, V.; Pohozaev, S.; Global solutions of higher-order semilinear parabolic equations in the supercritical range. Adv. Differ. Equat. 9 (2004), 1009- 1038. [20] Ellahi, R.; Hussain, F.; Ishtiaq, F.; et al.; Peristaltic transport of Jeffrey fluid in a rectangular duct through a porous medium under the effect of partial slip: An application to upgrade industrial sieves/filters. Pramana - J Phys 93 (2019), 34. DOI 10.1007/s12043-019-1781-8 [21] Enright, H.; Muir, P. H.; A Runge-Kutta type boundary value ODE solver with defect control. Teh. Rep. 267/93, University of Toronto, Dept. of Computer Sciences. Toronto. Canada, 1993. [22] Fisher, R. A.; The advance of advantageous genes. Ann. Eugenics, 7 (1937), 355–369. [23] Galaktionov, V. A.; Three types of self-similar blow-up for the fourth order p-Laplacian equation with source. Jour. Comp. and Appl. Math., 223 (2009), 326-355. [24] Galaktionov, V. A.; On a spectrum of blow-up patterns for a higher-order semilinear para- bolic equation. Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences. 2001. [25] Galaktionov, V.; Towards the KPP–Problem and Log-Front Shift for higher-Order Nonlinear PDEs I. Bi-Harmonic and Other Parabolic Equations. Cornwell University arXiv:1210.3513, 2012. [26] Galaktionov, V.; Shishkov, A.; Saint-Venant’s principle in blow-up for higher-order quasilin- ear parabolic equations. Proceedings of the Royal Society of Edinburgh: Section A Mathe- matics, 133(5) (2003), 1075-1119. DOI 10.1017/S0308210500002821. [27] Goldshtein, V.; Ukhlov, A.; Weighted Sobolev Spaces and embeddings Theorems. Transac- tions of the american mathematical society. 361 (2009), 3829-3850. [28] Hongjun, G.; Changchun, L.; Instabilities of traveling waves of the convective-diffusive Cahn- Hilliard equation. Chaos, Solitons and Fractals, 20 (2004), 253-258. [29] Kamin, S.; Vázquez, J. L.; Fundamental Solutions and Asymptotic Behaviour for the p- Laplacian Equation. Revista Matemática Iberoamericana. Vol 4 (1988), 2. [30] Keller, E. F.; Segel, L. A.; Traveling bands of chemotactic bacteria: a theoretical analysis. J. Theoret. Biol. 30 (1971), 235-248. 22 J. L. DÍAZ PALENCIA EJDE-2024/68 [31] Kesavan, S.; Topics in Functional Analysis and Applications, New Age International (for- merly Wiley-Eastern). 1989. [32] Kolmogorov, A .N.; Petrovskii, I. G.; Piskunov, N. S.; Study of the diffusion equation with growth of the quantity of matter and its application to a biological problem. Byull. Moskov. Gos. Univ., Sect. A, 1 (1937). [33] Ladyzhenskaya, O.; Some results on modifications of three-dimensional Navier-Stokes equa- tions. Nonlinear Analysis and Continuum Mechanics. (1998) pp. 73–84. [34] Li, Zhenbang; Liu, Changchun; On the nonlinear instability of traveling waves for a sixth- order parabolic equation. Abstr. Appl. Anal. Article ID 739156 (2012), 17 pages. DOI 10.1155/2012/739156 [35] Montaru, A.; Wellposedness and regularity for a degenerate parabolic equation arising in a model of chemotaxis with nonlinear sensitivity. Discrete and Continuous Dynamical Systems - B, 19,1 (2014), 231-256. [36] Niemela, J. J.; Ahlers, G.; Cannell, D. S.; Localized traveling-wave states in binary-fluid convection. Phys. Rev. Lett. 64 (1990), 1365–1368. DOI 10.1103/PhysRevLett.64.1365. [37] Okubo, A.; Levin, S. A.; The Basics of Diffusion. In: Diffusion and Ecological Problems: Modern Perspectives. Interdisciplinary Applied Mathematics, vol 14. Springer, New York, NY, 2001. DOI 10.1007/978-1-4757-4978-6. [38] Palencia, J. L. D.; Analysis of selfsimilar solutions and a comparison principle for an hetero- geneous diffusion cooperative system with advection and non-linear reaction. Comp. Appl. Math. 40 (2021), 302. DOI 10.1007/s40314-021-01689-y [39] Peletier, L. A.; Troy, W. C.; Spatial Patterns. higher order models in Physics and Mechanics. Progress in non linear differential equations and their applications. Volume 45. Universit’e Pierre et Marie Curie, 2001. [40] Rauprich, O.; Matsushita, M.; Weijer, C. J.; Siegert, F.; Esipov, S. E.; Shapiro, J. A.; Periodic phenomena in proteus mirabilis swarm colony development. J. Bacteriol. 178 (1996), 6525–6538. DOI 10.1128/jb.178.22.6525-6538.1996 [41] Rottschäfer, V.; Doelman, A.; On the transition from the Ginzburg-Landau equation to the extended Fisher-Kolmogorov equation. Physica D. 118 (1998), 261-292. [42] Shishkov, A. E.; Dead cores and instantaneous compactification of the supports of energy solutions of quasilinear parabolic equations at arbitrary order. Sb. Math., 190 (1999), 1843- 1869. [43] Strauss, W.; Wang, G.; Instabilities of travelling waves of the Kuramoto-Sivashinsky equation. Chin Ann Math B. 23 (2002), 267-276. [44] Tao, Y.; Winkler, M.; Effects of signal-dependent motilities in a keller–segel-type reaction- diffusion system. Math. Models Methods Appl. Sci., 27 (2017), pp. 16-45 [45] Yoon, C.; Kim, Y. J.; Global existence and aggregation in a keller–segel model with fokker- Planck diffusion. Acta Appl. Math., 149 (2016), pp. 101. José Luis D́ıaz Palencia Department of Mathematics and Education, Universidad a Distancia de Madrid, 28400 Madrid, Spain Email address: joseluis.diaz.p@udima.es 1. Problem description and objectives 2. Preliminaries 3. Boundedness of solutions 4. Travelling waves 4.1. Travelling wave instabilities 4.2. Exact travelling wave profiles and characteristic propagation speed 5. Scaling invariance and symmetry 6. Conclusions References