Scattering of SH waves by a scalene triangular hill with a shallow cavity in half-space

Sun Yingchao ,Yang Zailin ,Yang Yong and Wang Guangyi

1.College of Aerospace and Civil Engineering,Harbin Engineering University,Harbin 150001,China

2.AVIC-HONGDU,Nanchang 330024,China

3.Key Laboratory of Advanced Material of Ship and Mechanics,Ministry of Industry and Information Technology,Harbin Engineering University,Harbin 150001,China

4.Engineering Research Institute of Appraisal and Strengthening,Shandong Jianzhu University,Jinan 250101,China

5.Jiangxi Technology Innovation Center of Advanced Aviation Structure and System Integration,Nanchang 330024,China

Abstract: In this study,a theoretical approach is used to investigate the scattering problem of circular holes under a scalene triangle on the surface.The wave displacement function is obtained by solving the Helmholtz equation that meets the zero-stress boundary conditions by adopting the method of separation of variables.Based on the complex function,multi-polar coordinate method,and region-matching technique,algebraic equations are established at auxiliary boundaries and free boundaries conditions in a complex domain.The auxiliary circle is used to solve the singularity of the reflex angle at the triangle corner.Then,according to sample statistics,the least squares method is used instead of the Fourier expansion method to solve the undetermined coefficient of the algebraic equations by discrete boundary.Numerical results show that the continuity of the auxiliary boundaries and the accuracy of the zero-stress boundaries are adequate,and the displacement of the free surface and the stress of the circular hole are related to the shape of the triangle,the position of the circular hole,the direction of the incident wave,and the frequency content of the excitation.Finally,time-domain responses are calculated by FFT based on the frequency domain theory,and the results reveal the wave propagation mechanism in a complicated structure.

Keywords: earthquake ground motions;theoretical seismology;wave scattering and diffraction;wave propagation

1 Introduction

The scattering of waves by various convex or concave shapes on the surface has always been an important research topic in the field of wave motion.The approaches to solve these types of problems can be divided into two categories: theoretical analysis methods and numerical methods.

As to the concave topographies like canyons and alluvial valleys,many mature research results have been published.Pioneering work in this area was done in the early 1970s.The wave function expansion method was used to analyze the local semi-cylindrical canyon scattering in the half space under the action of incident SH waves (Trifunac,1973).The Hermite function and mapping function were applied to analyze scattering of SH waves with arbitrary shape depressions (Liu and Han,1991).Subsequently,the introduction of Graf’s addition theorem provided an analytical solution to the SH waves scattering by a cylindrical canyon of circular-arc cross-section (Yuan and Liao,1994),and the weighted residual method was applied to 2D canyons of arbitrary shape (Lee and Wu,1991),The degenerate kernels and Fourier series expansions were adopted in the null-field integral equation to analyze the surface motion of multiple alluvial valleys for incident plane SH waves (Chenet al.,2008,2017).In recent years,based on the region-matching technique,many more complex concave topographies have been theoretically studied,such as V-shaped canyon,circular sectorial canyon,and deep semi-elliptic canyon with a horizontal edge (e.g.,Tsaur and Chang,2008;Tsauret al.,2010,2018;Zhanget al.,2012,2015;Changet al.,2013;Liuet al.,2019).

Due to the multiple reflections of the incident wave on the convex surface,the scattering of the incident wave by the convex surface is more complicated than that on the concave surface.The two-dimensional model is divided into a closed circular region and a half space with a semi-circular concave topography by introducing a semi-circular auxiliary boundary,to study the scattering of plane SH waves by a semicylindrical hill in the half-space (e.g.,Leeet al.,2006;Yuan and Men,1992;Todorovskaet al.,2001;Yuanet al.,1996).The anti-plane response of an isosceles triangle,a scalene triangle,and an isosceles trapezoidal hill to SH waves is studied through the construction of complex displacement fields and complex coordinates(Qiu and Liu,2005;Linet al.,2010;Yanget al.,2019).Application of a series of expansion and construction complex displacement fields to solve scattering of plane SH waves by a circular-arc hill with a circular tunnel has been carried out (Lianget al.,2004).In addition,numerical methods have been used to solve the scattering of SH waves by concave and convex topography.A hybrid method based on a combination of transfinite interpolation and series expansion solves the problem of irregular surface scattering (Shyuet al.,2014).A hybrid method based on the combination of Lamb series and the finite element method was also used to investigate a dike with trapezoidal structure and a circular-arc foundation embedded in an elastic half-space (Shyuet al.,2017),and the boundary element solution was used to evaluate the geometric effects of the canyon site (Panjiet al.,2013;Tarinejadet al.,2019;Liuet al.,2020;Kumar and Chakraborty,2020).From the 1990s until now,the combination of terrain and underground structure has received much attention.Scattering of SH waves by a semicircular canyon above a subsurface unlined tunnel,circular-arc hill with a circular tunnel,and V shaped canyon with a circular structure was studied by Graf′s addition theorem (e.g.,Lee and Manoogian,1995;Leeet al.,1999;Lianget al.,2004;Gao and Chen,2018;Chen and Zhang,2019;Panji and Habibivand,2020).

In engineering,both surface motion and hole stress concentrations have increasingly drawn attention,especially for more complicated structures of the interaction between topography and hole (Chenet al.,2016),such as the problems in seismic design of subsea tunnels,earthquake-proof building problems in mountainous residential areas with drainage tunnels or underground caves,and vibration problems of irregular boundary plates with holes.For example,according to seismic design code,irregular buildings and highrise buildings need acceleration or displacement time history curves to evaluate strength by the time history analysis method.Terrain and underground structures have a significant amplification effect on the surface time history curve;therefore,the complex structure transfer function is used to obtain seismic time-domain spectroscopy,which provides indispensable research content for seismic analysis of building structures.

In this study,based on existing research results,the SH wave scattering problem of shallow circular holes under scalene triangles is studied.A more flexible nonsemicircle region division method is used to solve the shallow hole problem and an auxiliary circle is used to solve the singularity of the reflex angle at the triangle corner,which was proposed by Achenbach (1970).The wave displacement function is obtained by solving the Helmholtz equation that meets the scalene triangle zero-stress boundary conditions by the separation of variables method.Based on complex coordinates and the symmetric method,suitable wave functions that satisfy the governing equation and zero-stress free boundary condition are constructed.Since there are many wave equations of sub-regions which correspond to coordinate systems,multiple complex coordinate systems are involved to transform the coordinate system of each region wave field expressions.In particular,for the difference in the range of wave function series and multiple auxiliary boundary continuous conditions,a more effective least squares method which only need discrete boundary is proposed.The discrete point spacing and boundary equation amplitude adjustments are used to coordinate Euclidean distance weight.After numerical simulation,the high accuracy of the auxiliary boundary continuity and comparisons with FEM results prove the accuracy of region-matching technique,wave function equation and least squares method.Finally,the effects of different angles of incidence,frequency content of the excitation,slopes of triangular edges,and positions of hole on surface motion and hole stress are discussed in the frequency domain,and the process of wave propagation and scattering around the triangle and shallow circle are shown in the time domain.This study has proposed using an auxiliary circle to solve the singularity of the reflex angle to improve convergence,and a direct and effective discrete boundary least squares method,which solves the precision problem of complex boundary conditions of scattered waves.

2 Methodology

2.1 Model

The model,shown in Fig.1,consists of a scalene triangle on a half space and a hole under a free surface,whereOstands for the triangular peak,C1andC2symbolize hypotenuse with the gradients of 1:n1and 1:n2,O2is the circle center with rounded edgeD2,andSrepresents a half space free flat surface.The SH-wave propagates in the linear elastic half space with shear modulusμand densityρat an incident angleα.

Fig.1 Model of scalene triangle with a hole

2.2 Equation of motion

To solve the singularity of the reflex angle at a triangle corner and obtain the global displacement function which satisfies the Helmholtz equation and complicated boundary conditions,a region-matching technique and auxiliary circle are used.The space can be divided into regionby auxiliary boundaryD1with a circle centerO1,auxiliary boundaryD3with a circle centerO3and so on.P,P1andP2are the projection ofO,O1andO2on the flat surface in region.The coordinate systems are established as shown in the figure,where the x-axis bisects the angle of the triangle corner.More details about each angle are presented in the Appendix.

In homogeneous,isotropic,and linearly elastic medium,the steady-state out-of-plane motions are required to satisfy the governing Helmholtz equation:

The time-harmonic factor e-iωtis omitted,and the equation can be simplified into the following form:

Correspondingly,radial stress and hoop stress have the forms of

2.3 Wave function in region

There is only standing waveWD()in the closed regionit needs satisfy governing Eq.(2) and the free hypotenuse condition.

Zero-stress condition of free surface hypotenuseC1andC2can be expressed as:

In polar coordinates,the wave equation is

The separated variable method,introducingw=R(r) Θ(θ),Eq.(7) can be solved:

A standing-wave function which satisfies the zerostress condition at the triangular edges is written as:

In the complex plane,corresponding to coordinate systemoxy,according to Eq.(10),the standing waveWD3(1)satisfying governing Eq.(2) and boundary condition Eq.(6) can be written as

whereW0is the displacement amplitude,which is supposed to be unity in this study.Cmis a coefficient to be determined.Jmp0() is the Bessel functions with mp0th order,p0=π /(α1+α2).InWD3(1)characters,superscripts(1)means region,superscriptsD3 represent auxiliary boundaryD3.K1is the shear wave number of region.The following symbol marking method is similar.

According to moving coordinates,Z0′can be expressed as

whereq0=(α2-α1)/2 is the angle betweenoxand the vertical line,andb03=H-H3+(Htan(α1)-H3tan(α4))i is the complex coordinates ofO3with origin atO.

In the complex plane (z,),Eq.(11) can be rewritten as

The corresponding shear stresses are:

The corresponding shear stresses are:

2.4 Wave function in region

Based on complex coordinates and the symmetric method,suitable wave functions are constructed to satisfy the governing equation and zero-stress free boundary condition.A semi-infinite space scattering field that satisfies boundary conditions according to symmetry is also constructed.The scattering wave equation with two symmetrical holes is

According to Eq.(23),the equation of scattered waveWD1(2)generated by boundaryDand satisfying governing Eq.(2) and free boundary co1nditionSin the complex plane (zj,) can be written as

Fig.2 Circular hole symmetry in semi-infinite space

Similarly,equation of scattered waveWD2(2)generated by boundaryD2and satisfying governing Eq.(2)anzd free boundary conditionSin the complex plane(j,) can be written as

andb71=H1-L2i,b72=H2.

In the above formula,W0is the displacement amplitude.DmandEmare coefficients to be determined.() is the Hankel functions of the first kind withmth order.

The corresponding shear stresses are:

In the above stress formula,q=0,and see Appendix in details.

Incident and reflected waves with incidence angleα,can be represented in the Cartesian coordinate systemo6x6y6,and be rewritten in the complex plane (z6,)

The corresponding shear stresses are:

2.5 Wave function in region

The corresponding shear stresses are:

2.6 Boundary conditions and solution method

Based on the continuity conditions of displacement and stress at the auxiliary boundaryD1,D3,D4,D5,and radial zero-stress at the hole edgeD2,a system of equations is established for solving the unknown complex coefficients.

Currently,the Fourier expansion method is commonly used to solve the undetermined coefficients of algebraic equations,and is an average approximation of the entire boundary conditions.Due to the wave field high gradient of the triangle edge and the auxiliary boundary corner point,the Fourier expansion method convergence speed is slow,which makes it difficult to solve the problem of the scalene triangle.Therefore,this study proposes the least squares method with direct discrete boundary conditions.Discrete points are taken according to a set distance on the boundary,and the displacement and stress on the two sides of the discrete points are equal,as shown in the figure.An infinite number of pointsncan be taken on the boundary to form an infinite number of equations to solve the undetermined coefficientsCm,Dm,Em,….To minimize the error of the undetermined coefficient of the finite term,a large number of sample pointsn(n ＞＞ m) are approximated to the true solution by the least squares method.This study uses nonequidistant discrete points and stress terms divided byμk(not reflected in the formula) to coordinate the weights of Euclidean distance.The matrix is expressed as

See the Appendix for details.

2.7 Surface displacement amplitude and hole stress

Equations (48) to (52) can also be expressed as

where |Wj| is the displacement amplitude,andφjis the phase angle ofWj

The dimensionless frequency of incident waves can be expressed as

Fig.3 Discrete points of auxiliary boundary and hole edge

wherek=K1=K2=K3=K4=K5for homogeneous medium,andμis the same for different regions.λis the wavelength of the incident waves.It is well known that the effect of elastic waves on surface displacement and hole stress highly relies on wavelength.As can be seen from Eq.(52),the dimensionless frequencyηis introduced to represent the ratio of the radius (r2) of the hole to the wavelength,and indirectly represents the magnitude of the wave number.

The dimensionless hoop stress is

whereτ0=μ2k2W0.

3 Numerical examples and discussions

Typical structural dimensions studied in this chapter are shown in Table 1.

Table 1 Shape parameters of scalene triangle and hole

3.1 Precision discussion

An important method to verify the theory is to verify convergence of the seriesm,and fit the continuity of the auxiliary boundary and the zero-stress condition of the circular hole edge and free surface.The convergences of the displacementWat surface point 2y/Lare already very good form=25 after trial calculations from Fig.4 andFig.5.The numerical results of the auxiliary boundary and zero-stress condition are given for frequency of incidentη=2.0, incident angleα=0°,90°,triangular edge slopen1=0.5,n2=1.5,m=30.As shown in Fig.4 and Fig.5,the displacementWand stressτrzcontinuity of the auxiliary boundaryDare good,and the stress of the circular hole boundaryD2and free surface are close to 0,which also indicates that the wave function and the least squares method are effective.

Fig.4 Continuity of auxiliary boundary and free edge zero-stress at α=0°

Fig.5 Continuity of auxiliary boundary and free edge zero-stress at α=90°

Another important method to verify the theory is to compare the solution results of the finite element method as shown in Fig.6,which displays the free surface displacement amplitudes |Wj|,hole edge stressτθz*and displacement cloud at a certain time.The FEM results are obtained by the commercial software LS-DYNA.The geometric model is meshed by a shell element with an edge length of 0.1 and grid only with out-of-plane translational degrees of freedom;the mesh area is large enough to eliminate the effects of boundary reflections.Sine excitation is applied to the bottom or left edge of the analysis area,corresponding to incident angle 0° or 90°,and the calculation time is long enough to ensure that it is in a steady state.The surface displacement magnitude is measured from the displacement of the surface element nodes,and the stress of the finite element hole edge comes from the stress of the nearest element of the hole edge,which represents the stress in the hole edge area rather than the hole edge.Therefore,the FEM stress results may be somewhat inaccurate.The results by the proposed method shows good agreement with those obtained by the finite element method from Fig.6.

Fig.6 Comparison of proposed solution results with FEM results at η=0.5,n1=0.5

3.2 Parameter studies in the frequency domain

Each position of free surface can be expressed by dimensionlessy/(L/2) in the Cartesian coordinate systemo6x6y6,where-1 represents the left triangle foot point,(n2-n1)/(n2+n1) is the triangle vertex,and 1 represents the right triangle foot point.

Free surface displacement amplitude of the triangle without hole and with hole is given in Fig.7.It shows that the hole has a significant impact for free surface displacement.Free surface displacement amplitude and circular hole hoop stress distribution at different incident frequencies,angle and slope are shown in Figs.8 and 9.When the incident wave has a horizontal component(α≠ 0°),the amplitude and oscillation frequency of the right-side surface (2y/L＞ 1.0) both decrease,while they increase on the incident wave side (2y/L＜-1.0).It appears that the triangle and the circle hole become filters and amplifiers,and are related to the angle of incidence and the slope of the triangle.This is due to the multiple reflections of the incident wave at different angles on different slopes.Both the amplitude and oscillation frequency of the triangular area (-1.0 ＜ 2y/L＜ 1.0) are improved to varying degrees.The circular hole hoop stress increases significantly with increasing incident frequency,and the distribution is closely related to the triangular slope.In particular,the incident angleα=0,due to the symmetry of the structure,the displacement amplitude and hole hoop stress distribution show symmetry,including the symmetry slope ofn1=0.5 andn1=1.5.In addition,the displacement amplitude and hole stress on the slope triangular edge side are larger than the other side.

Fig.7 Free surface displacement amplitudes |Wj| of triangle without hole and with hole

Fig.8 Free surface displacement amplitudes |Wj| and hole edge stress τθz＊ at η=0.5

Fig.9 Free surface displacement amplitudes |Wj| and hole edge stress τθz＊ at η=2.0

Overall,the free surface displacement decreases with the increase of the angle of incidence.As the angle of incidence increases,the amount of reflection on the triangular edge decreases.The hole concentrated stress is distributed on both sides of the wave propagation direction,and shows an increasing trend with the increase of the incident angle.This is because the surface boundary effect increases on both sides of the wave propagation direction.

Free surface displacement amplitude and circular hole hoop stress distribution at different incident frequencies,angle and hole depth is shown in Fig.10 and Fig.11.When the incident angleα=0°,the amplitude of the surface displacement decreases as the hole depth increases,while the displacement amplitude of the triangular area increases.This is because when the hole is deeper,more wave energy enters triangle area.When the incident angleα=90°,the change of the hole stress is small.This is because the surface reflection wave that is applied to the hole is smaller.

Fig.10 Free surface displacement amplitudes |Wj| and hole edge stress τθz＊ for various H2 at η=0.5,n1=0.5

Fig.11 Free surface displacement amplitudes |Wj| and the hole edge stress τθz＊ for various H2 at η=2.0,n1=0.5

To better explain the effect of hole offset,only isosceles triangles with symmetrical structures are analyzed.Free surface displacement amplitude and circular hole hoop stress distribution at different incident frequencies,angle and distance of hole offset are shown in Figs.12 and 13.The hole offset has a significant effect on the displacement of the horizontal surface and triangular slope.Even when the horizontal projection of the hole center is located at the triangle hypotenuse bottom corner,the displacement amplitude of the triangular area is still relatively large.

To reveal the influence of dimensionless frequencies on free surface displacement and hole stress,the first row of graphs in Fig.14 and Fig.15 give the displacement amplitudes as a function of 2y/Landηat various angles of incidence (α=0°,45°,90°) and slope (n1=0.5,1.0),and the second column of graphs give the hole stress as a function ofθandη.It shows that the number of wave peaks in the triangular region increases as the wave dimensionless frequencies increases,and the peak and oscillation frequency increase on one side of the incident wave,while it decreases on the other side.At the same time,the peak and oscillation frequency of the area near the larger triangle slope significantly increase,but when the angle of incidence changes is large,the increase moves to the side of the wave incoming direction.The hole concentrated stress is distributed on both sides of the wave propagation direction,and the shear stress near the free boundary is greater than that on the infinite space.This is due to the superposition of the incident wave and the free boundary reflection wave.The displacement amplitude of the triangular surface is the peak shape,while the displacement amplitude of the free flat surface is the mountain shape.The ridge fluctuations are similar to the results obtained in numerical simulations by Shyu(2014).

Fig.12 Free surface displacement amplitudes |Wj| and hole edge stress τθz＊ for various L2 at η=0.5

Fig.13 Free surface displacement amplitudes |Wj| and hole edge stress τθz＊ for various L2 at η=2.0

Fig.14 3D plots of surface displacement amplitudes |Wj| and hole edge stress τθz＊ vs η at n1=1.0

3.3 Time domain response

The transient response is obtained from the frequency domain results through the inverse fast Fourier transform(FFT) algorithm.The incident time signal is a Ricker wavelet

with the characteristic frequencyfc=0.5 Hz.

The calculated frequencies range from 0 to 2.0 Hz with 1/33 Hz intervals.The transfer function for every position is deduced in the previous section for a particular frequencyω(or wave numberk).Then the time domain results can be synthesized by using the inverse FFT,and the shear wave propagates with the velocitycs=3.In Fig.16 the reference point is set to be (x,y)=(8,-16)fort=0 s,in Fig.17 and Fig.18,the reference point is set to be (x,y)=(20,-15) fort=0.The (x,y) reference coordinate system iso6x6y6.

Fig.16 3D plots of surface displacement amplitudes |Wj| vs time at n1=0.5

Fig.17 Snapshots for α=0° at nine specified times

Fig.18 Snapshots for α=90° at nine specified times

Figure 16 shows the synthetic displacement contour with half-space betweeny=-12 to 12 and contains 800 discrete positions located along the surface of the hill.Because the reference point is (8,-16),the vertical Ricker wave reaches the flat surface (x=0) att=2.7 s,after the vertical wave moves away from the flat surface,several scattered waves appear one after another with obviously different amplitudes.The horizontal incident Ricker wave reaches the flat surface position (y=-12)at=1.3 s.When the horizontal wave reaches the hill,several scattered waves appear one after another,also with obviously different amplitudes.

In Figs.17 and 18,snapshots for nodes with equal distance of 0.03 at 0° and 90° incident angles are computed.The snapshots show the wave fields at nine specified times to illustrate the process of the wave propagation and scattering around a triangle shape and shallow circle.In Fig.17,the incident wave passes through the hole,and the existence of the hole produces a scattered wave propagating in the opposite direction(t=6 s).When the incident wave reaches the free surface,a reflected wave is produced,and the triangular area produces a circular scattered wave (t=9 s-11 s);the appearance of displacement also corresponds to Fig.16.In Fig.18,the incident wave passes through the hole,and produces the scattered wave propagating in the opposite direction (t=6 s).The triangular area produces a less obvious circular scattered wave (t=8 s-10 s),when the wave reaches it.The amplitude and range of influence on the left side of the triangle are significantly higher than those on the right side,which also corresponds to Fig.16.In addition,through the time domain results of various points,it can be used for the transient response analysis of underground structures or surface structures to provide support for strength design.

4 Conclusions

The SH wave scattering problem of shallow circular holes under scalene triangles is studied.The wave displacement function,which satisfies the triangle zerostress boundary,is obtained by separation of variables and the wave displacement function which satisfies the half space surface zero-stress boundary is constructed by the complex coordinates and symmetric method.By applying the complex function method and regionmatching technique,a system of equations is established for solving the unknown complex coefficients.Finally,the least squares method is used to solve the undetermined coefficient of the algebraic equations by discrete boundary,and frequency domain and time domain results are obtained.

After numerical simulation,the following conclusions can be drawn:

(1) The validity of the wave function and the least squares method are verified by the continuity of displacement and stress at the auxiliary boundary,the holes free edge are close to zero,and FEM results.

(2) The surface displacement and circular stress are closely related to the position of the circular hole,and the hole center position,which ten times the hole radius depth and three times the hole radius offset,still has significant effects.

(3) The displacement amplitude of the triangular surface is the peak shape,while the displacement amplitude of the free flat surface is the undulating shape,which is similar to the numerical simulation results from Shyu (2014).The number of wave peaks in the triangular region increases as the wave dimensionless frequencies increases,and the peak and oscillation frequency increase on one side of the incident wave,while it decreases on the other side.At the same time,the peak and oscillation frequency of the area near the larger triangle slope significantly increase,but when the angle of incidence becomes larger,the increase moves to the side of the wave incoming direction.The overall surface displacement decreases with increasing incident angle.The hole concentrated stress is distributed on both sides of the wave propagation direction,and the shear stress near the free boundary is greater than that on the infinite space.

(4) The snapshots show the process of the wave propagation and scattering around the triangle and shallow circle in the time domain.In addition,the time domain results of various points can be used for the transient response analysis of underground structures or surface structures to provide support for structural strength design.

(5) Using the auxiliary circle to solve the singularity of the reflex angle at the triangle corner and the theory method in this study is not limited by the range (α1+α2≥60°) of the vertex angle (Linet al.,2010).

Acknowledgment

This work was supported by the National Natural Science Foundation of China (Grant No.12072085),Natural Science Foundation of Heilongjiang Province of China (No.ZD2021A001),the Opening Fund of Acoustics Science and Technology Laboratory (Grant No.SSKF2020011) and the Fundamental Research Funds for the Central Universities (No.3072021CF0206),a program for Innovative Research Teams in the China Earthquake Administration.

Appendix A: Expressions of each angle in the model shown in Fig.1

Appendix B: Expressions of functions

Heres′=s-Hj,`s=s+Hj,Hjis the depth of the corresponding circle center from the surface,and takes a negative value if it is above the horizontal plane.Hstands for the Bessel functionJ(),Hankel function of first kindH1() or Hankel function of second kindH2().