Seismic response of a layered soil deposit with inclusions

It is known that soil massifs can amplify or weaken seismic waves generated by earthquakes. Therefore, the problem of studying the impact of soil deposits on the passage of seismic waves is important in terms of the facilities in operation and the design of new earthquake-resistant objects. Soil deposits, which are allotted for building, are mainly layered. In addition, the materials in these layers are also significantly heterogeneous. To describe the dynamics of inhomogeneous soil massif, the model of an elastic continuum with oscillating non-interacting inclusions is used. Within the framework of this model, the resonant properties of multi-layered soil deposit are analyzed at the conditions of harmonic perturbations applied to the bedrock. On the basis of the solution to the boundary value problem concerning oscillations of the system subjected to the free surface and conjugation conditions on the boundaries between layers, it is derived the transfer function which characterizes the amplification of shear displacements by the layered system. Within the framework of problems on the oscillations of twoand five layered systems, the analytical studies were confirmed by numerical evaluations of transfer functions. In particular, using the built-in functions of the system «Mathematica», it is developed the numerical procedure for evaluating the frequency dependencies of amplification factor for layered Kelvin—Voigt media and media with oscillating inclusions. Moreover, for the two-layered system, it is analyzed the effect on the transfer function for the ratio of layers’ shear moduli and the ratio of the inclusions’ natural frequencies. It is also shown that the maxima in the transfer function correspond to the eigenfrequencies of the boundary value problem.The obtained results and the proposed approach to the study of the response of the layered inhomogeneous medium to vibrational perturbations can serve as a theoretical basis for earthquake-resistant design and construction.

Analysis of the destructive effects of past and recent earthquakes [Ishihara,1996;Adimoolam, Banerjee, 2019] has shown that the pattern of damage during earthquakes is mainly determined by the reaction of local soils to seismic loading. The near-surface layers of soil strata act as a filter that amplifies/weakens the amplitude of seismic waves.
Since the amplification of seismic oscillations by soil deposit can be very significant [Kausel, Roёsset, 1984;Pratt et al., 2017;Adimoolam, Banerjee,2019;Kumar et al., 2020], the analysis of the soil response became one of the most important tasks of engineering seismology [Gutowski, Dym, 1976]. The soil reaction modeling allows one to determine the dominant frequencies of a local ground site, to estimate the enhancement of seismic oscillations by local soils [With, Bodare, 2007;Hosseini, Pajouh, 2010;Kumar et al., 2020], and to obtain accelerograms and response spectra of soil oscillations on the surface. Information on the probabilistic characteristics of seismic oscillations on the soil surface is required for evaluating the critical dynamic stresses and strains that cause loss of soil DOI: https://doi.org/10. 24028/gzh.v43i2.230186 stability and, consequently, the destruction of buildings and structures [Gutowski, Dym, 1976;Takemiya, Yamada, 1981;Mandal et al., 2012] located on such ground.
To assess the seismic ground response of the local area , it is used the well-known [Sarma, 1994;Ishihara, 1996;Kramer, 1996] approach based on the shear waves which propagate vertically upwards from the underlying bedrock through the layers of ground. In engineering practice, as a rule, linear or equivalent linear methods are used with the involvement of SHAKE [Schnabel et al., 1972], EERA [Bardet et al., 2000], Deepsoil [Hashash, 2012], and others. The algorithms [Gutowski, Dym, 1976;Sarma, 1994;Hosseini, Pajouh, 2010] incorporated in these programs are based on the assumption that all boundaries between the layers are horizontal and each layer is homogeneous.
Another way to improve the description of soil state is related to the modification of equations of motion by means of introducing the auxiliary equations for the additional degrees of freedom. The bright example of this approach is the Cosserat model where the rotational degrees of freedom are incorporated [Erofeev, 2003;Green, Rivlin, 1964]. In this paper, we use the generalization of classical Lame equations [Slepjan, 1967;Palmov, 1969;Mishuris et al., 2019] taking into account the oscillating degrees of freedom. To do this, two continua are considered. The carrying medium is classical, while the other consists of a set of partial oscillators representing the medium's inclusions. Such models and their applications for the investigation of wave dynamics were studied in the papers [Danilenko, Skurativskyy, 2008;Danylenko, Skurativskyi, 2012, 2016.
The aim of the present studies is to evaluate the resonant properties of multilayer soil systems using a generalized model with oscillating inclusions.
Analytical estimation of soil response to the harmonic disturbances. The assessment for seismic response of heterogeneous soil deposit (single layer) on the basis of the model of mutually penetrating continua has been developed in [Kendzera et al., 2020 b]. The equations of motion for the one-dimensional shear deformation of the single layer soil deposit read as follows where u (z, t) and w (z, t) are the horizontal displacements of carrying medium and attached oscillator respectively, and G is the shear modulus. The quantity Ω stands for the natural frequency of a partial oscillator. Assume also that the relaxation time τ is constant for all oscillators. The function m (x) describing the distribution of oscillators over their natural frequencies is chosen as follows is the Dirac delta function and the parameter m coincides with the ratio of densities of continua.
It should be noted that in the case of delta distribution, model (1) has no integral term and is reduced to the model with oscillators possessing the identical natural frequencies.
In this work, we deal with the generalization of the approach developed in [Kendzera et al., 2020 b] to the case of layered soil deposit.
Thus, we consider the layered medium of thickness H (Fig. 1), which is composed of N layers and rests on a rigid base. Each layer is of thickness h s and contains inclusions (internal where u (z = H ; t) is the deformation at the depth H, and coincides with the value of displacement at the ground surface.
In contrast to the problem on the single soil layer, the problem on the multilayered medium contains auxiliary conditions related to the conjugation of displacements and stresses at the boundaries between layers. In this work, we use the physically motivated conditions of smoothness of deformations and stresses [Kramer, 1996;Kundu et al., 2019;Mondal et al., 2020] ( ) Applying the approach described in works [Kramer, 1996;Kendzera et al., 2020 b], the dependency of amplification factor of layered soil deposit on the frequency ω is constructed. To do this, let us derive the solution of model (1)-(3), which relates to the regime of deposit's oscillations with the constant amplitude and frequency ω.
In each layer, the partial solution can be found in the following form . From the condition (4) it follows that oscillators). In what follows, in each s-th layer the medium is described by the equation of motion (1), so that As shown in Fig. 1, each layer is subjected to the own reference frame with the origin placed at the upper layer boundary and the axis directed downwards. The typical assumption [Sarma, 1994;Kramer, 1996;Kundu et al., 2019;Mondal et al., 2020] about the absence of stresses at the free ground surface is used: We also assume that the bedrock is deformed according to the harmonic law where Q and ω are the amplitude and frequency of external disturbances respectively. Since we are interested in resonant properties of soil deposit, then it is suitable to assume that disturbances at the bottom edge are of unity amplitude, i. e. Q = 1. The value of amplification factor F for the deposits defined as a ratio of displacement amplitudes at its edges [Kramer, 1996;Rivin, 2003;With, Bodare, 2007] For comparison, it's convenient to take the Kelvin-Voigt medium [Sarma, 1994;Kramer, 1996] obeying the equation of motion where / η = ξ ω 2 stands for viscosity. The amplification factor F is defined by the same formula (8) , but in this case , we should put , s = 1, 2 . Numerical simulation of soil response to the harmonic disturbances. To check the analytic studies and develop the procedures for the more general form of functional dependencies of soil characteristics on depth, the numerical approach for evaluation of amplification factor is utilized. During numerical simulations tools incorporated into «Mathematica» system are used. To apply numerical methods, the system (2)-(4) should be adapted. At first, let us fix the reference frame (the point O is the origin, Oz is the axis) as it is shown in Fig. 1. We replace the set of systems with the conjugation conditions at the boundary between layers by the single system with discontinuous coefficients: Using the solutions (4), the conditions at the discontinuities can be written in the form of algebraic equation [Kramer, 1996] Solving the together with relation (7) and , one obtains the following values of coefficients A s and B s [Kramer, 1996]: 1 . It's worth noting that from relations obtained it follows that the coefficients are the function of the quantity A only, i. e. A s = A s (A) and B s = = B s (A). In turn, equation (7) allows one to derive the value of A. Thus, the evaluation of amplification factor F can be realized via the formula The resonant curve for the two-layered media. The main features of the analytical treatments presented above can be outlined for the two-layered soil deposit. In this case, it is easy to derive the analytic expression for the quantity .
The functions describing the soil characteristics with depth are defined as follows where ( ) Θ  is the Heaviside step function. Next, the discontinuous function G (z) should be smoothed out. To do this, we apply the following procedure: where the parameters G 1, 2 define the asymptotic values of function ( ) G z  , h represents the coordinate of discontinuity, ε is the rate of smoothness (the function ( ) G z  tends to a step function as ε → ∞ ). When the function G (z) possesses several discontinuities, the similar procedure can be used. Let h j , j = 1, ... K be the coordinates of discontinuities, then ( ) So, the system (10) with prescribed boundary and zero initial conditions is solved by the command NDSolve [...]. The solution is evaluated at a long time interval to provide that all transient processes are passed and the system reaches the steady oscillating regime. We are interested in the solution at free surface (z = 0). Namely, it is evaluated the amplitude of free surface oscillations when the steady mode is realized. This amplitude coincides with the amplification factor F due to its definition we are used. Finally, varying the frequency ω and deriving the corresponding values of F, the resonant curve is constructed.
To check the properties of numerical pro-cedure, we consider the Kelvin-Voigt twolayer (N = 2) model (9) with the boundary condition (3), (4). The parameters are fixed as follows: and, according to relation ( 11 ) , the function . The Fig. 2, a represents the results of the numerical derivation of resonant curve (points). The analytical expression (8) of amplification factor F is drawn by solid curve which fits the points well. The same numerical simulation was repeated for the five layered soil deposit which is characterized by the same power dependencies (13) but now N = 5. The resulting resonant curve is depicted in Fig. 2 The similar numerical procedure is used for the estimation of F in the case of twolayered medium with oscillating inclusions. At first, we consider the layers when only their shear module G are different. So, assume that For comparison , the analytically derived function F is drawn by the solid line. It is evidence of the perfect fitting of both approaches. ] ω ∈ 0.2, 4 , and the amplification factor for the first resonant frequency ω 1 is much higher than for the second resonant frequency ω 2 . Analytical expression ( 8 ) for the amplification factor in the case of two-layer application allows us to analyze how it is affected by the relationship between the parameters of the layers. The contours of the surface ( ) F F = γ ω , shown in Fig. 4 indicate that the arctan resonant frequencies depend on the ratio of the shear moduli for the second and first layers / G G γ = 2 1 : when this ratio γ increases, the resonant frequencies shift towards their increase. In addition, in the case of the ratio γ = 1, the amplification factor F for the first resonant frequency is much higher than for subsequent resonant frequencies. As the ratio γ increases, the amplification factor becomes more uniform for all resonant frequencies. It should be noted that only cases when γ ≥ 1 are analyzed, because the second layer, being under the first layer, is denser, and hence G G > 2 1 . It is also important to analyze the influence of the ratio of natural frequencies of par-tial oscillators / θ = Ω Ω 2 1 on the dependence ( ) F F = ω . Fig. 5 shows that there is the range [ ] θ ∈ 1, 2 that separates the two areas with the different nature of the dependencies ( ) F F = ω . In the area [ ] θ ∈ 0, 1 the maxima of the function F (ω) are shifted toward the low frequencies relative to the maxima in the area [ ] θ ∈ 2, 11 . In the area [ ] θ ∈ 1, 2 there is almost none except one natural frequency. It should also be noted that in the aria [ ] θ ∈ 2, 11 there is no shift of the maxima of the function F (ω) with increasing the ratio / θ = Ω Ω 2 1 . The construction of the amplification factor is performed also for the five-layered soil deposit with oscillating inclusions. Used parameters are as follows:  is observed as well (see Fig. 6).
Relation to the boundary value problem. It is worth noting that the estimation of factor F and solving the boundary value problem are related tasks. To show this, consider the system (10) written for the two-layer soil deposit. It is easy to show that the boundary conditions can be reduced to homogenous (zero) by the change of variable are the step functions and using the solution (6), we derive the solution in the first layer  cies ω 1 and ω 2 gives the values 0.698 and 2.826, which coincide perfectly with the maxima in the profile of F (see Fig. 3, a). The similar procedure can be developed for the estimation of eigenfrequencies of multilayered deposit. Concluding remark. The presented research is concerned with the assessment of resonant properties of multi-layered non-homogenous soil deposits. To describe the dynamics of these media, the mathematical model in the form of mutually penetrating continua was used. The procedure for amplification factor evaluation was adapted on the basis of this model. The good agreement of numerical and analytical results on the amplification factor evaluation was obtained. For a two-layer deposit, the influence of the ratios of the shear moduli and the natural frequencies of the oscillatory inclusions for the two layers on the amplification factor is analyzed. Moreover, we have shown that maxima in the amplification factor profile correspond to the eigenfrequencies of boundary value problem.
The obtained findings contribute to the improvement of methods for response assessment of soil deposits to seismic waves. These results and the proposed theoretical approaches provide increased protection of artificially created objects from seismically dangerous influences [With, 2007]. The related branches of science [Li, 2000;Smolyakov et al., 2010;Mondal et al., 2020] can be interested in the results as well.