New Analytical Solutions for Time-Fractional Stochastic (3+1)-Dimensional Equations for Fluids with Gas Bubbles and Hydrodynamics

: This paper explores the effects of spatial fractional derivatives and the multiplicative Wiener process on the analytical solutions for (3+1)-dimensional fractional stochastic equations for ﬂuids with gas bubbles. We study the bifurcation of the analytical solutions and introduce new fractional stochastic solutions. We also discuss how the solutions differ depending on the initial conditions. The new solutions are notably more beneﬁcial and impactful for understanding various, signiﬁcant, and incredibly hard physical phenomena due to the signiﬁcance of the modiﬁed fractional stochastic (3+1)-dimensional equations for ﬂuids with gas bubbles and hydrodynamics. We also discuss the effects of the fractional order and the Wiener process on the obtained analytical solutions.


Introduction
Partial differential equations are used to model the most complicated natural phenomena in diverse branches of science, including fluid dynamics, condensed matter physics, optical fibers, geochemistry, mathematical biology, oceanography, chemical physics, etc.By analyzing the exact solutions, we are able to obtain insights into the specific physical systems that underly the dynamical behavior of these equations.Therefore, finding precise solutions to these nonlinear mathematical physical equations is essential for comprehending nonlinear phenomena.Nevertheless, there are no unified methods that can provide all the exact solutions to nonlinear partial differential equations because nonlinear partial differential equations have different forms and characteristics that present distinct challenges when trying to detect their exact analytical solutions.Numerous efficient methods have been reported that can find exact solutions using diverse properties.For instance, inverse scattering methods [1-3], Hirota's bilinear operators [4][5][6], the variational iteration method [7], Lie group methods [8,9], the homogeneous balance method [10,11], the first integral method [12,13], and bifurcation methods [14][15][16][17][18] have been used, among others that are outlined in the references of those studies.Recently, partial differential equations have been extended or generalized in three directions: fractional partial differential equations, which take into account the fractional order instead of the integer order; stochastic partial differential equations, which take into account random effects in the modeling of complicated phenomena; fractional stochastic partial differential equations, which take into account both fractional order and random effects.Many branches of science, such as engineering and physics, have focused on the importance of taking random effects into account during the analysis and modeling of complex phenomena, as the noise effect may explain more important statistical properties of those phenomena that cannot be neglected [19][20][21][22][23].
It is worth mentioning that when stochastic terms are introduced into partial differential equations, the analytical solutions are more difficult to derive than those of deterministic partial differential equations.
Scientists have become increasingly interested in liquids with gas bubbles, which are commonly found in nature, science, engineering, and physics [24][25][26].Bubble-liquid mixture equations have been found to weakly expound nonlinear waves that propagate in liquids with gas bubbles [25,26].In the current work, we consider a space fractional stochastic (3+1)-dimensional wave equation for fluids with gas bubbles with the following form: where U is the mixture velocity, which is a function of the scaled spatial coordinates x, y, z, characterizing the transverse y and z perturbations of waves propagating in the x direction; t is the temporal coordinate; D α indicates the conformable fractional-order derivative of the fractional order α, where 0 < α ≤ 1; σ refers to the noise intensity; G(t) is the standard Brownian motion.According to the authors of [27], the enclosed parameters a i are interpreted as a 1 , a 2 , and a 3 , which stand for nonlinearity, dispersion, and perturbed effect, respectively, while a 4 and a 5 indicate the troubling wave velocities along the y and z directions, respectively.

Motivations
We had many motivations to examine Equation (1), including the following: (a) Due to the difficulties in considering every aspect of the problem, it is useful to consider the stochastic perturbations of the nonlinear (3+1)-dimensional wave equation that describes liquids with gas bubbles, which takes the following form [28] of Equation (1) with space fractional derivatives.Hence, Equation (1) provides a good description for the liquids with gas bubbles; (b) Equation (1) has not previously been considered in the literature.Consequently, it is a good model for investigation.Consequently, the obtained results are new; (c) Equation ( 1) is considered as an extension to Equation (2) by taking into account the effect of fractional order only, random effects, and the combined influence of both fractional order and random effects.Hence, the obtained results for Equation (1) can be employed to restore previous results or introduce new solutions to the equations, which can be obtained as special case forms of Equation (1) (as outlined in Table 1).(1, 1, 0, 0, 0) (2+1)-dimensional Kadomtsev-Petviashvili [32] with constant coefficients 3.
(6, 1, 0, 0.5, 0.5) (6, 1, 0, 0, 0) (1+1)-dimensional KdV equation that describes [1] the fluid dynamics in some complicated phenomena, such as plasma physics and lattice dynamics This work is prepared as follows.Section 2 introduces some facts about fractional derivatives and the Wiener process that are applied throughout the paper.Section 3 presents a mathematical analysis of Equation (1) to reduce it to a planar dynamic system, as well as the application of the qualitative theory to the dynamics system, including the properties of the phase portrait.Section 4 offers some new solutions to the fractional stochastic equation in Equation (1).Section 5 considers the influence of the Brownian motion, fractional-order derivatives, and their combined effects on the derived solutions.Section 6 summarizes the obtained results.

Preliminaries
Since this work deals with time-fractional stochastic partial differential equations, we present some facts about fractional derivatives and the Wiener process in order to produce a self-contained article.
The application of fractional partial differential equation models has increased in many research fields over recent years, including fluid mechanics, plasma, nuclear physics, the mechanics of materials, biology, chemistry and finance.
To model and best analyze real-world problems, some appropriate operators had to be proposed.Recently, novel operators of fractional calculus have been proposed and developed, such as Riesz, Liouville, Hilfer, and conformable fractional.According to the definitions of the operators of conformable fractional, some restrictions of other fractional operators are resolved and basic rules of classical calculus, such as the chain rule, the product of two functions and the mean value theorem have been preserved.It is very simple and efficient to apply conformable derivatives, and it improves the study of how physical phenomena behave [34].

Definition 1 ([34]
).The conformable fractional derivative of order ε for the function P : R + → R is defined, for any t ∈ R + and 0 < ε ≤ 1, as Next, we state some of basic rules of conformable derivatives that are helpful in the current work, assuming that P 1 and P 2 are ε− conformable differential functions for t > 0 and that c and d are two real constants: , where P 2 (t) is included in the range of P 1 .
Next, we present the definition of one-dimensional Brownian motion [35].
Definition 2. A one-dimensional continuous time stochastic process {J(t)} t≥0 is termed Brownian motion if the following are verified:

Mathematical Analysis
The previous studies concerning the construction of wave solutions to Equation (2) assumed that the solutions were in the form where b i (i = 1, 2, 3) and c are constants (see, e.g., [36]).Solution (4) usually aims to convert Equation ( 2) into an ordinary differential equation.After extending Equation (2) into the fractional stochastic partial differential Equation ( 1), the solutions must also be imposed in extended forms, i.e., in the form Notice that when σ = 0, α → 1, Solutions ( 4) and (5) become identical.Hence, the fractional stochastic forms of Solution (5) can be reduced to solutions for Equation (1) when σ = 0, α → 1.After some calculations, we obtain where v = dv dζ .Substituting the expressions in (6) into Solution (5), we obtain According to our assumption, G(t) is normally distributed, so E(e σG(t) ) = e 1 2 σ 2 t .Hence, the expectation of both sides of Equation ( 7) leads to The double integration of both sides with respect to ζ implies where γ and β are introduced instead of the original parameters for convenience.These are given by The dynamical system corresponding to Equation ( 9) System ( 11) is conservative, since div(v , w ) = 0 and Hamiltonian, where The Hamiltonian system corresponding to Equation ( 12) is a conserved quantity, i.e., where s is a constant that represents the total energy of the system.The equilibrium points for System (11) are O(0, 0) and E −γ β , 0 .The qualitative theory of dynamical systems [37] states that when γ > 0, the equilibrium points O and E are saddle and center points, respectively, and when γ < 0, they are converse points.The phase portraits of System (11) are shown in Figures 1 and 2. The total energy at the equilibrium points is denoted by The rest of this paper is devoted to analyzing System (11), as the solutions to this system lead to the solutions of Equation (1).We begin with the following clear result.According to Figures 1 and 2, the phase portraits of System (11) can be summarized by the following proposition.Proposition 1.For γ > 0 and β > 0, it follows: (a) For s < s 1 , System (11) has a family of unbounded orbits, as illustrated by the brown lines in Figure 1a; (b) For s = s 1 , System (11) has an unbounded orbit that separates the preceding family and the succeeding family, as shown by the black lines in Figure 1a; (c) For s 1 < s < s 2 , System (11) has two families of orbits: a family of periodic orbits around the center E and a family of unbounded orbits to the right of the saddle O, as shown by the green lines in Figure 1a; (d) For s = s 2 , System (11) has a homoclinic orbit passing through the saddle O and surrounding the family of periodic orbits and two unbounded orbits passing through the saddle O, as shown by the red lines in Figure 1a; (e) For s > s 2 , System (11) has a family of unbounded orbits outside the homoclinic orbit, as shown by the blue lines in Figure 1a.
The phase portraits illustrated in Figures 1b and 2a,b for the other cases of γ and β can be similarly described using Lemma 2; consequently, we do not consider them here.

Solution Construction
In this section, we derive some novel solutions for Equation (1).We only focus on establishing wave solutions when γ > 0 and β > 0 and do not study any other cases due to their similarities (either in terms of signs, displacements, or both), according to Lemma 2. By substituting the first equation from System (11) into the conserved quantity, we obtain where • For s 2 < s < s 1 , there is a family of orbits consisting of periodic orbits and unbounded orbits, as shown by the green lines in Figure 1a.Since all orbits intersect the v− axis at three different points, F (v) takes the form , where v 0 < v 1 < 0 and v 2 > 0. We then focus on the regions that permit real propagation, i.e., the shaded regions in Figure 3a that are restricted by ]v 0 , v 1 [ and ]v 2 , ∞[.Therefore, for Equation (15), using (v, w) ∈ R × R : 3 2 w 2 = F (v), v 0 < v < v 1 , produces the following periodic solution: , produces the following singular solution: Consequently, Equation (1) has the following solutions: and • For s = s 1 , there is a homoclinic orbit passing through the equilibrium point O around E and two unbounded orbits passing through the equilibrium point O, as shown by the red lines in Figure 1a.These orbits intersect the v− axis in two points; therefore, F (v) has two real zeros: one is single (e.g., v 3 ) and the other (at O) is double.Hence, F (v) = βv 2 (v − v 3 ).We then focus on the regions that permit real propagation, i.e., the shaded regions in Figure 3b that are restricted by ]v 3 , 0[ and ]0, ∞[.Therefore, for Equation ( 15), using (v, w) ∈ R × R : 3 2 w 2 = F (v), v 3 < v < 0 produces the following one-soliton solution: Consequently, Equation ( 1) has the following solutions: and • For s > s 1 , there is a family of unbounded orbits outside the homoclinic orbit and its two unbounded orbits, as shown by the blue lines in Figure 1a.All of these orbits intersect the v− axis on the left-hand side of equilibrium point E; therefore, F (v) has one real zero v 4 < 0 and two conjugate complex pair zeros z and z.Hence, We then focus on the region that permit real propagation, i.e., the shaded regions in Figure 3c that are restricted by ]v 4 , ∞[.Hence, for Equation (15), produces the following solution: where A 2 = (e 1 − v 14 ) 2 + e 2 ; e 1 = e and e 2 = e.Consequently, Equation (1) has the following solution:

Physical Interpretation
This section aims to explore the effects of noise intensity and fractional derivative order on the solutions, both individually and combined.We demonstrate these effects on the two most important types of solutions, namely Solutions ( 19) and (23), which reduce to the periodic and the solitary solutions in the deterministic state (i.e., when σ = 0 and α = 1).It is worth mentioning that the selected values for the β, γ, and s parameters are based on the bifurcation analysis illustrated in Section 3. The next graphic representations are based on finite-difference methods, with the help of Mathematica software.

Deviation from the periodic solution
Let us assume that β = 1, γ = 2, and (a) The noise effect.Figure 4a shows the 2D representation of Solution (27) for diverse noise intensity σ, where the fractional differential order α is 1.The disturbance around the blue line, which represents the deterministic state, increases with an increase in noise strength.Note that the wave amplitude increases while the wavelength remains almost the same.It is evident from Figure 4b that the 3D solution represents a smooth periodic surface in the deterministic state (σ = 0) but becomes a rough surface in the non-deterministic status.Moreover, it loses its periodicity for larger values of noise intensity, so the surface becomes flat.(b) The fractional order effect.Figure 5a illustrates the impact of the fractional order α in the deterministic status (σ = 0) through the deviations surrounding the the blue line.
The wavelength of Solution (27) increases as the fractional order α moves away from 1, while the amplitude remains almost the same.The corresponding 3D representation of the solution is shown in Figure 5b.
(c) The combined influence.Figure 6a shows that the amplitude and wavelength of Solution ( 27) increase when both the noise intensity and the fractional order move away from 0 and 1, respectively.The surface characterized by Solution (27) that is smooth for σ = 0, γ → 1 becomes a rough surface, as shown in Figure 6b.This discussion can be summarized by the following lemma.
Lemma 3.For Solution (19), the further σ is from 0, the noisier the curve.On the other hand, the further σ is from 1 and the closer it is to 0, the more the curve loses its periodicity by gradually increasing the wavelength.

Deviation from the soliton solution
For the values β = 1 and γ = 2, we choose s = 0. Hence, Solution (23), becomes (a) The noise effect.Figure 7a illustrates the 2D representation of Solution (28) for different noise intensity values when the fractional order α = 1.The wave height increases while the width remains almost the same as the noise intensity increases.Figure 7b shows the surface characterized by this solution, which loses its smoothness as the noise intensity increases.Further, as the noise intensity increases, the surface flattens.(b) The fractional order effect.Figure 8a shows the effect of the fractional order derivative α on Solution (28) when σ = 0.As the fractional order moves away from 1, Solution (28) gradually loses part of its wave structure and smoothness.The 3D representation is shown in Figure 8b.(c) The combined influence.Figure 9a shows the combined effect of σ and α.Solution (28) loses part of its wave structure and smoothness.The same conclusion can be drawn for the corresponding 3D representation and the surface gradually becomes rougher, as shown in Figure 9b.The following lemma summarizes this discussion.
Lemma 4. For solution (23), the further σ is from 0, the noisier the curve.On the other hand, the further σ is from 1 and the closer it is to 0, the more the curve loses part of its structure and gradually loses its smoothness.

Conclusions
The current work aimed to explore the impact of spatial fractional derivatives and the multiplicative Wiener process on the analytical solutions for (3+1)-dimensional fractional stochastic equations for fluids with gas bubbles.The bifurcation of the solutions was studied and some new solutions were introduced.These solutions were also new to the deterministic version of Equation ( 1), with derivatives that were either integers or fractions.Additionally, they were new for the stochastic integer-order derivative version of Equation (1).We concentrated on only forming real solutions for particular regions, which are called regions of real propagation.This type of solution is significant because it is applicable to real-world scenarios.On the other hand, the obtained exact solutions were significant because of their ability to describe diverse crucial challenges in weakly nonlinear waves propagating in liquids with gas bubbles, which are described by stochastic Equation (1).We graphically illustrated the individual impact of noise intensity and fractional order derivative, as well as the combined effect of both.These effects could be clearly observed in the amplitude, wavelength, and periodicity of the solutions.Furthermore, the smooth surface represented by the solution in the deterministic case became rough due to the noise term and lost its periodicity due to the fractional-order term.

Figure 4 .
Figure 4.The effect of noise intensity on Solution (27) when the fractional order α is 1.

Figure 5 .Figure 6 .
Figure 5.The effect of the fractional order derivative α on Solution (27) in the deterministic case.

Figure 7 .
Figure 7.The effect of noise intensity on Solution (28) when the fractional order α is one.

Figure 8 .
Figure 8.The effect of the fractional order derivative α on Solution (28) in the deterministic case.

Figure 9 .
Figure 9.The combined effect of noise intensity and fractional order derivative on Solution (28).