Introduction

The study of seismic waves propagation caused by an earthquake phenomena are crucial to examine the interior of Earth as well as to analyze the seismic activity at the different margins of the Earth. The Earth consists of various kind of layers which are subjected to seismic waves due to the occurrence of earthquake. The propagation characteristic of seismic waves in layered structures is of primary importance due to its applications in geophysical prospective, civil engineering, mechanical engineering, etc. Ewing et al1. elaborated the different type of Earth structure according to the dispersive property of elastic layered media. Many investigations2,3have been done over the propagation behaviour of elastic waves in layered anisotropic media. A number of literature4,5provide the explanation on characterization of elastic waves in elastic layered media. Love waves are one of the devastating seismic surface waves produced by the earthquake. It may have an enormous effect on engineering designs, reservoirs and other man made structures6,7. In order to minimize the catastrophic effects of Love waves on human life and engineering structures, analysis of its propagation characteristics on different types of earth structures is crucial. Chattopadhyay et al8. analysed the propagation behaviour of Love waves in an isotropic elastic media lying over a semi-infinite elastic media and concluded a uniform decrease in the Love waves phase velocity with increasing non-dimensional wave number. Palermo et al9. investigated the control of Love wave dispersion through resonator parameter variation. Saha et al10. incorporated a numerical modelling of SH-wave propagation in initially-stressed multilayered composite structures.

In view of the problems with elastic stability for anisotropic media, the investigation of wave propagation in fluid saturated porous layer has been found to be very crucial. A porous medium is a type of material medium made up of the solid skeleton including with pore space that is often filled with fluid, is regarded as the poroelastic medium11. As the layer of the Earth is made up of various kind of rocks such as sandstone, lime stones, shale, etc., exhibiting the poroelastic nature. These poroelastic rocks are subjected to seismic waves during the occurrence of earthquake phenomena. Due to this, the investigation of the seismic waves propagation in poroelastic media is a matter of great interest for geoscientists and engineers. In view of above-mentioned reasons, many researchers made a lot of efforts to study the influence of properties of the poroelastic rock structures on the seismic wave’s propagation. The basic concepts regarding the waves propagation in poroelastic rock medium have been given by Biot12. Konczak13analysed the nature of Love wave propagation in the anisotropic porous layer. Also investigations have been done over the propagation of SH- wave in a poroelastic layer and on the dynamic interfacial crack effect in poroelastic strip14,15.

Generally the Earth medium are not homogeneous. Some kind of heterogeneity are found due to the variation of elastic properties in vertical depth. The heterogeneity can occur due to the phenomena of formation, growth, and coalescence of microcracks inside the rock medium16. The literature gives substantial evidence that the rock structures in the earth medium contains several sorts of heterogeneity, which can be expressed mathematically by different functions17,18. Based upon this idea, a lot of investigations have performed on the impact of various forms of heterogeneity viz. linear, quadratic, exponential or any other non-linear function of depths on the SH-wave’s propagation in an isotropic layered structure19,20. It has been examined by many studies that the Love wave’s propagation is strongly influenced by the heterogeneity of the structures. Kumari et al21. and Kundu et al22. investigated the effect of heterogeneity parameter on the Love wave’s propagation in a layered isotropic medium and heterogeneous micropolar media respectively. Gupta et al23. analysed the influence of propagation behaviour of Love wave on heterogeneity anisotropic porous media. Saha et al24. analyzed the impact of inhomogeneity on SH-type wave propagation in an initially stressed composite structure. Also, Sattari et al25. analyzed the body waves propagation through heterogeneous and discontinuous media. The problem having non-homogeneous boundary conditions may be homogenized with the use of some variable changes and will be able to provide solutions. Such type of non-homogeneous boundary value problem can be easily solved by using Green’s function methodology. Covert26was the researcher, who evolve a fresh method for finding green function for developed bodies. Also, Saha et al27. and Mario et al28. analyzed the behaviour of SH-waves propagation subjected to the curved boundary and imperfect boundary conditions respectively. Further, Venkatesan29 studied the Love-type wave propagation due to an impulsive point source by using the Green’s function technique.

Earthquakes are triggered by body forces, modeled as space-time dependent impulsive line sources, typically represented by the Dirac delta function. The resulting elastic displacement in layered structures is determined using the Green’s function technique, a proven method for solving elastodynamic problems involving impulsive sources. Several investigations30,31have been carried out over the propagation behaviour of Love waves utilizing Green’s function approach to handle impulsive point source problem in a heterogeneous medium, viscoelastic FGM half-space media respectively. Gupta et al32. analysed the effect of point source on propagation of SH-wave through an orthotropic crustal layer. Further, Gupta et al33. incorporated a comparative analysis due to the effect of point source on generation of SH wave. Hoop and Hijden34investigated the acoustic wave motion which depends on time-space are proliferated due to the impulsive line source. Also, Chattopadhyay et al35. investigated the propagation behaviour of SH-wave generated by impulsive line source in elastic media by considering Green’s function approach. Recently, Mahanty et al36. investigated the shear wave propagation in a heterogeneous poroelastic layered structure subjected to an impulsive line source. In many porous rocks, interconnected pores allow fluid to flow freely. However, fluid can sometimes become trapped, and prolonged entrapment can lead to mineral dissolution, creating patches within the rock. This process can form unique rock types where spherical inclusions are embedded within a porous host rock. As a result, two types of porosity emerge, storage porosity and transport porosity37. These rock formations, known as double porous structures, are commonly found near reservoir regions38. The formation of a large reservoir can cause the occurrence of an earthquake in the surrounding area39,40. Therefore, it is vital to investigate the behaviour of Love wave’s propagation in these double porous rock structures. Many investigations41,42have been carried out over the propagation behaviour of Love waves in various kind of double porosity structure. Recently, Gupta et al43. analyzed the propagation characteristic of Love wave in heterogeneous porosity media having discontinuity.

A critical review of the available literature reveals a gap in knowledge regarding the impact of impulsive line sources on Love wave propagation in heterogeneous double-porosity composite rock structures. Notably, no prior studies have investigated the propagation behavior of Love waves in a layered composite double-porosity structure consisting of a transversely isotropic double-porosity (TIDP) layer overlying a heterogeneous isotropic double-porosity (IDP) half-space influenced by an impulsive line source. It motivates the authors for performing the present investigation, which aims to formulate a problem that more closely reflects realistic scenarios by incorporating an impulsive line source in a double-porosity composite structure. This study assumes a non-linear function of vertically depths and stiffness to represent the material’s heterogeneity. To achieve an efficient analytical solution, the dispersion equation is derived in closed form by employing Green’s function technique in conjunction with Fourier transforms. The investigation is motivated with following research questions which aroused after carrying out a detail literature survey.

Key questions

  • What are the influence of heterogeneity parameter (linear as well as non linear quadratic heterogeneity) of the IDP half-space on the phase velocity of Love waves in the double porous rock structure?

  • How the porosity parameter (corresponds to both upper TIDP layer and lower IDP half-space) can effect the phase velocity of the Love waves in the double porous rock structure?

  • What are the effect of anisotropic parameter of TIDP layer on the phase velocity of Love waves in the double porous rock structure?

Mathematical model of dynamical problem

In the present study, we assume the Love wave’s propagation caused by an impulsive line source (S) in a double porous rock structure comprises an upper transversely isotropic double porous layer overlying the lower heterogeneous isotropic double porous half-space medium. Let H is the thickness of layer measured from the origin O. Here, the half-space medium is not homogeneous because of the consideration of varying stiffness and density. The heterogeneity in half-space can be expressed as the non-linear function of vertically downward depth. The Cartesian co-ordinate system containing (xz) axes has been selected from this perspective that the waves propagation is considered in x-direction while z-axis is considered to be placed in the positive direction of vertical depth. The location of a line source S (impulsive in nature) is considered at a distance H below from the origin. It also represents the interface (\(z=H\)) between the upper layer and lower half-space medium as depicted in the Figure 1. The particle displacements in x, y, z direction are \(u_j\), \(v_j\) and \(w_j\) respectively, as for the layer j=1 and for the half-space medium j=2 respectively.

Fig. 1
figure 1

Schematic diagram.

In view of the properties of Love waves in xz-plane causing the particle displacement in y-axis only, the displacement components for Love wave’s propagation can be represented as

$$\begin{aligned} \left. \begin{array}{l} {u_j}= 0,\,\,\ w_j=0,\,\ v_j= v_j{(x,z,t)}\, \frac{\partial {(.)}}{\partial y} \equiv 0\\ U_k^{(l)}= 0, \,W_k^{(l)}= 0,\, V_k^{(l)}= {V_k^{(l)}{(x,z,t)}} \end{array} \right\} \end{aligned}$$
(1)

where \(u_k^{(l)}, v_k^{(l)}, w_k^{(l)}\) and \(U_k^{(l)}, V_k^{(l)}, W_k^{(l)}\) are defined in Nomenclature respectively.

Dynamics of the upper transversely isotropic homogeneous double porous (TIDP) layer influenced by the line source S:

The constitutive equations governing the TIDP layer can be expressed as44

$$\begin{aligned} \left. \begin{array}{l} \tau _{xx}^{(1)} = (\lambda _1 + 2N_1) \frac{\partial u_1}{\partial x} + \lambda _1 \frac{\partial v_1}{\partial y} +\mu _1 \frac{\partial w_1}{\partial z} + S_1(\xi _1^1 + \phi _1^2 \zeta _1) + S_2(\xi _1^2-\phi _1^1\zeta _1) \\ \tau _{yy}^{(1)}=\lambda _1 \frac{\partial u_1}{\partial x} + (\lambda _1 + 2N_1)\frac{\partial v_1}{\partial y} +\mu _1 \frac{\partial w_1}{\partial z} + S_1(\xi _1^1 + \phi _1^2 \zeta _1) + S_2(\xi _1^2-\phi _1^1\zeta _1) \\ \tau _{zz}^{(1)}=\mu _1 \frac{\partial u_1}{\partial x} + \mu _1 \frac{\partial v_1}{\partial y} +(I_1 + 2L_1) \frac{\partial w_1}{\partial z} + S_3(\xi _1^1 + \phi _1^2 \zeta _1) + S_4(\xi _1^2-\phi _1^1\zeta _1)\\ \tau _{zx}^{(1)}=L_1\left( \frac{\partial u_1}{\partial z}+\frac{\partial w_1}{\partial x} \right) \\ \tau _{zy}^{(1)} =L_1\left( \frac{\partial v_1}{\partial z}+\frac{\partial w_1}{\partial y} \right) \\ \tau _{xy}^{(1)} =N_1\left( \frac{\partial u_1}{\partial y}+\frac{\partial v_1}{\partial x} \right) \\ -p_1^{f_1}= S_1\left( \frac{\partial u_1}{\partial x}+\frac{\partial v_1}{\partial y} \right) +S_3\frac{\partial w_1}{\partial z} +K_1(\xi _1^1+\phi _1^2\zeta _1)\\ -p_1^{f_2}=S_2\left( \frac{\partial u_1}{\partial x}+\frac{\partial v_1}{\partial y} \right) +S_4\frac{\partial w_1}{\partial z} +K_2(\xi _1^2+\phi _1^1\zeta _1) \end{array} \right\} \end{aligned}$$
(2)

Considering equations (1) and (2), the non-vanishing equations of motion describing the Love wave’s propagation within a fluid-saturated TIDP layer can be stated as41

$$\begin{aligned} \left. \begin{array}{l} \frac{\partial \tau _{yx}^{(1)}}{\partial x}+\frac{\partial \tau _{yz}^{(1)}}{\partial z}=\rho _{00}^{(1)}\frac{\partial ^2v_1}{\partial t^2}+\rho _{01}^{(1)}\frac{\partial ^2V_1^{(1)}}{\partial t^2} + \rho _{02}^{(1)}\frac{\partial ^2 V_1^{(2)}}{\partial t^2}\\ \frac{\partial ^2V_1^{(1)}}{\partial t^2}= -\frac{\rho _{01}^{(1)}}{\rho _{11}^{(1)}}\frac{\partial ^2v_1}{\partial t^2}\\ \frac{\partial ^2V_1^{(2)}}{\partial t^2}= -\frac{\rho _{02}^{(1)}}{\rho _{22}^{(1)}}\frac{\partial ^2v_1}{\partial t^2} \end{array} \right\} \end{aligned}$$
(3)

Now, introducing the function associated with the distribution of force density \(F_1(r,t)\) caused by an impulsive line source at the interface between TIDP layer and IDP half-space acting inside the medium. The existing equation of motion describing Love wave propagation in the TIDP layer caused by an impulsive line source derived from equations (13) may be implemented as

$$\begin{aligned} \frac{\partial }{\partial x}\left( N_1 \frac{\partial v_1}{\partial x}\right) + \frac{\partial }{\partial z}\left( \text{L } _1 \frac{\partial v_{1}}{\partial z}\right) + 4 \pi F_1{(r,t)}&= \rho ^{(1)}\frac{\partial ^2 v_1}{\partial t^2}, \end{aligned}$$
(4)

where, r represents the distance between origin and coordinate points where the force is applied and

$$\begin{aligned} \rho ^{(1)}&=\rho _{00}^{(1)}-\frac{(\rho _{01}^{2})^{(1)}}{\rho _{11}^{(1)}}-\frac{(\rho _{02}^{2})^{(1)}}{\rho _{22}^{(1)}}. \end{aligned}$$
(5)

In view of harmonic nature of Love wave’s propagation, the following transformation can be taken

$$\begin{aligned} F_1{(r,t)} = e^{i\omega t} {\bar{F}}_1{(r)} \,\,\, \text {and} \,\,\, v_1 {(x,z,t)} = e^{i\omega t} \bar{v}_1{(x,z)}. \end{aligned}$$
(6)

Using equation (6), equation (4) yields

$$\begin{aligned} \frac{\partial ^2\bar{v}_1}{\partial x^2} +\left( \frac{L_1}{N_1}\right) \frac{\partial ^2\bar{v}_1}{\partial z^2}+\frac{\rho ^{(1)}}{N_1}\omega ^2\bar{v}_1&=\frac{4\pi {\bar{F}}_1(r)}{N_1}. \end{aligned}$$
(7)

We can now express the force density distribution function, \(F_1(r)\), caused by the impulsive line source disturbance by using the Dirac delta function.

$$\begin{aligned} F_1{(r)}= \delta {(x)}\delta {(z-H)} \end{aligned}$$
(8)

Using equation (8) in equation (7), the only existing equation of motion impacted by the line source for the TIDP layer can be represented as

$$\begin{aligned} \frac{\partial ^2\bar{v}_1}{\partial x^2} +\left( \frac{L_1}{N_1}\right) \frac{\partial ^2\bar{v}_1}{\partial z^2}+\frac{\rho ^{(1)}}{N_1}\omega ^2\bar{v}_1=\frac{4\pi \delta {(x)}\delta {(z-H)}}{N_1}. \end{aligned}$$
(9)

Considering the following Fourier transform of \(\bar{v}_1{(x,z)}\) as

$$\begin{aligned} V_1{(\Theta ,z)}=\frac{1}{2\pi }\int _{-\infty }^{\infty } \,\bar{v}_1{(x,z)}e^{i\Theta x}dx, \end{aligned}$$
(10)

and the related inverse Fourier transform as

$$\begin{aligned} \bar{v}_1{(x,z)}=\frac{1}{2\pi }\int _{-\infty }^{\infty } \, V_1{(\Theta ,z)}e^{-i\Theta x}d\Theta , \end{aligned}$$
(11)

where \(\Theta\) denotes the transformation parameter.

Using equation (10) in equation (9), we get

$$\begin{aligned} \frac{1}{\gamma ^2}\frac{d^2V_1}{dz^2}-\alpha _1^2V_1=4\pi F_1{(z)}, \end{aligned}$$
(12)

where \(4\pi F_1{(z)}\) = \(2\delta {(z-H)}/N_1\) and \(\alpha _1^2\)=\(\Theta ^2-{(\omega /\beta _1)^2}\) such that

$$\begin{aligned} \beta _1=\sqrt{\frac{N_1}{\rho ^{(1)}}}, \end{aligned}$$
(13)

denotes the velocity of Love waves corresponding to the TIDP layer and \(\gamma ^2\)=\({\frac{L_1}{N_1}}\) denotes the anisotropic parameter of TIDP layer.

Dynamics of isotropic heterogeneous double porous (IDP) half-space

The constitutive equations governing the IDP half-space is being expressed by45,

$$\begin{aligned} \left. \begin{array}{l} \tau _{xx}^{(2)} = (\lambda _2 + 2N_2) \frac{\partial u_2}{\partial x} + \lambda _2 \frac{\partial v_2}{\partial y} +\lambda _2 \frac{\partial w_2}{\partial z} + S_5(\xi _2^1 + \phi _2^2 \zeta _2) + S_6(\xi _2^2-\phi _2^1\zeta _2) \\ \tau _{yy}^{(2)}=\lambda _2 \frac{\partial u_2}{\partial x} + (\lambda _2 + 2N_2)\frac{\partial v_2}{\partial y} +\lambda _2 \frac{\partial w_2}{\partial z} + S_5(\xi _2^1 + \phi _2^2 \zeta _2) + S_6(\xi _2^2-\phi _2^1\zeta _2) \\ \tau _{zz}^{(2)}=\lambda _2 \frac{\partial u_2}{\partial x} + \lambda _2 \frac{\partial v_2}{\partial y} + (\lambda _2 + 2N_2) \frac{\partial w_2}{\partial z} + S_5(\xi _2^1 + \phi _2^2 \zeta _2) + S_6(\xi _2^2-\phi _2^1\zeta _2)\\ \tau _{zx}^{(2)}=N_2\left( \frac{\partial u_2}{\partial z}+\frac{\partial w_2}{\partial x} \right) \\ \tau _{zy}^{(2)} =N_2\left( \frac{\partial v_2}{\partial z}+\frac{\partial w_2}{\partial y} \right) \\ \tau _{xy}^{(2)} =N_2\left( \frac{\partial u_2}{\partial y}+\frac{\partial v_2}{\partial x} \right) \\ -p_2^{f_1}= S_5\left( \frac{\partial u_2}{\partial x}+\frac{\partial v_2}{\partial y} +\frac{\partial w_2}{\partial z}\right) +K_3(\xi _2^1+\phi _2^2\zeta _2)\\ -p_2^{f_2}=S_2\left( \frac{\partial u_2}{\partial x}+\frac{\partial v_2}{\partial y} \right) +\frac{\partial w_2}{\partial z} +K_4(\xi _2^2+\phi _2^1\zeta _2) \end{array} \right\} \end{aligned}$$
(14)

Considering equations (1) and (14), the non-vanishing equations of motion describing the Love wave’s propagation for IDP half-space can be represented as follows41

$$\begin{aligned} \left. \begin{array}{l} \frac{\partial \tau _{yx}^{(2)}}{\partial x}+\frac{\partial \tau _{yz}^{(2)}}{\partial z}=\rho _{00}^{(2)}\frac{\partial ^2v_2}{\partial t^2}+\rho _{01}^{(2)}\frac{\partial ^2V_2^{(1)}}{\partial t^2} + \rho _{02}^{(2)}\frac{\partial ^2 V_2^{(2)}}{\partial t^2}\\ \frac{\partial ^2V_2^{(1)}}{\partial t^2}= -\frac{\rho _{01}^{(2)}}{\rho _{11}^{(2)}}\frac{\partial ^2v_2}{\partial t^2}\\ \frac{\partial ^2V_2^{(2)}}{\partial t^2}= -\frac{\rho _{02}^{(2)}}{\rho _{22}^{(2)}}\frac{\partial ^2v_2}{\partial t^2} \end{array} \right\} \end{aligned}$$
(15)

The presence of material heterogeneity can significantly impact the behaviour of the medium’s elastic constants and density. On account of, considering the heterogeneity in the form of linear as well as the quadratic function of zin vertical depths and stiffness as19

$$\begin{aligned} {N_2= N_{02} + R}, \end{aligned}$$
(16)

where \(R ={\epsilon _1(z-H) + \epsilon _2(z-H)^2}\).

Here, \(\epsilon _1\) and \(\epsilon _2\) denote the linear and quadratic heterogeneity parameters respectively. \(\epsilon _1\) and \(\epsilon _1\) are very small such that \(\epsilon _1<<1\) and \(\epsilon _1<<1\). \(N_{02}\) denotes the stiffness for the isotropic double porous half-space.

On account of equations (1), (14), (15) and (16), the only existing equation of motion for the Love waves propagation in the IDP half-space may be stated as

$$\begin{aligned} \frac{\partial ^2{v}_2}{\partial x^2} + \frac{\partial ^2{v}_2}{\partial z^2}+\frac{\rho ^{(2)}}{N_{02}}\omega ^2{v}_2&=\frac{-1}{N_{02}}\left[ R \frac{\partial ^2{v}_2}{\partial x^2}+\{\epsilon _1+2\epsilon _2{(z-h)\}\frac{\partial {v_2}}{\partial z}} +R \frac{\partial ^2{v}_2}{\partial z^2}\right] . \end{aligned}$$
(17)

Now, using the following transformation (due to the harmonic nature of waves), we have

$$\begin{aligned} v_2 {(x,z,t)} = e^{i\omega t} \bar{v}_2{(x,z)}. \end{aligned}$$
(18)

Using equation (18), equation (17) yields

$$\begin{aligned} \frac{\partial ^2\bar{v}_2}{\partial x^2} + \frac{\partial ^2\bar{v}_2}{\partial z^2}+\frac{\rho ^{(2)}}{N_{02}}\omega ^2\bar{v}_2&=\frac{-1}{N_{02}}\left[ R \frac{\partial ^2\bar{v}_2}{\partial x^2}+\{\epsilon _1+2\epsilon _2{(z-h)\}\frac{\partial \bar{v_2}}{\partial z}} +R \frac{\partial ^2\bar{v}_2}{\partial z^2}\right] . \end{aligned}$$
(19)

On using the following Fourier transform for \(\bar{v}_2{(x,z)}\)

$$\begin{aligned} V_2{(\Theta ,z)}=\frac{1}{2\pi }\int _{-\infty }^{\infty } \,\bar{v}_2{(x,z)}e^{i\Theta x}dx, \end{aligned}$$
(20)

and the related inverse Fourier transform can be given as

$$\begin{aligned} \bar{v}_2{(x,z)}=\frac{1}{2\pi }\int _{-\infty }^{\infty } \, V_2{(\Theta ,z)}e^{-i\Theta x}d\Theta . \end{aligned}$$
(21)

Now, using the equation (20) in equation (19), we get

$$\begin{aligned} \frac{d^2V_2}{dz^2}-\alpha _2^2V_2=4\pi F_2{(z)}, \end{aligned}$$
(22)

where \(\alpha _2^2\) = \(\Theta ^2-{(\omega /\beta _2)^2}\) in which

$$\begin{aligned}\beta _2=\sqrt{\frac{N_{02}}{\rho ^{(2)}}}\end{aligned}$$

is represented as the velocity of the love waves corresponding the half-space medium, where

$$\begin{aligned} \rho ^{(2)}=\rho _{00}^{(2)}-\frac{(\rho _{01}^{2})^{(2)}}{\rho _{11}^{(2)}}-\frac{(\rho _{02}^{2})^{(2)}}{\rho _{22}^{(2)}}, \end{aligned}$$
(23)

and

$$\begin{aligned} 4\pi F_2(z)&=\frac{-1}{N_{02}}\left[ R\frac{d^2{V}_2}{d x^2}+\{\epsilon _1+2\epsilon _2{(z-h)\}\frac{d{V_2}}{dz}}-R\Theta ^2{V_2}\right] . \end{aligned}$$
(24)

Boundary conditions

On considering the schematic diagram 1 of present problem and equations (25101416), the following possible boundary conditions can be provided

  • The upper transversely isotropic double porous layer has traction free upper surface, i.e.

    $$\begin{aligned} \frac{dV_1}{dz}=0\,\,\,\ \text {at}\,\,\,\ z=0 \end{aligned}$$
    (25)
  • At the common interface between the double porous layer and the half-space, there are continuous components for both stress and displacement. i.e, at \(z=H\),

    $$\begin{aligned} L_1\frac{dV_1}{dz}&= N_{02}\frac{dV_2}{dz}\,\,\,\ \end{aligned}$$
    (26)
    $$\begin{aligned} V_1&=V_2\,\,\,\,\,\,\, \end{aligned}$$
    (27)

Solution of problem

We will now use the Green’s function technique for solving equations (12) and (22) with the aid of the necessary boundary conditions (25-27).

Let us introduce the Green’s function \(G_1{(z,z_0)}\) for any arbitrary \(z_0\) for the upper TIDP layer holding the boundary condition \(\frac{dG_1}{dz}\)=0 at \(z=0\) and at \(z=H\). Therefore, \(G_1{(z,z_0)}\) satisfies the equation

$$\begin{aligned} \frac{1}{\gamma ^2}\frac{d^2G_1{(z,z_0)}}{dz^2}-\alpha _1^2G_1{(z,z_0)}=\delta {(z-z_0)}, \end{aligned}$$
(28)

where \(z_0\) locates any where in the TIDP layer.

By multiplying equation (12) and equation (28) by \(G_1{(z,z_0)}\) and \(V_1(z_0)\) respectively, subtracting the resultant equations, and integrating within the range z=0 to \(z=H\), we get

$$\begin{aligned} V_1(z_0)&= \frac{2}{N_1}G_1(H,z_0)-\frac{1}{\gamma ^2}G_1(H,z_0)\left[ \frac{dV_1}{dz}\right] _{z=H} \end{aligned}$$
(29)

Replacing \(z_0\) with z and utilizing Green’s function symmetry \(G_1(H,z) = G_1(z,H)\), equation (29) yields

$$\begin{aligned} V_1(z)&= \frac{2}{N_1}G_1(z,H)-\frac{1}{\gamma ^2}G_1(z,H)\left[ \frac{dV_1}{dz}\right] _{z=H} \end{aligned}$$
(30)

Similarly, for IDP half-space, let \(G_2{(z,z_0)}\) be the Green’s function with fulfilling the following equations:

$$\begin{aligned} \frac{d^2G_2{(z,z_0)}}{dz^2}-\alpha _2^2G_2{(z,z_0)}=\delta {(z-z_0)}, \end{aligned}$$
(31)

where \(z_0\) locates any where in the IDP half-space medium.

The similar procedure as discussed for the case of upper TIDP layer can be used for the lower IDP half-space. The equations (22) and (31) are multiplied by \(G_2{(z,z_0)}\) and \(V_2\) respectively; subtracting the resultant equations, and then integrating within the range z=H to \(z=\infty\), and using Green’s function property, we get

$$\begin{aligned} V_2(z)&= G_2(z,H)\left[ \frac{dV_2(z)}{dz}\right] _{z=H}+\int _{H}^{\infty } T G_2{(z,z_0)} dz_0, \end{aligned}$$
(32)

where \(T = 4\pi F_2(z_0)\).

On using equations (30), (32) and boundary conditions (26-27), we get

$$\begin{aligned} \left[ \frac{dV_1(z)}{dz}\right] _{z=H}&=\frac{1}{A}\left[ \frac{2}{N_1}G_1(H,H) - \int _{H}^{\infty } T G_2{(H,z_0)} dz_0\right] , \end{aligned}$$
(33)

where \(A = \frac{1}{\gamma ^2}G_1{(H,H)}+\frac{L_1}{N_{02}} G_2{(H,H)}\).

Using the equations (33) in equation (30) and putting the value of \(4\pi F_2{(z_0)}\) from equation (24), we get

$$\begin{aligned} V_1(z)&=\frac{2}{N_1}G_1(z,H) \frac{L_1}{AN_{02}}G_2(H,H)-\frac{1}{A\gamma ^2}G_1(z,H)\nonumber \\ &\phantom { =}\left[ \int _{H}^{\infty }\frac{1}{N_{o2}}\left( R'\frac{d^2V_2}{dz_0^2}+\{\epsilon _1 +2\epsilon _1(z_0-H)\}\frac{d^2V_2}{dz_0}-R'\Theta ^2V_2\right) G_2(H,z_0) dz_0 \right] , \end{aligned}$$
(34)

where

$$\begin{aligned}R' ={\epsilon _1(z_0-H) + \epsilon _2(z_0-H)^2}.\end{aligned}$$

Furthermore, using the equations (33) and boundary condition(26), equation (32) leads to the following form

$$\begin{aligned} V_2{(z)}&=\frac{G_2(z,H)}{A}\frac{L_1}{N_{02}}\left[ \frac{2}{N_1}G_1(H,H) - \int _{H}^{\infty } T G_2{(H,z_0)} dz_0\right] +\int _{H}^{\infty } T G_2{(z,z_0)} dz_0. \end{aligned}$$
(35)

Now, putting the value of T in the above equation, we get

$$\begin{aligned} V_2{(z)}&=\frac{G_2(z,H)}{A}\frac{L_1}{N_{02}}\left[ \frac{2}{N_1}G_1(H,H) +\frac{1}{N_{02}}\int _{H}^{\infty }I G_2(H,z_0)dz_0\right] - \frac{1}{N_{02}}\int _{H}^{\infty }I G_2(z,z_0)dz_0, \end{aligned}$$
(36)

where

$$\begin{aligned}I =R' \frac{d^2V_2}{dz_0^2}+\{\epsilon _1 +2\epsilon _2(z_0-H)\}\frac{d^2V_2}{dz_0}-R' \Theta ^2V_2.\end{aligned}$$

Now, \({V_2{(z)}}\) in equation (36) can be simplified by using the procedure of approximation. The obtained expression can be further used in the determination of \(V_1(z)\). In view of this, the first order approximation for \(V_2(z)\) can be derived by avoiding the higher powers of both parameters \(\epsilon _1\) and \(\epsilon _2\) present in equation (36) as

$$\begin{aligned} V_2{(z)}&=\left( \frac{G_2(z,H)}{A}\frac{L_1}{N_{02}}\right) \frac{2}{N_1}G_1(H,H). \end{aligned}$$
(37)

Equation (37) mathematically describes the displacement \(V_2(z)\) for any arbitrary position in the considered IDP half-space medium.

On account of equations (24), (34) and (37), the following expressions of \(V_1{(z)}\) can be obtained

$$\begin{aligned} V_1(z)&=\frac{2}{N_1}G_1(z,H) \frac{L_1}{AN_{02}}G_2(H,H)-\frac{2L_1}{A^2\gamma ^2N_{02}^2}G_1(z,H)G_1(H,H)\left( \int _{H}^{\infty }I G_2(H,z_0) dz_0 \right) . \end{aligned}$$
(38)

In order to find the final value of \(V_1{(z)}\), it is needed to calculate the value of \(G_1\) and \(G_2\), which is accomplished by solving the following differential equation

$$\begin{aligned} \frac{d^2S}{dz^2}-\alpha _1^2S=0,\,\,\,\,\ z \ne z_0 \end{aligned}$$
(39)

Let \(S_1=e^{\alpha _1z}\) and \(S_2=e^{-\alpha _1z}\) represent the two distinct solutions to the differential equation (39) which vanish at \(z=\infty\) and \(z=-\infty\), respectively.

Therefore, the solution of equation (39) for an infinite medium can be taken as \({(e^{-\alpha _1|z-z_0|}/2\alpha _1)}\).

Consequently, \(S_1\) and \(S_2\)are combined to represent Green’s function for upper homogeneous TIDP layer as46

$$\begin{aligned} G_1{(z,z_0)}=\frac{e^{-\alpha _1|z-z_0|}}{2\alpha _1} +Pe^{\alpha _1z}+Qe^{-\alpha _1z}, \end{aligned}$$
(40)

where P and Q are undetermined coefficients which can be found by following conditions

$$\begin{aligned} \frac{dG_1{(z,z_0)}}{dz}=0 \,\,\,\,\,\,\,\text {at both } z=0 \text { and } z=H.\,\ \end{aligned}$$
(41)

On using equation (41) into equation (40), we get

$$\begin{aligned} G_1{(z,z_0)}=-\frac{1}{2\alpha _1}[e^{-\alpha _1|z-z_0|}+\frac{e^{-\alpha _1z}(e^{\alpha _1(H+z_0)}+e^{-\alpha _1(H+z_0)})}{e^{\alpha _1 h}-e^{-\alpha _1 h}} + \frac{e^{\alpha _1z}(e^{\alpha _1(H-z_0)}+e^{-\alpha _1(H+z_0)})}{e^{-\alpha _1 h}-e^{-\alpha _1 h}}] \end{aligned}$$
(42)

Using the similar approach, the value of \(G_2{(z,z_0)}\) can be calculated from equation (31), we get

$$\begin{aligned} G_2{(z,z_0)}=-\frac{1}{2\alpha _2}[e^{-\alpha _2|z-z_0|} +e^{-\alpha _2(z+z_0-2H)}] \end{aligned}$$
(43)

Dispersion relation

The substitution of \(G_1(z,z_0)\) and \(G_2(z,z_0)\) from equations (42-43) into the equation (38) provides

$$\begin{aligned} V_1{(z)}&=\frac{L_1(e^{\alpha _1z}+e^{-\alpha _1z})}{\alpha _1\alpha _2AN_1N_{02}\sinh {(\alpha _1H)}}\left[ 1-\frac{\alpha _2}{A\gamma ^2N_{02}}\left( \frac{-1}{\alpha _1\tanh {(\alpha _1H)}}\right) \{\frac{1}{4\alpha _2^5}{(\alpha _2\epsilon _1+\epsilon _2)}{(\alpha _2^2+\Theta ^2)}\}\right] \end{aligned}$$
(44)

Since, we have \(\epsilon _1<<1\) and \(\epsilon _2<<1\), therefore, the higher powers of \(\epsilon _1\) and \(\epsilon _2\) can be ignored. In view of this, the above equation (44) can be further expressed as

$$\begin{aligned} V_1{(z)}&=\frac{L_1(e^{\alpha _1z}+e^{-\alpha _1z})}{\alpha _1\alpha _2AN_1N_{02}\sinh {(\alpha _1H)}}\left[ 1-\frac{\alpha _2}{A\gamma ^2N_{02}}\left( \frac{-1}{\alpha _1\tanh {(\alpha _1H)}}\right) \{\frac{1}{4\alpha _2^5}{(\alpha _2\epsilon _1+\epsilon _2)}{(\alpha _2^2+\Theta ^2)}\}\right] \end{aligned}$$
(45)

By applying the inverse Fourier transform mentioned in equation (11) into equation (45), we may further determine the displacement component at any position in TIDP layer as

$$\begin{aligned} {\bar{v}}_1{(x,z)}&=\int _{-\infty }^{\infty }\frac{L_1(e^{\alpha _1z}+e^{-\alpha _1z})}{\alpha _1\alpha _2AN_1N_{02} \sinh {(\alpha _1H)}}e^{-i\Theta x}d\Theta \nonumber \\ &\phantom { =}\left[ 1-\frac{\alpha _2}{A\gamma ^2N_{02}}\left( \frac{-1}{\alpha _1\tanh {(\alpha _1H)}}\right) \{\frac{1}{4\alpha _2^5}{(\alpha _2\epsilon _1+\epsilon _2)}{(\alpha _2^2+\Theta ^2)}\}\right] \end{aligned}$$
(46)

The integral value of equation (46) is completely influenced by the integrand’s poles, which are situated at the equation’s roots

$$\begin{aligned} 1+\frac{\alpha _2}{AN_{02}\alpha _1\gamma ^2}\left( \frac{-1}{\tanh {(\alpha _1H)}}\right) \left[ \frac{1}{4\alpha _2^5}{(\alpha _2\epsilon _1+\epsilon _2)}{(\alpha _2^2+\Theta ^2)}\right]&=0. \end{aligned}$$
(47)

On replacing \(\Theta\) by k and \(\alpha _1\) by \(i \alpha _1^{\prime}\), equation (47) yields

$$\begin{aligned} \tan {(\alpha _1^{\prime}H)}&= -\frac{{(\alpha _2\epsilon _1+\epsilon _2)}{(\alpha _2^2+k^2)}}{4\alpha _1^{\prime}\alpha _2^4A\gamma ^2N_{02}}, \end{aligned}$$
(48)

where

$$\begin{aligned}A= -\frac{\alpha _2 N_{02}+L_1\gamma ^2\alpha _1^{\prime}\tanh {(\alpha _1^{\prime}H)}}{\alpha _1^{\prime}\alpha _2N_{02}\gamma ^2 \tanh {(\alpha _1^{\prime}H)}},\,\,\ \alpha _1^{\prime}=k\sqrt{\frac{c^2}{\beta _1^2}-1} ,\,\,\,\,\ \alpha _2=k\sqrt{1-\frac{c^2}{\beta _2^2}}\end{aligned}$$

Equation (48) characterizes the dispersion relation for Love waves propagating in a homogeneous TIDP layer overlying a heterogeneous IDP half-space.

Particular cases

Case (1)

When the TIDP layer becomes the single poroelastic layer and IDP half-space becomes the single poro-elastic half space i.e (\(\frac{(\rho _{02}^{2})^{(l)}}{\rho _{22}^{(l)}}\rightarrow 0\)) in the considered structure, the expression for dispersion relation in (48) reduces to

$$\tan {(\alpha_1{^{{\prime}{\prime}}H})} = \frac{{N_{02}}{N_1{{{\alpha}_{2}^{\prime}}}} + ({{\alpha}_{2}^{\prime}}{\epsilon}_1 + {\epsilon}_2)({\alpha}_2^{{\prime}2} + k^2)/4{\alpha}_2^{{\prime}3}}{L_1^2{\alpha}_1^{{\prime}{\prime}}}$$
(49)

where

$$\alpha _1^{\prime}=\sqrt{k^2-(\omega ^2/\beta _1{'}^2)},\quad \beta _1^{\prime}= \sqrt{\frac{N_{1}}{\rho ^{(1)}{'}}},\quad \rho ^{(1)'}=\rho _{00}^{(1)}-\frac{(\rho _{01}^{2})^{(1)}}{\rho _{11}^{(1)}},$$
(50)

and

$$\alpha _2^{\prime}=\sqrt{k^2-(\omega ^2/\beta _2{'}^2)},\quad \beta _2^{\prime}= \sqrt{\frac{N_{02}}{\rho ^{(2)}{'}}},\quad \rho ^{(2)'}=\rho _{00}^{(2)}-\frac{(\rho _{01}^{2})^{(2)}}{\rho _{11}^{(2)}},$$
(51)

where (49) express the dispersion relation for love waves which propagates in a layered heterogeneous single poroelastic structure.

Case (2)

When the TIDP layer and IDP half-space in the considered layered double porous structure becomes the isotropic elastic layer i.e   \(L_1\) = \(N_1\), (\(\frac{(\rho _{0k}^{2})^{(1)}}{\rho _{kk}^{(1)}}\rightarrow 0\)) and the isotropic elastic half-space (\(\frac{(\rho _{0k}^{2})^{(2)}}{\rho _{kk}^{(2)}}\rightarrow 0\)) respectively, then the dispersion relation (48) becomes

$$\tan {(\alpha_1{^{'''}H})} = \frac{{N_{02}}{\alpha}_2^{''} + ({\alpha}_2^{''}{\epsilon}_1 + {\epsilon}_2)({\alpha}_2^{{''}2} + k^2)/4{\alpha}_2^{{''}3}}{N_1{\alpha}_1^{'''}}$$
(52)

where

$$\begin{aligned} \alpha _1{'''}=\sqrt{k^2-(\omega ^2/\beta _1{''}^2)},\,\,\,\,\ \,\,\ \beta {_1{''}}= \sqrt{\frac{N_{1}}{\rho ^{(1)}{''}}},\,\,\,\,\,\,\,\,\,\, {\rho ^{(1)}{''}}=\rho _{00}^{(1)}, \end{aligned}$$
(53)

and

$$\begin{aligned} \alpha _2{''}=\sqrt{k^2-(\omega ^2/\beta _2{''}^2)},\,\,\,\,\ \,\,\ \beta _2{''}= \sqrt{\frac{N_{02}}{\rho ^{(2)}{''}}},\,\,\,\,\,\,\,\,\,\, {\rho ^{(2)}{''}}=\rho _{(00)}^{(2)}. \end{aligned}$$
(54)

Equation (52) provides the dispersion relation for the Love wave’s propagation caused by the line source in the layered rock structure made up of isotropic elastic layer overlying the isotropic elastic half-space. The result obtained in (52) matches with the result (dispersion relation) obtained by Kumar et al19. for the case when the uppermost layer is removed \((H\,\rightarrow \,0\), H is the thickness of the uppermost layer) and the structure is comprised of only isotropic elastic layer overlying isotropic half-space.

Case (3)

When the considered layered double porous structure becomes the isotropic layered elastic structure i.e \(N_1=L_1\), (\(\frac{(\rho _{0k}^{2})^{(l)}}{\rho _{kk}^{(l)}}\rightarrow 0\)) and the isotropic elastic half-space is free from linear heterogeneity i.e \(\epsilon _1=0\), then the expression for dispersion relation (48) becomes

$$\tan {(\alpha_1{^{'''}H})} = \frac{{N_{02}}{\alpha}_2^{''} + {\epsilon}_2({\alpha}_2^{{''}2} + k^2)/4{\alpha}_2^{{''}3}}{N_1{\alpha}_1^{'''}}$$
(55)

Equation (55) provides the dispersion relation for the Love wave’s propagation caused by the line source in the layered rock structure made up of transversely isotropic double porous layer lying over the isotropic double porous half-space having quadratic heterogeneity only. The result obtained in equation (55) matches with the result (dispersion relation) obtained by Chattopadhyay et al35. after ignoring viscosity.

Case (4)

When the considered layered double porous structure becomes the isotropic layered elastic structure and the isotropic elastic half-space i.e , \(N_1=L_1\), (\(\frac{(\rho _{0k}^{2})^{(l)}}{\rho _{kk}^{(l)}}\rightarrow 0\)) is free from quadratic heterogeneity i.e \(\epsilon _2=0\), then the expression for dispersion relation in (48) becomes

$$\tan {(\alpha_1{^{'''}H})} = \frac{{N_{02}}{\alpha}_2^{''} + {\epsilon}_1({\alpha}_2^{{''}2} + k^2)/4{\alpha}_2^{{''}2}}{N_1{\alpha}_1^{'''}}$$
(56)

Equation (56) provides the dispersion relation for the Love wave’s propagation caused by the line source in the layered rock structure made up of isotropic elastic layer overlying the isotropic elastic half-space having linear heterogeneity only. The result obtained in equation (56) matches comparatively with the result (dispersion relation) obtained by Chattopadhyay et al47..

Case (5)

When the considered layered double porous structure becomes the isotropic layered elastic structure and the isotropic elastic homogeneous half-space i.e \(\epsilon _1=0\) and \(\epsilon _2=0\), (\(\frac{(\rho _{0k}^{2})^{(l)}}{\rho _{kk}^{(l)}}\rightarrow 0\)), then the expression for dispersion relation in (48) becomes

$$\begin{aligned} \tan \left( kH\sqrt{\frac{c^2}{\beta {_1{''}^2}}-1}\right) =\frac{N{_{02}}\sqrt{1-\frac{c^2}{\beta {_2{''}^2}}}}{N_1\sqrt{\frac{c^2}{\beta {_1{''}^2}}-1}}, \end{aligned}$$
(57)

where

$$\begin{aligned} \beta _1{''}= \sqrt{\frac{N_1}{\rho _1{''}}}\,\,\,\, \& \,\,\,\, \beta _2{''}= \sqrt{\frac{N_{02}}{\rho _2{''}}}. \end{aligned}$$
(58)

The expression obtained in equation (57) represents the propagation of Love wave in layered isotropic elastic structure made up of isotropic elastic upper layer overlying the isotropic elastic half-space. The expression obtained in equation (57) matches exactly with the standard Love wave’s equation1.

Numerical results and discussion

For the purpose of numerical simulations and graphical depiction, the data of the following physio mechanical properties have been considered:

For Upper Transversely Isotropic Double-Porous Layer44,45

$$\begin{aligned} N_1= & 8.3368\times 10^9 N/m^2,\,\,\,\,\,\ L_1 = 5.5211\times 10^9 N/m^2,\,\,\,\,\,\,\ \rho _{00}^{(1)} = 2282.4 kg/m^3,\,\,\,\,\ \rho _{01}^{(1)} = -291.2 kg/m^3\\ \rho _{11}^{(1)}= & 540.8 kg/m^3,\,\,\,\,\,\ \rho _{02}^{(1)} = -83.2 kg/m^3,\,\,\,\,\,\ \rho _{22}^{(1)} = 124.8 kg/m^3\end{aligned}$$

For Lower Isotropic Double-Porous Half-Space44,45

$$\begin{aligned} N_2= & 18.937\times 10^9 N/m^2,\,\,\,\,\,\ \rho _{00}^{(2)} = 2.8295\times 10^3 kg/m^3,\,\,\,\,\ \rho _{01}^{(2)} = -450.6840 kg/m^3\\ \rho _{11}^{(2)}= & 550.8360 kg/m^3,\,\,\,\,\ \rho _{02}^{(2)} = -13.4680 kg/m^3,\,\,\,\,\ \rho _{22}^{(2)} = 25.012 kg/m^3\end{aligned}$$

Note: Aggregate density \(d^{(l)}\) = \(\rho _{00}^{(l)}\)+\(2\left( \rho _{01}^{(l)}+\rho _{02}^{(l)}\right)\) + \(\rho _{11}^{(l)}+ \rho _{22}^{(l)}\), as for the layer \(l=1\) and for the half-space \(l=2\) respectively.

Fig. 2
figure 2

Scaled phase velocity (\(\frac{c_{ph}}{\beta _1})\) of Love waves against scaled wave number (kH) under the influence of linear heterogeneity parameter \((\nu _1=\frac{\epsilon _1 H}{N_{02}})\).

Fig. 3
figure 3

Scaled phase velocity (\(\frac{c_{ph}}{\beta _1})\) of Love waves against scaled wave number (kH) under the influence of quadratic heterogeneity parameter \((\nu _2=\frac{\epsilon _2 H^2}{N_{02}})\).

Fig. 4
figure 4

Scaled phase velocity (\(\frac{c_{ph}}{\beta _1}\)) of Love waves against scaled wave number (kH) under the influence of porosity parameter \(d_l=\frac{\rho ^{(1)}}{d^{(1)}}\) of the double porous upper layer.

Fig. 5
figure 5

Scaled phase velocity (\(\frac{c_{ph}}{\beta _1}\)) of Love waves against scaled wave number (kH) under the influence of porosity parameter \(d_h=\frac{\rho ^{(2)}}{d^{(2)}}\) of the double porous half-space.

Fig. 6
figure 6

Scaled phase velocity (\(\frac{c_{ph}}{\beta _1}\)) of Love waves against scaled wave number (kH) under the influence of anisotropic parameter (\(\gamma\)) of the double porous upper layer.

Fig. 7
figure 7

Comparison of case-wise phase velocity (\(\frac{c_{ph}}{\beta _1}\)) of Love waves against non-dimensional wave number (kH).

Figure 2 represents the effect of linear heterogeneity parameter \((\nu _1=\frac{\epsilon _1 H}{N_{02}})\) of a IDP half-space on the phase velocity \((\frac{c_{ph}}{\beta _1})\) of Love waves in the layered double porous structures containing TIDP layer lying over IDP half-space. It is noted from Fig. 2 that as the linear heterogeneity parameter \((\nu _1)\) increases, the Love wave’s phase velocity decreases in the aforementioned structure.

Figure 3 shows the influence of quadratic heterogeneity parameter \((\nu _2=\frac{\epsilon _2 H^2}{N_{02}})\) of a IDP half space on the phase velocity \((\frac{c_{ph}}{\beta _1})\) of Love wave. It can be seen from the Fig. 3 that as the quadratic heterogeneity parameter \((\nu _2)\) of IDP half-space rises, the Love wave’s phase velocity in the structure under consideration falls.

It is observed from both figures that the heterogeneity parameters (both linear and quadratic) have adverse impact on the phase velocity of Love waves in the considered double porous structures. It may be due to the reason that as the heterogeneity (associated with the heterogeneity parameters) increases, the size of the grains of the rock materials varies. It may lead to the gapping at some positions which further may results into the decrease of phase velocity of Love waves.

Figure 4 demonstrate the impact of porosity parameter of TIDP layer \((d_l=\frac{\rho ^(1)}{d^(1)})\) on the phase velocity \((\frac{c_{ph}}{\beta _1})\) of Love wave. It’s been noticed from Fig. 4 that as the porosity parameter of the TIDP layer increases (i.e., porosity decreases), the Love wave’s phase velocity increases significantly in the studied double porous rock structure. Now, since most of the Love wave’s energy is centred on the finite layer, so, in the examined rock formation, the phase velocity of the Love waves increases as the upper layer becomes more compressed (i.e. porosity decreases).

Figure 5 explains the effectiveness of porosity parameter of a IDP half-space \((d_h=\frac{\rho ^(2)}{d^(2)})\) on the phase velocity \((\frac{c_{ph}}{\beta _1})\) of Love waves. It has been noticed from Fig. 5 that as the porosity parameter of the IDP half-space increases (i.e., porosity decreases), a significant reduction has occurred in the Love wave’s phase velocity. Probably because, due to the heterogeneity in the half-space, a feeble compaction has been found comparatively in the lower half-space. Also, as the energy of Love waves centred on the finite layer, only a minor portion of it propagates in the lower half-space which could be the reason behind the drop in phase velocity of Love waves.

Figure  6 exhibits the impact of anisotropic parameter (\(\gamma\)) of TIDP layer on the phase velocity \((\frac{c_{ph}}{\beta _1})\) of Love waves. It has been revealed From the Fig. 6 that as the anisotropic parameter (\(\gamma\)) increases, the Love wave’s phase velocity decreases considerably in the regarded double-porous rock structure, it happens because as the ratio of anisotropic parameter (\(L_1\)/\(N_1\)) associated with the TIDP layer increases, the energy of Love waves gets dispersed in many direction, causing the Love wave’s phase velocity reduces.

Figure 7 exhibits a comparative study of the effect to the phase velocity \((\frac{c_{ph}}{\beta _1})\) of Love wave in different kind of scenario. The figure demonstrates that the phase velocity is maximum for case (5) i.e for the case of isotropic layered elastic structure and the isotropic elastic homogeneous half-space and minimum for case (3) i.e for the case of isotropic layered elastic structure and the isotropic elastic half-space is free from linear heterogeneity.

Conclusion

This study investigates Love wave propagation in a heterogeneous double-porosity structure influenced by an impulsive line source. The considered structure comprises a homogeneous transversely isotropic double porous (TIDP) layer overlying a heterogeneous isotropic double porous (IDP) half-space with depth-dependent stiffness. Green’s function technique yields the Love wave dispersion relation, and the influence of heterogeneity, porosity, and anisotropy on phase velocity is analyzed. Numerical simulations and graphical representations are included. The single-porosity case is derived as a special instance. Validation is achieved through comparisons with existing literature. Key findings are presented subsequently.

  1. 1.

    The increase in the linear heterogeneity parameter of the IDP half-space causes reduction in phase velocity of Love wave’s in the double porous rock structure.

  2. 2.

    The IDP half-space’s non-linear quadratic heterogeneity parameter adversely affects Love wave phase velocity in the double porous rock structure.

  3. 3.

    The augmentation in the porosity parameter (i.e decrease in the porosity) of the TIDP layer causes the upsurge in the phase velocity of Love waves in the studied double porous structure. On the contrary, the increment in the porosity parameter of the IDP half space significantly resist the Love wave’s phase velocity.

  4. 4.

    The phase velocity of Love wave’s propagation decreases with the enhancement of anisotropic parameter (associated with the extent of anisotropy) of the TIDP layer in the double porous rock structure.

Engineering applications

The regions near the reservoir sometimes encounter the seismic activity. The impulsive line source in the rock materials can significantly affect the behavior of seismic wave propagation. The results obtained from the current study can have significant possible applications to examine the dynamic strength of the rock structure that can be found near the reservoirs. From the current analysis, the medium parameters that are responsible for the rise and fall of phase velocity can be examined. It can be useful to select the appropriate rocks for the construction of more safe engineering structure in the earthquake prone area nearby pools, pond, and reservoirs regions.

It has been examined from the obtained results that higher the porosity parameter (i.e. lower the porosity) of upper TIDP layer causes the increment in the phase velocity of Love waves, while an increase in the anisotropy of the upper TIDP layer is associated with a decrease in the phase velocity of Love waves. On the other side, if the porosity is less and the heterogeneity is more in the lower half-space medium then the Love wave’s velocity reduces in the considered double porous structure. In view of above, the double porous rock blocks for the construction purpose have been selected in such a way that its upper portion is of higher porosity and higher anisotropy while the lower bigger portion is of lower porosity and higher heterogeneity. Due to this, the velocity of Love waves (i.e the damaging effect of Love waves) can be reduced. In this way, the engineering structures influenced by Love waves near the reservoir region can be constructed with the more efficient and safe design.