Transient Propagation of Longitudinal and Transverse Waves in Cancellous Bone: Application of Biot Theory and Fractional Calculus

: In this paper, the inﬂuence of the transverse wave on sound propagation in a porous medium with a ﬂexible structure is considered. The study is carried out in the time domain using the modiﬁed Biot theory obtained by the symmetry of the Lagrangian (invariance by translation and rotation). The viscous exchanges between the ﬂuid and the structure are described by fractional calculus. When a sound pulse arrives at normal incidence on a porous material with a ﬂexible structure, the transverse waves interfere with the longitudinal waves during propagation because of the viscous interactions that appear between the ﬂuid and the structure. By performing a calculation in the Laplace domain, the reﬂection and transmission operators are derived. Their time domain expressions depend on the Green functions of the longitudinal and transverse waves. In order to study the effects of the transverse wave on the transmitted longitudinal waves, numerical simulations of the transmitted waves in the time domain by varying the characteristic parameters of the medium are realized whether the transverse wave is considered or not.


Introduction
Osteoporosis is defined as a sickness that affects the bones, and to understand it, an overview of the structure of bone tissues is necessary. Bones are in perpetual renewal; this phenomenon is called bone remodeling: constant construction of new tissue and destruction of old tissue. In osteoporosis, biochemical and hormonal changes occur and affect this balance. The bone is, therefore, considerably weakened by reductions in its tissue and mass [1], a fragility that is increased tenfold when the contribution of bone minerals is reduced [2][3][4]. Therefore, there is a growing need for improved diagnosis and management of osteoporosis. Indeed, the sooner the disease is detected, the more effective the care can be.
As it is known, bone is an inhomogeneous porous material [5][6][7][8][9] and the determination of its characteristics by means of ultrasound techniques can be a game-changer in the diagnosis of osteoporosis since its elastic properties affect acoustic propagation. Within this context, several studies have been conducted on cancelous and trabecular bone [10][11][12][13]. There are ultrasonic methods that rely on the speed of sound and the broadband attenuation of a single wave to estimate the different stages of osteoporosis, but the poroelastic nature of the bone generates two types of waves during the propagation of the acoustic wave in the bone. The most appropriate theory for the study of acoustic wave propagation in porous media such as trabecular bone is the Biot theory [14][15][16][17][18][19][20][21][22][23][24][25][26][27]. This theory [14][15][16][17][18][19][20][21][22][23][24][25][26][27] has been proven to be a very advantageous model for the description of wave propagation in porous media saturated with fluids; it takes into account the movements of the solid and fluid phases, and through the mass, elastic, and viscous couplings, the interaction between these two phases is considered. Its main success lies in the prediction of the existence of three waves in a porous material subjected to external stress: two longitudinal waves and a transverse wave. The first compressional wave has the highest speed and is called the "fast wave", while the second compressional wave has a lower speed and is called the "slow wave".
Because of its very general and rather fundamental character, the Biot theory has been used widely in the field of petroleum exploration and geoscientific testing, but many have also applied it to the study of bone structures [14,15,[28][29][30][31]. McKelvie [28,31] qualitatively predicted that, for cancelous bone, attenuation is related to ultrasound frequency. The second compressional wave with a lower velocity was reported in the studies by Lakes, Yoon, and Katz [32] on wet human and bovine cortical bone. The connection between porosity and longitudinal wave velocities and the relevance of the presence of these two waves to bone characterization were demonstrated by Hosokawa and Otani [14] in their studies of cancelous bone. By implementing the concepts of dynamic tortuosity, permeability, and viscous characteristic length, in 1987, Johnson et al. [33] made a major input to the study of porous media.
The modified Biot theory has been used to study the transmission and reflection of acoustic waves for different media in the frequency domain [34][35][36][37][38][39][40][41][42]. In order to offer an alternative to the frequency methods, Fellah and Depollier [43,44] have proposed a temporal approach [43][44][45][46][47][48][49][50][51][52] that allows the phenomena that follow the propagation of acoustic waves in porous media to be seen from a different angle. The viscous fluid/structure interactions are described using fractional calculus, which has been proven to be a well-adapted tool to describe the dispersion effects that develop when a wave passes in a saturated porous medium. Fractional derivatives and integrals [53][54][55][56][57][58][59][60][61][62][63][64][65] have long been the focus of research in mathematics, but their use outside mathematics is only recent. Indeed, the fractional derivative was only used towards the end of the 1960s by Caputo [57][58][59] to describe dissipation phenomena in solids and, later, by Bagley and Torvik [60] to treat viscoelastic materials, but the most significant contribution in the field of acoustics remains the one made by Matignon [61]. The reflection and transmission operators have been established as a result to external solicitation of a saturated porous medium by a viscous fluid in the frequency domain for normal incidence by considering the three Biot waves instead of only the longitudinal waves, by Hodaei, Maghoul, and Popplewell [66].
The purpose of this article is to study, from a new perspective, the behavior of a porous medium subjected to an external load when the three Biot waves (i.e., the transverse wave in addition to the longitudinal waves) are considered. For that, we establish the expressions of the diffusion operators in the time domain while including the effect of the transverse wave; the fractional calculus as well as the Biot model modified by Johnson et al. [33] will be used. Numerical simulations are given for the transmitted signals, with shear being taken into account or not while applying a variation of the characteristic parameters of the porous medium. The effect of transverse waves as well as the sensitivity of the parameters are studied. The advantage of this model is that it presents a temporal modeling of the propagation by taking into account both longitudinal and transverse waves in normal incidence.

Modified Biot Theory
The displacements of the solid phase u and the fluid phase U are modeled via the Biot theory, where the equations of motion are given by the following: N is the shear modulus of the solid structure. P, Q, and R are the elastic coefficients of Biot; they are functions of the bulk modulus of the pore fluid K f , the bulk modulus of the elastic solid K s , the bulk modulus of the porous skeletal frame K b as well as of the porosity φ of the porous medium and are defined by the relations below: K f , K s , and K b depend on the classical mechanical parameters, which are Young's moduli, Poisson's coefficients of the solid E s and ν s , and Poisson's coefficients of the skeleton E b and ν b as follows: The coefficients ρ ij involved in the motion equations are mass coefficients. They have a density dimension and are linked to the skeleton-constituting solid density ρ s , which constitutes the skeleton, and to the fluid density ρ f by: The term ρ 12 expresses the fluid-structure inertial coupling. It appears as a contribution to the densities of the fluid and the solid because of the interactions between the fluid and the structure.
α ∞ is the tortuosity of the medium. When an acoustic wave crosses a fluid-saturated porous medium, it is attenuated. This attenuation is mainly due to the viscous exchanges that occur between the structure and the saturating fluid because of relative fluid/solid displacements, leading to a modification of the densities. The tortuosity becomes a function of the frequency and is called dynamic tortuosity α(ω): In order to estimate the fluid domain affected by these exchanges, the thickness of the viscous boundary layer δ = 2η ωρ 0 1/2 (η fluid viscosity and ω angular frequency) is compared to a microscopic dimension characteristic of the medium such as the pore size "a". In the high-frequency mode, the viscous boundary layer thickness is small compared to the pore size "a" i.e., if ( δ a 1). In this mode, the viscous effects take place only in a small fluid thickness, at the walls where the velocities are affected by the fluid/solid displacements. In the rest of the fluid volume, the fluid behaves like a perfect fluid.
In this case, the dynamic tortuosity takes the following form: This model was proposed by Johnson et al. [33]. It involves the tortuosity α ∞ and the viscous characteristic length Λ, a geometric parameter indicative of small pores, in the high-frequency limit and represents the thickness where viscous effects are important.
By replacing the new expressions for the mass coefficients in (6) where α(ω) is given by (7) in the system of Equation (1), the motion equations becomes the following: This system is written neither in the time domain nor in the frequency domain, since the densities depend on the frequency and the derivatives are in the time domain. In the following section, we try to describe the above problem in the time domain.

Temporal Formulation of the Modified Biot Theory
In order to study the waves transmitted in the time domain through a porous medium, we focus on the temporal formulation of the Biot theory. To achieve this, we introduce the concept of the fractional derivative of order ν, the definition of which [56] is given by the following: where ν is a real number and Γ is the Gamma function. The high-frequency expression of the dynamic tortuosity includes the term 1 jω , of which the time equivalent is a fractional derivative of order 1/2. Following the definition in (9), the expression of the response factor α takes the following expression in the time domain [43,49]: where δ(t) is the Dirac function. With this new expression of the dynamic tortuosity α(t) in the system of Equation (1), the equations of motion in the time domain take the following shape: where the temporal operators ρ 11 , ρ 12 , and ρ 22 are obtained by replacing α(ω) in Equation (6) by its temporal equivalent α given by Equation (10).

Solution of the Motion Equations
According to Helmholtz's theorem, a solution to the system in (12) has the following form: where Φ s and Φ f are the scalar displacement potentials of the solid u and fluid U, respectively, and Ψ s and Ψ f are the vector potentials of u and U, respectively. The Biot equations of motion are linear, and it it is possible to look for the scalar and vector solutions separately.

Longitudinal Waves
In the case of the longitudinal wave, the contribution of scalar displacement potentials are obtained by substituting u(x, t) = ∇Φ s (x, t) and U = ∇Φ f (x, t) in the equations of motion, which gives where a = 2 . To solve this system of equations, we use the Laplace transform.
The Biot equations are then written as follows: where Φ s and Φ f are, respectively, the Laplace transform of Φ s and Φ f . We note that Φ s (x, s) is the Laplace transform of Φ s (x, t) defined by the following: The matrix M is given by the following: with The system is solved by finding the eigenvalues of the matrix M: This corresponds to the following equations: The solutions to the above are in the following form: The displacements of the solid and fluid are the linear combinations of the displacements due to the fast wave (first-type wave) characterized by the scalar potential Φ 1 (x, s) and the slow wave (second-type wave) characterized by the scalar potential Φ 2 (x, s).
The expressions for the eigenvalues λ 1 (s) and λ 2 (s)) and the components V 1 (s) and V 2 (s) of the eigenvectors are given in Appendix A.

Transverse Wave
Fo the wave, the displacements of the solid and fluid phases can be written as follows: The equations of motion in (12) are then A relation of proportionality between the solid and fluid displacement potentials is obtained from Equation (24) by where This last relation allows us to write the system of Equation (12) in the following form: Using the Laplace transform, where Ψ s (x, s) is the Laplace transform of Ψ s (x, t) as defined previously in (15) gives the following: Posing The solution in the Laplace domain takes the following form The proportionality relation in (25) in the Laplace domain becomes with and The following section focuses on the determination of the diffusion operators in the time domain; the expressions of the longitudinal and transverse waves in (19), (20), and (29); the expressions of the pressure and velocity fields in the different media; as well as the boundary conditions at the interfaces.

The Scattering Operators
In a two-dimensional problem (X-Z), a poroelastic material of thickness L, which is supposed to be homogeneous and isotropic, is saturated by a viscous fluid and found to be disturbed by a high-frequency sound wave that arrives in normal incidence on the left boundary of the material (x = 0); the geometry of the problem is represented in the Figure 1. A part of the incident wave is reflected in the fluid and the other one is transmitted through the porous medium in the shape of the three waves of Biot: the two compressional waves and a shear wave. By pointing out the fact that the porous medium occupies the region (0 ≤ x ≤ L), we can see that the expression of the pressure field in the region (I), where x < 0, is the sum of the incident p i and the reflected p r fields.
is the total pressure field at the left of the medium. In the limit x > L, the pressure field results from the transmission of the incident wave through the porous material and is given by the following: where p t is the pressure of the transmitted wave. The expressions of the reflected and transmitted pressure fields can be established by integrating the product of the scattering operators (R and T) by the incident wave, which gives the following: The fact that the lower limit of the integrals of these expressions is zero indicates that it is exactly at t = 0 that the incident wavefront hits the material for the first time.
Accordingly, p 1 (x, t) and p 3 (x, t) are expressed by the following equations: T 1 (t) and c 1 are the transmission operator and the propagation speed of the fast wave, respectively. Similarly, T 2 (t) and c 2 are the transmission operator and the propagation speed of the slow wave, respectively.
To simplify the analysis, as previously achieved, we use the Laplace transform, which is more adequate for our problem. We note that, following the definition in (15), P(x, s) is the Laplace transform of p(x, t). The pressure fields on both sides of the medium in the Laplace domain are given by the following: with P 1 (x, s) and P 3 (x, s) are, respectively, the Laplace transform of the field at the left and right of the material; P(s) is the Laplace transform of the incident field; and finally, R(s) and T(s) are the Laplace transforms of the reflection and transmission operators.
The stress-strain equations in the porous medium are given by where σ a ij are the components of the stress tensor that act on the solid (a = s) or fluid (a = f) phase and ij = 1 2 (u ij + u ji ) is the skeletal strain tensor. Since the pressure and stress fields at the boundaries of the medium are continuous [39], the relations between the pressure field and the stresses in the Laplace domain at the interfaces x = 0 and x = L are given by the following: σ s and σ f are the normal stresses in the skeleton and fluid, respectively, and σ s is the shear stress in the skeleton; their expressions in the Laplace domain take the following form: where Φ s Φ f and Ψ s are the Laplace transforms of φ s , φ f , and Ψ s . Hence, using the boundary conditions in (38), the expressions of the scalar functions To determine the reflection and transmission coefficients, we consider the conservation conditions of the flow at the limits of the medium (x = 0 and x = L): In those expressions, V 1 (x, s) and V 3 (x, s) are, respectively, the acoustic velocity fields in the regions x < 0 and x > L, while V s (x, s) and V f (x, s) are the velocities of the solid and fluid in the porous medium, respectively.
The acoustic velocity fields are obtained by considering the field pressure of the surrounding fluid in (34) and (35), and the Euler equation: The expressions of V s and V f are found using the following expressions: Using the boundary conditions in (38) and the expressions for the fluid and solid velocity fields, V s (x, s)), V f (x, s), V 1 (x, s), and V 3 (x, s) , we obtain the expressions for the transmission and reflection coefficients in the Laplace domain: where The coefficients Z i (s) are given by The scattering operators are obtained in the Laplace domain. The purpose of our work is to determine these operators in the time domain, but using the inverse Laplace transform to go directly to the time domain is a very difficult task to perform because of the complexity of these operators. For this reason, modifications and simplifications are applied (Appendix C) to the relations in (40) and (41), which allows us to find the time equivalent of the reflection and transmission operators: where and The expressions of P i (t) and G i (t, L) are given by G i (t), i = 1, 2, 3, are the Green functions of compressional (first and second species) and shear waves, respectively; they describe the propagation of each of these waves inside the porous medium, and their expressions are given in Appendix D.
G i (t), i = 1, 2, 3, are the Green functions of the longitudinal (fast and slow) and transverse waves, respectively. They describe the propagation of each of these waves inside the porous medium, and their expressions are given in Appendix D.
Knowing the strong attenuation that the acoustic wave undergoes in a porous media saturated by a viscous fluid makes multiple reflections negligible; thus, we take only the first reflections at the interfaces x = 0 and x = L in the expressions of the operators in (42) and (45).
The first term of the second member of Equation (43) corresponds to the instantaneous reflection of the wave by the first interface x = 0. The reflected wave does not have time to propagate in the medium, which gives it the advantage of not being dispersive but only attenuated by the term 1+U 1−U . The second term of r(t) expresses the attenuation and dispersion process since it depends on the fractional derivative operator, which takes into account the viscousfluidstructure interaction effects during the wave propagation. This memory term expresses the phenomenon of relaxation of the wave reflected by the first interface.
The term of the second member of Equation (45) corresponds to the reflection by the second interface x = L. This term depends on the Green function of the longitudinal and transverse waves that describe the propagation and dispersion of these waves inside the porous material just as R(t) the temporal expression of the transmission operator in (45) depends on the Green functions of the three Biot waves, highlighting the contribution of each of these waves to the total transmitted wave.
From the temporal expressions of the scattering operators obtained in (42) and (45), the effect of the shear wave plus the effect of the compressional waves (fast and slow) on the porous medium's response to an external load can be examined. It is also possible to study the influence of the characteristic parameters of the porous medium on the amplitudes and speed of these waves.

Numerical Validation
In order to test the validity of our analytical development and to show the influence that the shear waves can have on the response of bone-type porous materials, we proceed to compare transmitted signals obtained during this study, the shear waves being taken into account, and the transmitted signals, where only the longitudinal waves are observed because the shear waves are very weak in normal incidence; however, whether or not they are taken into account is visible on the longitudinal waveforms, as shown in the following study. The characteristics of the bone samples used in this study were taken from work carried out previously by Fellah [67] and are shown in Table 1. The transmitted signals are computed by performing the convolution product of the transmission operator given by (45) with the incident signal, as in the relation in (33).

Porosity Variation
The two waves observed in the numerical simulations ( Figure 2) correspond to the slow and fast longitudinal waves. The shear wave is not noticeable because of its low amplitude at normal incidence; however, its effect is visible on the velocities and amplitudes of the longitudinal waves. This is the subject of the following sensitivity study.  The transmitted signals are simulated by taking or not taking into account the transverse wave for two distinct porosities: φ = 0.2 and φ = 0.85. Figure 2 depicts the impact of the shear wave on different values of the porosity. When the porosity value is equal to 0.2, the shear effect is not felt in the slow longitudinal wave, and its amplitude and speed do not change, unlike the amplitude and speed of the fast wave. For the value of the porosity of 0.85, we notice that, whether the shear wave is considered is much more important than in the previous case (when the porosity was equal to 0.2) for the amplitudes of the fast and slow longitudinal waves. We can therefore conclude at this stage that the effect of shear is more visible as the value of porosity is high and at the level of the amplitudes of both longitudinal waves. The velocities of the longitudinal waves do not seem to be affected by whether or not shear is taken into account.    Readers can notice that the amplitudes of the two compressional waves (the fast and slow) are different for the two values of tortuosity, even though it is more apparent for α ∞ = 1.3 than it is for α ∞ = 1. The velocities of both compressional waves remain unchanged, regardless of whether the effect of the transverse waves is included or not.

Tortuosity Variation
Tortuosity is an important parameter to describe the inertial effects, which occur between the fluid and solid parts of a porous material because it provides information on the geometry of the pores of the porous medium, which are rarely straight and have variable sections. Like porosity variations, tortuosity variations cause changes in the fast and slow waves and their velocities.
Comparing the figures for two values of the tortuosity, we notice that an increase in the tortuosity induces an increase in the amplitude of the fast wave and a decrease in the amplitude of the slow-wave. A temporal shift of the two Biot waves occurs (decrease in the velocities), which is more distinct on the slow wave than it is on the fast wave as the tortuosity increases.

Variation of the Viscous Characteristic Length Λ
The effect of whether transverse waves are considered along with longitudinal waves on the transmitted signals for two different viscous characteristic lengths is illustrated in Figure 4. In the comparison of the simulated transmitted signals, the effect of shear can be observed on the amplitudes of both the fast and slow waves. This effect leads to an increase in the amplitudes of the fast and slow waves regardless of the value of the viscous characteristic length. The speeds of the fast and slow waves remain unchanged in all cases. The comparison between the figures corresponding to two values of the viscous characteristic length provides information on the effect of the viscous characteristic length on the transmitted signals. This parameter was first introduced by Johnson et al. [33] in 1987 and is considered to be a geometrical parameter indicative of small pore sizes; it represents the thickness where viscous effects are important.  Through our comparison, we can see that, by increasing the viscous characteristic length from 5 µm to 9 µm, the amplitude of the slow wave increases while the velocities of the slow and fast waves as well as the amplitude of the fast wave remain practically the same.
The viscous characteristic length plays a less important role than the tortuosity, since the variation in its value barely affects the transmitted signals whether the shear wave is considered or not.

Solid Density ρ s Variation
The difference between the transmitted signals when shear effects are considered and not considered along with longitudinal waves for two distinct values of the solid density ρ s = 1000 kg/m 3 and ρ s = 2400 kg/m 3 is given in Figure 5. Unlike for the previous case of the effect of the viscous characteristic length, we notice that the only significant difference is for the amplitudes of the fast waves that increase when the shear waves are considered for both values of the solid density; the slow waves are practically the same for both signals in Figure 5. It is worth mentioning that, for a high value of the density ρ s = 2400 kg/m 3 , the fast wave with shear taken into account arrives faster than the one with shear not taken into account, which means that its speed is greater without than with shear taken into account, which is not observed for the slow wave.  We can also notice from these figures that the fast waves (amplitude and speed) are sensitive to the variation of the solid density; indeed, a decrease in the amplitudes of the fast waves as well as a decrease in their speed appear when the solid density goes from 1000 (kg/m 3 ) to 2400 (kg/m 3 ).

Bulk Modulus of the Elastic Solid K s Variation
The simulated transmitted signals with and without considering the effect of transverse waves in addition to longitudinal waves are compared in Figure 6 for two values of the bulk modulus of the elastic solid K s . We notice that the effect of the shear wave is visible mainly on the amplitude of the fast wave for the value of K s = 2.45 GPa, and the amplitude of the slow wave varies but in a lesser way. For the value of K s = 10 GPa, it is the amplitude of the slow wave that is most sensitive to the incorporation of shear in the model, with the change in the amplitude of the fast wave being, however, smaller for this value of K s . Shear affects the amplitudes of the fast and slow waves differently depending on the values of K s .  The speed of the fast wave also seems to be sensitive to the shear; in fact, when decreasing the value of K s , the fast wave arrives delayed when the shear is not taken into account. In other words, taking into account the shear increases the speed of the fast wave, which makes it arrive faster in comparison with the fast wave obtained without shear. This phenomenon is not observable for the speed of the fast wave when K s = 10 GPa or for the speeds of the slow waves, regardless of the value of K s .

Bulk Modulus of the Porous Skeletal Frame K b Variation
The sensitivity study of the bulk modulus of the porous skeletal frame K b taking or not taking into account the shear shows quite interesting results compared with the previous cases (see Figure 7). For example, we notice as an example that, for the high value of K b = 3.6 GPa, the amplitudes of the slow and fast waves are very affected by whether shear is considered. The amplitude of the slow wave in this case can be doubled when the shear is taken into account, while the amplitude of the slow wave decreases sharply. For the value of K b = 1.1 GPa, the variation in the slow wave amplitude is very small while that of the fast wave is much larger. We also notice generally that the velocities of the waves do not practically change with the taking into account of the shear, except for that of the fast wave for K b = 1.1 GPa, where we see that taking shear into account of the shear decreases the velocity of the fast wave, which slows down compared with the case where the shear is taken into account.

The Shear Modulus N Variation
The simulated transmitted signals including and excluding the effect of the shear wave in addition to the longitudinal waves for two values of the shear modulus are shown in Figure 8. The comparison between the two signals for each value of N allows us to say that the inclusion of the shear wave in the model affects the fast wave more than the slow wave and more particularly for large values of N (N = 8.9 GPa). However, the effect of shear is more visible on the speed of the fast wave for low values of N (N = 5.8 GPa). However, the slow wave velocities remain unchanged whether shear is considered or not.

Discussion and Conclusions
This study allows us to examine the influence of the shear wave in the modeling of ultrasound propagation in porous media. This effect is studied by comparing the transmitted waveforms by varying each physical parameter of the porous material (porosity, tortuosity, viscous characteristic length, density of the solid, bulk modulus of the elastic solid, bulk modulus of the porous skeletal frame, and shear modulus). The amplitudes and velocities of the transmitted longitudinal waves (slow and fast) are observed and compared in order to see the impact of taking into account the shear wave in the model. In our previous work [19,62], we did not take into account the shear wave at normal incidence, but recent studies [62,67] have shown that the shear wave must be taken into account because of the viscous exchanges between fluid and structure. At normal incidence, the shear wave is not as visible as the fast and slow longitudinal waves because its amplitude is very low; however, its effect is felt in the amplitudes and speeds of longitudinal waves, and this effect varies depending on the values of the physical parameters of the porous medium, as shown in this study.
The originality of our approach is that it makes it possible not only to model in the time domain acoustic propagation in a porous medium using fractional calculus but also to obtain analytical solutions to the response (reflection and transmission operators) of the porous medium to an acoustic excitation (incident wave). This analytical solution in the time domain expressed by the Green functions of longitudinal (fast and slow) and transverse waves using fractional calculation is the strong point (original) of our work. Indeed, the the existing studies do not provide analytical solutions but only propose numerical solutions.
The study of the sensitivity of the physical parameters of the longitudinal waves by taking or not taking into account the shear wave in the model allows us to see and quantify the importance of the shear wave in the modeling and its impact on the longitudinal waves. It appears that the variation in the the values of some parameters such as the bulk modulus of the porous skeletal frame reflects a significant change in the amplitude of the fast and slow longitudinal waves. The velocity of the fast wave is also modified, contrary to that of the fast wave, which is unchanged, and thus independent of whether shear is taken into account. The variation in other parameters (such as porosity, tortuosity, viscous characteristic length, solid density, and shear modulus) modifies the amplitude of the fast wave more than that of the slow wave. Practically for all the parameters, their variation has no impact on the speed of the slow wave, only on that of the fast wave.
The results of this study clearly show that the inclusion of the shear wave in the model has a considerable impact on the amplitude of the slow and fast waves as well as on the speed of the fast wave. The fractional calculus have proven to be very useful in obtaining the analytical solutions in the time domain. We can then conclude that the shear wave at normal incidence should no longer be ignored as it was the case in the previous studies and that it is necessary to consider the shear effects, as this will have consequences for the acoustic and mechanical characterization of porous materials by solving the inverse problems using experimental data of transient waves.

Conflicts of Interest:
The authors declare no conflicts of interest.

Appendix A. The Eigenvalues and Eigenvectors of the Matrix M
The matrix M is given in the Laplace domain by the following: The eigenvalues of the matrix M are written as follows: where Tra(M) and Det(M) are, respectively, the trace and determinant of this same matrix: Posing we obtain the expressions of the eigenvalues: The expressions of the components V 1,2 of the eigenvectors are given by . with: The stress-strain relation for the solid in the Laplace domain is in the following form: where the solid displacement vectors are given by The normal component of the stress tensor acting on the solid part becomes the following: Considering the expressions for Φ 1 (x, s) and Φ 2 (x, s), we obtain the following: with Z 1 (s) = (V 1 (s)Q + P)λ 1 (s), Z 2 (s) = (V 2 (s)Q + P)λ 2 (s).
For shear wave, the stress-strain relationship is written as follows: Using the fact that the displacement vector of the solid is the shear stress becomes The continuity conditions for solid and fluid stresses at x = 0 and x = L are as follows: By using the expressions of the constraints and those of the pressure fields in (34) and (35), we obtain sinh(l χ(s)) , ).

Appendix C. The Scattering Operators in the Time Domain
The expressions for the scattering operators in the Laplace domain are given by the following: These expressions can be decomposed into simple elements: Using the expressions of F 4 (s) and F 5 (s) given by (42), the calculus of (F 4 (s) + F 5 (s)) as well as that of (F 4 (s) − F 5 (s)) gives the following: Posing , we obtain where e −l √ λ i (s) are assumed to be very small for large values of s because of the high-frequency condition assumed in our study, which allows us to use the following development: By considering only the small powers of W i (s),we obtain the following approximations: For large values of s, we obtain the expressions of the reflection and transmission operators in the Laplace domain: (A14) To obtain the Laplace inverse transform for these operators, it would be easier to explicitly define X(s), Y(s) and Z(s). Knowing that (A15) Using the high-frequency approximation, s λ 1 (s) we obtain X(s) The same steps for Y(s) and Z(s) lead us to the following relations, with the following notations: and by posing where Then, by putting U = x 1 + x 2 + x 3 , V = y 1 + y 2 + y 3 , .

Appendix D. Green Functions of Longitudinal and Transverse Waves
The following relations are used for the Green functions G i (t, K) of the porous material for the fast wave (i = 1), the slow (i = 2), and the rotational (i = 3) waves.
where h i (τ, ζ) has the following form, and where