A spectral unit modeling method for seismic wave field simulation suitable for complex terrain

Abstract

The present invention provides a crustal spectral unit optimization modeling method suitable for simulating strong earthquake wave fields in complex terrain conditions, which relates to the field of earthquake engineering. The method comprises determining the terrain parameters of the study area; determining the grid refinement factor based on the wave velocity structure parameters of the study area, and obtaining the minimum distance of the wave velocity isosurface from the bottom; determining the location of the crustal spectral unit refinement transition layer; determining the total number of three-dimensional units before grid refinement based on the location of the crustal spectral unit refinement transition layer; using hexahedral units to mesh the entire area; refining the shallow crust model based on the number of crustal spectral unit layers above the transition layer; calculating the grid quality, counting the number of crustal spectral units to obtain the expected value of the crustal spectral units, and comparing it with the estimated value. The present invention can solve the quality problem of spectral units when modeling crustal structure in areas with complex terrain conditions such as high mountains and canyons, while minimizing the computational degrees of freedom, and achieving stable and low-cost simulation of earthquake motion based on the spectral element method.

Description

Technical Field

[0001] The present invention belongs to the technical field of earthquake engineering, and in particular relates to a seismic wave field simulation spectrum unit modeling method suitable for complex terrain. Background Art

[0002] In structural seismic design, reasonable seismic motion parameters are a prerequisite for ensuring structural safety. For the seismic design of critical structures such as high dams, large reservoirs, and nuclear power facilities, seismic motion time histories are required as input. However, for specific projects, real records meeting similar seismic and geological conditions are often difficult to obtain, necessitating the use of artificial simulation methods to supplement these.

[0003] Currently, commonly used methods for simulating earthquake motion include random methods and physics-based direct simulation methods. Physics-based direct simulation methods (such as the finite difference method, the finite element method, and the spectral element method) can fully consider various influencing factors related to the earthquake source, propagation path, and site, and are naturally reasonable based on the basic physical mechanism of wave propagation. Among them, the spectral element method combines the flexibility of the finite element method with the accuracy and fast convergence of the pseudo-spectral method. It is particularly good at simulating the propagation of earthquake motion in inhomogeneous media and complex terrain conditions.

[0004] Compared to stochastic methods, direct physics-based simulations are often computationally intensive. In spectral element seismic simulations, the maximum simulated frequency is proportional to the wave velocity in the medium and inversely proportional to the distance between element nodes. Therefore, shallow mesh refinement based on the crustal wave velocity structure is often required to maximize computational efficiency. However, in areas with complex terrain, such as high mountains and canyons, sudden changes in steep and complex terrain can degrade mesh quality, leading to non-convergence in the simulation process. Summary of the Invention

[0005] In response to the above-mentioned deficiencies in the prior art, the present invention provides a spectral unit modeling method for seismic wave field simulation suitable for complex terrain. It can solve the quality problem of spectral units when modeling crustal structure in areas with complex terrain conditions such as high mountains and canyons, while minimizing the computational freedom and realizing stable and low-cost simulation of seismic motion based on the spectral element method.

[0006] In order to achieve the above objectives, the present invention adopts a technical solution: a seismic wave field simulation spectrum unit modeling method suitable for complex terrain, comprising the following steps:

[0007] S1. Determine the topographic parameters of the study area;

[0008] S2. Determine the grid refinement factor based on the wave velocity structure parameters of the study area and obtain the minimum distance between the wave velocity isosurface and the bottom;

[0009] S3. Determine the position of the crustal spectrum unit refinement conversion layer based on the determined topographic parameters and the minimum distance between the wave velocity isosurface and the bottom;

[0010] S4. Refine the position of the conversion layer according to the crustal spectrum unit and determine the total number of three-dimensional units before grid refinement;

[0011] S5. Based on the determined topographic parameters of the study area and the total number of three-dimensional elements before mesh refinement, hexahedral elements are used to mesh the entire area.

[0012] S6. Refine the shallow crust model according to the number of crustal spectral units above the transition layer;

[0013] S7. Based on the refinement results, the grid quality is calculated, the number of crustal spectrum units is counted to obtain the expected value of the crustal spectrum unit, and the expected value is compared with the estimated value to complete the optimization modeling of the crustal spectrum unit.

[0014] The beneficial effect of the present invention is that in areas with complex terrain conditions such as high mountains and canyons, steep and complex terrain mutations will deteriorate the quality of the spectral unit grid, resulting in non-convergence of the spectral element method simulation process. The present invention proposes a spectral unit modeling method for seismic wave field simulation suitable for complex terrain that comprehensively considers the spectral unit quality and the computational efficiency of the spectral element method. By calculating the sizes of the upper and lower parts of the conversion layer, and using hexahedral units to perform grid division in the entire domain, and calculating the number of crustal spectral units, the optimized modeling processing of the crustal spectral units is achieved. The present invention can simultaneously solve the spectral unit quality problem and minimize the computational freedom, thereby achieving stable and low-cost simulation of seismic motion using the spectral element method.

[0015] Furthermore, the step S1 includes the following steps:

[0016] S101, determining a minimum simulation altitude according to the maximum fault depth, and determining a simulation length L and a simulation width W according to the relative position of the site and the fault;

[0017] S102, modeling the terrain structure of the study area, and calculating the simulated height H from the highest point of the terrain to the minimum simulated altitude;

[0018] S103. Calculate the elevation difference ΔH between the lowest point and the highest point in the study area to complete parameter determination.

[0019] Furthermore, step S2 includes the following steps:

[0020] S201. Conduct three-dimensional modeling of the crustal wave velocity structure;

[0021] S202, based on the modeling results, by obtaining the surface shear wave velocity V s,0 and regional maximum wave velocity V s,m , calculate the mesh refinement multiple;

[0022] S203, calculate the minimum distance between the velocity isosurface and the bottom according to the grid refinement multiple, where the minimum distance is the height H of the velocity isosurface. r .

[0023] The beneficial effect of the above further solution is: according to the surface shear wave velocity V s,0 and regional maximum wave velocity V s,m , calculating the mesh refinement multiple can ensure that all meshes meet the requirements of the maximum simulation frequency.

[0024] Furthermore, the expression of the wave velocity isosurface is as follows:

[0025] V sr =n r ·V s,0

[0026]

[0027] Among them, V sr represents the wave velocity isosurface, n r Indicates the grid refinement multiple, and ceil() indicates the rounding up function.

[0028] The beneficial effect of the above further solution is that the minimum mesh refinement factor is determined according to the wave velocity, thereby avoiding a reduction in simulation efficiency due to an overly small mesh.

[0029] Furthermore, the expression for refining the transition layer position of the crustal spectrum unit in step S3 is as follows:

[0030]

[0031] Among them, H2 represents the depth of the crustal spectrum unit refinement transition layer, ΔH represents the elevation difference between the lowest and highest points in the study area, and C r and C H All represent constants, n r Indicates the grid refinement multiple, and H indicates the simulated height from the highest point of the terrain to the minimum simulated altitude.

[0032] The beneficial effect of the above further solution is that the present invention can determine the depth of the conversion layer when the number of units is minimized through the above formula.

[0033] Furthermore, step S4 includes the following steps:

[0034] S401. Calculate the vertical dimension D of the upper part of the crustal spectrum unit conversion layer using the following formula: near :

[0035]

[0036] Among them, D near (f m ) indicates that the maximum frequency to be simulated is f m The vertical dimension of the upper part of the crustal spectrum unit conversion layer, f m Indicates the maximum frequency to be simulated, N gll Indicates the number of integration points of the crustal spectrum unit. D0 represents the reference size of the eight-node hexahedral unit when the number of integration points of the crustal spectrum unit is 5, the maximum frequency is 1 Hz, and the surface shear wave velocity is 1000 m / s. V s,0 Represents the surface shear wave velocity at the site, N msh Indicates the number of node segments on the edge of the crustal spectrum unit. For an eight-node hexahedral unit, N msh =1, 27-node hexahedral surface element N msh =2;

[0037] S402, according to D near (f m ), the number of crustal spectral units above the transition layer, n, is calculated using the following formula: z1 :

[0038]

[0039] Among them, H2 represents the depth of the crustal spectrum unit refinement conversion layer;

[0040] S403, according to the number of crustal spectrum unit layers n z1 , the horizontal dimension a1 of the upper part of the crustal spectrum unit transition layer is calculated using the following formula;

[0041]

[0042] Among them, H1 represents the depth from the lowest point of the regional terrain to the transition layer;

[0043] S404. Based on the horizontal dimension a1 of the upper portion of the crustal spectrum unit transition layer, the number of horizontal units in the upper portion of the transition layer is calculated using the following formula: and

[0044]

[0045] Where L represents the simulation length and W represents the simulation width;

[0046] S405, based on the horizontal dimension a1 of the upper portion of the crustal spectrum unit conversion layer, calculate the three-dimensional unit dimension a2 of the lower portion of the conversion layer using the following formula;

[0047] a2=a1n r

[0048] Among them, nr Indicates the mesh refinement multiple;

[0049] S406. Based on the three-dimensional unit size a2, the number of vertical units n at the bottom of the transfer layer is calculated using the following formula: z2 :

[0050]

[0051] Where H represents the simulated height from the highest point of terrain elevation to the minimum simulated elevation;

[0052] S407, according to the number of horizontal units n x1 、n y1 and the number of vertical units n z2 , the total number of three-dimensional elements n before mesh refinement is calculated using the following formulas: x 、n y and n z , complete the determination of the crustal spectrum unit size:

[0053]

[0054] Among them, ceil() represents the rounding up function.

[0055] The beneficial effect of the above further solution is that the present invention can conveniently calculate the total number of elements and estimate the total degrees of freedom of the model by confirming the total number of three-dimensional elements before mesh refinement.

[0056] Furthermore, the expression of the expected value of the crustal spectrum unit in step S7 is as follows:

[0057]

[0058] Where n represents the expected value of the crustal spectrum unit, C represents a constant, H2 represents the depth of the crustal spectrum unit refinement transition layer, ΔH represents the elevation difference between the lowest and highest points in the study area, and n r Indicates the grid refinement multiple, H indicates the simulated height from the highest point of terrain elevation to the minimum simulated elevation, L indicates the simulated length, W indicates the simulated width, and D near (f m ) indicates that the maximum frequency to be simulated is f m Vertical dimensions of the upper part of the transition layer of the crustal spectrum unit.

[0059] The above further solution has the beneficial effect of conveniently calculating the total number of elements for comparison with other projects and estimating the total degrees of freedom of the model. Based on the total number of degrees of freedom and the comparison with other projects, it is possible to allocate reasonable computing resources for earthquake simulation. BRIEF DESCRIPTION OF THE DRAWINGS

[0060] Figure 1 Flow chart of the method of the present invention.

[0061] Figure 2 It is a schematic diagram of the terrain profile in the present invention.

[0062] Figure 3 This is a schematic diagram of the decomposition of the global spectrum unit in the present invention.

[0063] Figure 4 This is a schematic diagram of the refinement of shallow crustal units in the present invention.

[0064] Figure 5 Schematic diagram of the overall skewness of the spectrum unit in this embodiment.

[0065] Figure 6 Schematic diagram of three-component earthquake motion obtained by simulating the crustal spectrum unit model using the technical solution of the present invention in this embodiment. DETAILED DESCRIPTION

[0066] The specific embodiments of the present invention are described below to facilitate understanding of the present invention by those skilled in the art. However, it should be clear that the present invention is not limited to the scope of the specific embodiments. For those skilled in the art, as long as various changes are within the spirit and scope of the present invention as defined and determined by the appended claims, these changes are obvious, and all inventions and creations utilizing the concepts of the present invention are protected.

[0067] Example

[0068] like Figure 1 As shown, the present invention provides a seismic wave field simulation spectrum unit modeling method suitable for complex terrain, and its implementation method is as follows:

[0069] S1. Determine the terrain parameters of the study area. The implementation method is as follows:

[0070] S101, determining a minimum simulation altitude according to the maximum fault depth, and determining a simulation length L and a simulation width W according to the relative position of the site and the fault;

[0071] S102, modeling the terrain structure of the study area, and calculating the simulated height H from the highest point of the terrain to the minimum simulated altitude;

[0072] S103. Calculate the elevation difference ΔH between the lowest point and the highest point in the study area to complete parameter determination.

[0073] S2. Determine the grid refinement factor based on the wave velocity structure parameters of the study area and obtain the minimum distance between the wave velocity isosurface and the bottom. The implementation method is as follows:

[0074] S201. Conduct three-dimensional modeling of the crustal wave velocity structure;

[0075] S202, based on the modeling results, by obtaining the surface shear wave velocity V s,0 and regional maximum wave velocity V s,m , calculate the mesh refinement factor:

[0076]

[0077] Among them, ceil() represents the upward rounding function;

[0078] S203, calculate the minimum distance between the velocity isosurface and the bottom according to the grid refinement multiple, where the minimum distance is the height H of the velocity isosurface. r ,like Figure 2 As shown;

[0079] V sr =n r ·V s,0 (2)

[0080] Among them, V sr represents the wave velocity isosurface, n r Indicates the mesh refinement multiple;

[0081] S3. Determine the position of the crustal spectrum unit refinement conversion layer based on the determined topographic parameters and the minimum distance between the wave velocity isosurface and the bottom;

[0082] In this embodiment, the key to the present invention is to r ) to create a horizontal transition zone within this region, making the spectral units relatively flat. If H2 is too small, thin, flake-like units will be generated near the lowest point of the terrain, and their unit quality will also be poor. If H2 is too large, the simulation efficiency will be reduced because the unit density in the upper part is much greater than that in the lower part. Therefore, choosing the depth of H2 to optimize the overall quality of the crustal spectral unit model while minimizing the simulation scale is the core of the present invention.

[0083] In this embodiment, the depth H2 of the conversion layer can be obtained by taking the derivative of the total number of upper and lower units n of the conversion layer with respect to H2. Minimizing the simulation scale is equivalent to solving the following equation:

[0084]

[0085] Among them, C r and C H Both represent constants, and ΔH represents the elevation difference between the lowest and highest points in the study area.

[0086]

[0087] Among them, n rIndicates the grid refinement multiple, and H indicates the simulated height from the highest point of the terrain to the minimum simulated altitude.

[0088] The following conditions must be met to solve equation (3):

[0089]

[0090] When n r When ≥2, the above formula always holds. From formula (3), we can get that when H2 takes the following value, the total number of units n has the minimum value:

[0091]

[0092] S4. Determine the total number of three-dimensional cells before grid refinement based on the position of the conversion layer of the crustal spectrum unit refinement. The implementation method is as follows:

[0093] S401. Calculate the vertical dimension D of the upper part of the crustal spectrum unit conversion layer using the following formula: near :

[0094]

[0095] Among them, D near (f m ) indicates that the maximum frequency to be simulated is f m The vertical dimension of the upper part of the crustal spectrum unit conversion layer, f m Indicates the maximum frequency to be simulated, N gll Indicates the number of integration points of the crustal spectrum unit. D0 represents the reference size of the eight-node hexahedral unit when the number of integration points of the crustal spectrum unit is 5, the maximum frequency is 1 Hz, and the surface shear wave velocity is 1000 m / s. V s,0 Represents the surface shear wave velocity at the site, N msh Indicates the number of node segments on the edge of the crustal spectrum unit. For an eight-node hexahedral unit, N msh =1, twenty-seven-node hexahedral surface element N msh =2;

[0096] S402, according to D near (f m ), the number of crustal spectral units above the transition layer, n, is calculated using the following formula: z1 :

[0097]

[0098] Among them, H2 represents the depth of the crustal spectrum unit refinement conversion layer;

[0099] S403, according to the number of crustal spectrum unit layers n z1The horizontal dimension a1 of the upper portion of the crustal spectrum unit conversion layer is calculated using the following formula. In this embodiment, the horizontal dimension a1 on the unit conversion layer is the same as the vertical dimension of the unit at the average terrain elevation, that is:

[0100]

[0101] Among them, H1 represents the depth from the lowest point of the regional terrain to the transition layer;

[0102] S404. Based on the horizontal dimension a1 of the upper portion of the crustal spectrum unit transition layer, the number of horizontal units in the upper portion of the transition layer is calculated using the following formula: and

[0103]

[0104] Where L represents the simulation length and W represents the simulation width;

[0105] S405, based on the horizontal dimension a1 of the upper portion of the crustal spectrum unit conversion layer, calculate the three-dimensional unit dimension a2 of the lower portion of the conversion layer using the following formula;

[0106] a2=a1n r (13)

[0107] Among them, n r Indicates the mesh refinement multiple;

[0108] S406. Based on the three-dimensional unit size a2, the number of vertical units n at the bottom of the transfer layer is calculated using the following formula: z2 :

[0109]

[0110] Where H represents the simulated height from the highest point of terrain elevation to the minimum simulated elevation;

[0111] S407, according to the number of horizontal units n x1 、n y1 and the number of vertical units n z2 , the total number of three-dimensional elements n before mesh refinement is calculated using the following formulas: x 、n y and n z , complete the determination of the crustal spectrum unit size:

[0112]

[0113] Among them, ceil() represents the rounding up function.

[0114] S5. Based on the determined topographic parameters of the study area and the total number of three-dimensional elements before mesh refinement, hexahedral elements are used to mesh the entire area.

[0115] In this embodiment, according to the determined terrain parameters and the determined total number of three-dimensional units, hexahedral units are used to perform mesh generation in the entire domain. The results are as follows: Figure 3 shown.

[0116] S6. Refine the shallow crust model according to the number of crustal spectral units above the transition layer;

[0117] In this embodiment, the shallow crust model is refined according to the determined number of upper unit layers of the transition layer, and the results are as follows: Figure 4 shown.

[0118] S7. Based on the refinement results, calculate the grid quality, count the number of crustal spectrum units to obtain the expected value of the crustal spectrum unit, and compare it with the estimated value to complete the optimization modeling of the crustal spectrum unit. The expression of the expected value of the crustal spectrum unit is as follows:

[0119]

[0120] Where n represents the expected value of the crustal spectrum unit, C represents a constant, H2 represents the depth of the crustal spectrum unit refinement transition layer, ΔH represents the elevation difference between the lowest and highest points in the study area, and n r Indicates the grid refinement multiple, H indicates the simulated height from the highest point of terrain elevation to the minimum simulated elevation, L indicates the simulated length, W indicates the simulated width, and D near (f m ) indicates that the maximum frequency to be simulated is f m Vertical dimensions of the upper part of the transition layer of the crustal spectrum unit.

[0121] In this embodiment, Figure 5 As shown in the figure, the skewness distribution of the spectral unit using the technical solution of the present invention under complex terrain conditions is given. It can be seen that the maximum value of the unit skewness under complex terrain conditions is 0.708, which is less than the maximum limit of the spectral element method software specfem3D. That is, the model can directly simulate seismic motion without any additional processing, ensuring the stability and efficiency of the simulation.

[0122] In this embodiment, the spectral element model obtained by the technical solution of the present invention has a relatively small number of elements. The total number of spectral elements in the technical solution of the present invention is 899,232. The existing solution can correctly perform spectral element method simulation with 50% more total elements than the technical solution of the present invention. The technical solution of the present invention significantly improves simulation efficiency and saves about half of the simulation cost.

[0123] In this embodiment, the model obtained by the technical solution of the present invention can be directly used for spectral element method simulation, such as Figure 6 As shown, Figure 6 The three-component seismic motion of the site obtained by spectral element method simulation using the crustal spectral unit model of the technical solution of the present invention is given. Its maximum effective frequency meets the set target value, and there is no non-convergence in the simulation process.

Claims

1. A spectral unit modeling method for seismic wave field simulation suitable for complex terrain, characterized in that: The following steps are involved: S1. Determine the topographic parameters of the study area; S2. Determine the grid refinement factor based on the wave velocity structure parameters of the study area and obtain the minimum distance between the wave velocity isosurface and the bottom; S3. Determine the position of the crustal spectrum unit refinement conversion layer based on the determined topographic parameters and the minimum distance between the wave velocity isosurface and the bottom; S4. According to the position of the conversion layer of the crustal spectrum unit refinement, the total number of three-dimensional units before grid refinement is determined, which is specifically: According to the horizontal size of the upper part of the crustal spectrum unit transition layer , the three-dimensional unit size of the lower part of the conversion layer is calculated using the following formula: ; in, Indicates the mesh refinement multiple; According to the three-dimensional element size , the number of vertical units below the transfer layer is calculated using the following formula: : in, Indicates the simulated height from the highest point of terrain elevation to the minimum simulated elevation. Indicates the depth of the crustal spectrum unit refinement conversion layer; According to the number of horizontal units 、 and the number of vertical units , the total number of three-dimensional elements before mesh refinement is calculated using the following formulas: 、 and , complete the determination of the crustal spectrum unit size: in, represents the ceiling function, Indicates the number of crustal spectral unit layers; S5. Based on the determined topographic parameters of the study area and the total number of three-dimensional elements before mesh refinement, hexahedral elements are used to mesh the entire area. S6. Refine the shallow crust model according to the number of crustal spectral units above the transition layer; S7. Based on the refinement results, the grid quality is calculated, the number of crustal spectrum units is counted to obtain the expected value of the crustal spectrum unit, and the expected value is compared with the estimated value to complete the optimization modeling of the crustal spectrum unit.

2. The seismic wave field simulation spectrum unit modeling method applicable to complex terrain according to claim 1, characterized in that: The step S1 comprises the following steps: S101. Determine the minimum simulated altitude based on the maximum fault depth, and determine the simulation length based on the relative position of the site and the fault. and simulated width ; S102. Model the terrain structure of the study area and calculate the simulated altitude from the highest point of the terrain to the minimum simulated altitude. ; S103. Calculate the elevation difference between the lowest and highest points in the study area , complete the parameter determination.

3. The seismic wave field simulation spectrum unit modeling method applicable to complex terrain according to claim 1, characterized in that: The step S2 comprises the following steps: S201. Conduct three-dimensional modeling of the crustal wave velocity structure; S202, based on the modeling results, by obtaining the surface shear wave velocity of the site and regional maximum wave speed , calculate the mesh refinement multiple; S203, according to the grid refinement multiple, calculate the minimum distance between the velocity isosurface and the bottom, where the minimum distance is the height of the constant velocity surface .

4. The seismic wave field simulation spectrum unit modeling method applicable to complex terrain according to claim 3, characterized in that: The expression of the wave velocity isosurface is as follows: in, represents the wave velocity isosurface, represents the mesh refinement factor, Represents the ceiling function.

5. The seismic wave field simulation spectrum unit modeling method applicable to complex terrain according to claim 1, characterized in that: The expression for refining the conversion layer position of the crustal spectrum unit in step S3 is as follows: in, Indicates the depth of the crustal spectrum unit refinement conversion layer, It represents the elevation difference between the lowest and highest points in the study area. and All represent constants, represents the mesh refinement factor, Indicates the simulated altitude from the highest point of the terrain to the minimum simulated altitude.

6. The method for modeling a seismic wave field simulation spectrum unit suitable for complex terrain according to claim 1, characterized in that: The step S4 comprises the following steps: S401. Calculate the vertical dimensions of the upper part of the crustal spectrum unit transition layer using the following formula: : in, The maximum frequency to be simulated is The vertical size of the upper part of the transition layer of the crustal spectrum unit, represents the maximum frequency to be simulated, represents the number of integration points of the crustal spectrum unit, The reference size of an eight-node hexahedral element when the number of integration points of the crustal spectrum element is 5, the maximum frequency is 1 Hz, and the surface shear wave velocity is 1000 m / s. represents the surface shear wave velocity at the site, Indicates the number of node segments on the edge of the crustal spectrum unit. For eight-node hexahedral units , 27-node hexahedral surface element ; S402, according to , the number of crustal spectral units above the transition layer is calculated using the following formula: : ； S403, according to the number of crustal spectrum units The horizontal dimension of the upper part of the crustal spectrum unit transition layer is calculated using the following formula: ; in, Indicates the depth from the lowest point of the regional terrain to the transition layer; S404, according to the horizontal dimension of the upper part of the crustal spectrum unit conversion layer , the number of horizontal units above the conversion layer is calculated using the following formula: and : in, Indicates the simulation length, Indicates the simulated width.

7. The method for modeling a seismic wave field simulation spectrum unit suitable for complex terrain according to claim 1, characterized in that: The expression of the expected value of the crustal spectrum unit in step S7 is as follows: in, represents the expected value of the crustal spectrum unit, represents a constant, Indicates the depth of the crustal spectrum unit refinement conversion layer, It represents the elevation difference between the lowest and highest points in the study area. represents the mesh refinement factor, Indicates the simulated height from the highest point of terrain elevation to the minimum simulated elevation. Indicates the simulation length, Indicates the simulated width, The maximum frequency to be simulated is Vertical dimensions of the upper part of the transition layer of the crustal spectrum unit.

Images