CN117763914A - Numerical analysis method, terminal and medium for stress field and earthquake dynamic response of near-fault tunnel - Google Patents

Numerical analysis method, terminal and medium for stress field and earthquake dynamic response of near-fault tunnel Download PDF

Info

Publication number
CN117763914A
CN117763914A CN202311813634.4A CN202311813634A CN117763914A CN 117763914 A CN117763914 A CN 117763914A CN 202311813634 A CN202311813634 A CN 202311813634A CN 117763914 A CN117763914 A CN 117763914A
Authority
CN
China
Prior art keywords
model
stress
tunnel
structural surface
numerical
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202311813634.4A
Other languages
Chinese (zh)
Other versions
CN117763914B (en
Inventor
吴苏
石崇
杨再宏
张聪
李汪洋
乔天蛟
杨泽峰
姚震
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Hohai University HHU
Original Assignee
Hohai University HHU
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Hohai University HHU filed Critical Hohai University HHU
Priority to CN202311813634.4A priority Critical patent/CN117763914B/en
Publication of CN117763914A publication Critical patent/CN117763914A/en
Application granted granted Critical
Publication of CN117763914B publication Critical patent/CN117763914B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Testing Of Devices, Machine Parts, Or Other Structures Thereof (AREA)

Abstract

The invention relates to a numerical analysis method, a terminal and a medium of a tunnel stress field and an earthquake dynamic response of a near fault, which belong to the technical field of earthquake analysis, and are characterized in that based on a finite difference numerical platform, numerical model parameters and indoor test parameters are matched and calibrated, a near fault tunnel model and a fault contact surface model are established, the secondary stress field of the near fault tunnel model after tunnel excavation is simulated, and stress concentration coefficients are calculated; synthesizing and designing a seismic acceleration time course according to the seismic acceleration response spectrum, and filtering and correcting a base line to obtain an artificial seismic acceleration response time course meeting calculation requirements; and according to the artificial seismic wave generated by the designed response spectrum as a dynamic boundary condition, calculating a stress field of surrounding rock of the near-fault tunnel under the seismic action and a seismic dynamic response. According to the invention, the mechanical behavior of the surrounding rock under the influence of faults of the tunnel is fully considered, and the stability of the surrounding rock of the tunnel and the efficiency of structural earthquake-resistant design are improved through analysis of stress fields and dynamic responses.

Description

Numerical analysis method, terminal and medium for stress field and earthquake dynamic response of near-fault tunnel
Technical Field
The invention relates to the technical field of earthquake analysis, in particular to a numerical analysis method, a terminal and a medium for stress field and earthquake dynamic response of a tunnel near fault.
Background
The fault is a common bad geological structure in the process of tunnel excavation and use, and the mechanical behavior and stability of surrounding rocks of the near-fault tunnel are obviously influenced by a structural stress field. The construction of tunnel engineering should be selected in more stable stratum at first, should first choose to dodge and dodge the distance and need accord with the requirement of design rule if meet the fault, perhaps make tunnel and fault's trend mutually perpendicular. However, in the actual tunnel engineering construction, due to the constraints of geological conditions, design schemes, construction methods and other factors, part of tunnel engineering is inevitably constructed at a place close to a fault, and is oblique or parallel to the fault trend, so that great challenges are caused to the control of surrounding rock stability.
Current approaches to investigate such problems have involved simulation of near-fault tunnel engineering by preparing similar materials, such as discontinuous faults. The theoretical basis of the model test is that a model meeting the similar condition is constructed according to the similar principle, namely, the model is required to be similar to an entity (prototype), and the model can reflect the condition of the entity. At this time, the model and the prototype, in addition to being geometrically similar, must also be proportionally similar to the same kind of physical quantities (e.g. stress, strain, displacement). These proportions must satisfy various mechanical conditions. However, due to the selection of materials and the manufacture of the model, the period of the model test is longer, the workload is large, the similar conditions are difficult to completely meet, and the seismic dynamic boundary conditions are difficult to be applied to the similar model according to the designed response spectrum.
Therefore, the model of the fault contact surface is required to be adopted for simulation calculation, a computer numerical platform is used for replacing a similar model in a real laboratory, the preparation period of the model test is shortened, and the workload of the model test is saved.
Disclosure of Invention
The invention aims to overcome the defects in the prior art, provides a numerical analysis method for stress field and earthquake dynamic response of a tunnel near fault, establishes a geomechanical model, simulates dynamic boundary conditions of the geomechanical model, utilizes the high efficiency of a numerical simulation platform and the adaptability to different projects, is easy to operate and realize, and is beneficial to improving the efficiency of earthquake-resistant design of a tunnel structure.
In order to achieve the above purpose, the present invention provides the following technical solutions:
in a first aspect, the invention provides a method for analyzing the stress field and the seismic dynamic response of a near-fault tunnel, which is characterized by comprising the following steps:
s1: based on a Flac3D6.0 platform, a numerical model of a complete rock sample without a structural surface and a numerical model of a rock sample with the structural surface are pre-constructed, the numerical model of the rock sample with the structural surface comprises the structural surface, deformation parameters of the structural surface are obtained by comparing the numerical model of the complete rock sample without the structural surface with the numerical model of the rock sample with the structural surface, a speed boundary condition is applied according to the deformation parameters of the structural surface, normal rigidity and shear rigidity parameters of the structural surface are calculated, a fault contact surface model is obtained by a Flac3D software model library or a finite difference method, and the normal rigidity and the shear rigidity of the structural surface are calibrated for the fault contact surface model;
S2: constructing a geomechanical model according to the normal stiffness of the structural surface and the shear stiffness of the structural surface, wherein the geomechanical model comprises a plurality of units, and each unit comprises a plurality of nodes and a plurality of surface units; assuming that surrounding rock materials in the geomechanical model obey a line elastic strength criterion, and that a fault contact surface model obeys a Moire-Coulomb strength criterion; giving deformation parameters of the structural surface to the geomechanical model, applying normal constraint conditions and ground stress boundary conditions according to boundary coordinates of the geomechanical model, performing simulation calculation on the geomechanical model, setting displacement and speed of nodes of the geomechanical model to be zero, and reserving calculation results of an initial ground stress field;
s3: s2, constructing a near-fault tunnel model and excavating a tunnel based on the calculation result of the initial ground stress field in the step S; simulating a secondary stress field of the near-fault tunnel model after tunnel excavation and calculating a stress concentration coefficient;
s4: based on a Matlab platform, according to the related requirements on tunnel earthquake-resistant design in highway tunnel earthquake-resistant design specification (JTG 2232-2019) and underground structure earthquake-resistant design standard (GB/T51336-2018), artificial earthquake waves meeting the requirements of the design specification are obtained by synthesizing and designing artificial earthquake acceleration time course according to an earthquake acceleration response spectrum and performing filtering and baseline correction;
S5: and (3) applying a viscoelastic artificial boundary and a damping coefficient to the geomechanical model based on the calculation result of the stress concentration coefficient of the secondary stress field in the step (S3), applying the artificial seismic wave to the geomechanical model boundary as a seismic power boundary condition according to the relative positions of the seismic source and the tunnel engineering, recording the peak acceleration of the unit and the peak stress concentration coefficient of the surrounding rock in the reaction process, and analyzing the dynamic response of the geomechanical model to obtain a stress field diagram of the distribution of the peak stress concentration coefficient of the near-fault tunnel and the horizontal acceleration reaction time course of each monitoring point.
Further, the construction of the rock sample numerical model in step S1 includes:
s1.1: generating a numerical model of a complete rock sample without a structural surface and a numerical model of a rock sample with a structural surface by adopting software, and setting the ratio of the height to the diameter of the numerical model of the complete rock sample without the structural surface and the numerical model of the rock sample with the structural surface to be 2 according to a uniaxial compressive strength test standard recommended by an engineering rock mass test method standard (GB/T50266-2013): 0;
s1.2: providing the rock sample material obeys the Moire-Coulomb strength criterion, comparing a numerical model of a complete rock sample without a structural surface with a numerical model of a rock sample with a structural surface to obtain deformation parameters of the structural surface, wherein the normal deformation of the numerical model of the complete rock sample without the structural surface is always smaller than that of the numerical model of the rock sample with the structural surface, the slope of a stress-deformation curve in the numerical model of the complete rock sample without the structural surface is larger than that of the numerical model of the rock sample with the structural surface, and applying normal constraint conditions to the boundary of the numerical model of the rock sample with the structural surface according to the boundary coordinates of the numerical model of the rock sample with the structural surface;
S1.3: simulating the loading rate of a uniaxial compression test, and respectively applying speed boundary conditions to a complete rock sample without a structural surface and a rock sample with the structural surface; the sample takes a coordinate axis z as a rotation axis, and a circular plane is positioned on an xoy coordinate plane; generating an upper boundary and a lower boundary, applying a rate of 10 to the upper boundary and the lower boundary -5 Speed boundary conditions of m/s; the method adopts the fish language to realize the uniaxial compression test, and comprises the following steps:
s1.3.1: generating a global mapping variable surface;
s1.3.2: the numerical model of the complete rock sample without the structural surface and the numerical model of the rock sample with the structural surface comprise node bodies, the node bodies are traversed, and if a node body gp is positioned at the upper boundary, the node bodies are saved to a surface;
s1.3.3: applying (1-10) x 10 to global map variable surface 5 m/s speed boundary condition, simulating the loading rate of a uniaxial compression test;
s1.3.4: traversing the node body, recording the unbalanced force in the z direction on the global mapping variable surface by using a function gp.force.unbal.z (gp) of the maximum unbalanced force, and storing the unbalanced force in the z direction to a variable zforce;
s1.3.5: calculate normal stress axiamaress = zdoece/pi r 2 ；
Wherein AxialStress represents normal stress; zforce represents the z-direction imbalance force; r represents the radius of the numerical model of the complete rock sample without structural face and the numerical model of the rock sample with structural face; pi represents the circumference ratio;
s1.3.6: setting artificial convergence conditions, and calculating convergence by the model when the displacement of the numerical model of the complete rock sample without the structural surface and the numerical model of the rock sample with the structural surface reaches the termination strain;
let the deformation of the complete rock sample be DeltaV r The deformation of the rock sample containing the structural plane was DeltaV t Normal closing deformation DeltaV of structural surface j The method comprises the following steps:
ΔV j ＝ΔV t -ΔV r (1)
s1.4: fitting σ using the equation n -ΔV j Stress-deformation relationship curve:
v in m Representing the maximum normal closing deformation of the structural surface, sigma n -ΔV j An asymptote of the stress-deformation relationship curve; sigma (sigma) n Representing normal stress;representing the normal stiffness of the structural plane fitted by the equation;
s1.5: solving the tangential slope K of the asymptote equation at the origin according to the equation fitted by equation (2) ni I.e. the initial normal stiffness of the structural face:
s1.6: calibrating the normal stress sigma according to the rock sample parameters containing the structural surface n Under the action, the normal rigidity of the structural surface and the shear rigidity of the structural surface of the fault contact surface model:
Wherein L is the trend length of the fault; JRC is the roughness coefficient of the fault structural plane; JCS is the wall rock strength of the structural face; phi (phi) r Residual friction angle of the fault structural surface; k (K) n Is the normal rigidity of the structural surface; k (K) s Is the shear stiffness of the structural face.
Further, in step S2, the construction of the tunnel model and the fault interface model includes:
s2.1: constructing a tunnel model, namely constructing the model by adopting a structured grid for a tunnel section with a simpler shape, such as a circular or semicircular arch tunnel; for tunnel sections with relatively complex shapes, such as straight wall three-center arches and horseshoe-shaped tunnels, adopting unstructured grids to build a model;
establishing a xoz plane model of the tunnel section according to the real size, and enabling the geometric center of the model to be located at the origin of coordinates; setting the thickness of the model to be 1.0;
s2.2: constructing a surrounding rock model, setting surrounding rocks with the range of 4-6 times of tunnel radius around the built tunnel model, and respectively grouping and naming an upper disc and a lower disc of the surrounding rocks;
s2.3: constructing a fault contact surface model, and building the contact surface model according to the occurrence and the property of faults based on a finite difference numerical platform; for the structured grid, creating a surface unit between the upper disc and the lower disc of the fault to generate a contact surface model; for unstructured grids, directly generating a contact surface model between the upper disc grouping and the lower disc grouping; giving K to the contact surface model n 、K s Mechanical parameters of the iso-fault.
And S2.4, applying a normal constraint condition and a ground stress boundary condition according to the boundary coordinates of the geomechanical model, setting the normal speed of a model boundary surface unit to be 0, and applying horizontal stress and vertical stress to the model.
Further, in the step S3, simulation of surrounding rock stress fields after tunnel excavation is realized by adopting a fish language:
s3.1: traversing the model unit, extracting the stress state of each unit zp, and realizing by using a fish language:
middle sigma x The horizontal stress of the unit body in the numerical calculation model is calculated; sigma (sigma) z The vertical stress of the unit body in the numerical calculation model is calculated; τ xz The shear stress of the unit body in the numerical calculation model is calculated; zp is the pointer variable of each unit;
s3.2: calculating the distance R from the center of each unit body to the origin of coordinates/the geometric center point of the tunnel:
s3.3: taking a straight wall semicircular arch as an example without losing generality, calculating sine and cosine of a connecting line from the center of each unit body to the origin of coordinates according to the shape of the section of the tunnel, and realizing by using a fish language:
s3.3.1: the unit body is positioned at the arc part, namely zone. Pos.z (zp) >0.0
S3.3.2: the unit body is positioned at the straight wall part, namely-1.2 < zone.pos.z (zp) <0.0
S3.3.3: the unit body is positioned at the corner part, namely zone.pos.z (zp) < -R, and math.abs (zone.pos.x (zp)) > R
S3.3.4: the unit body is located at the bottom plate position, i.e., zone.pos.z (zp) < -R, and math.abs (zone.pos.x (zp)) < R
S3.4: calculating radial stress sigma of surrounding rock around tunnel according to rotating shaft formula r And hoop stress sigma θ ：
S3.5: defining stress concentration coefficient K as hoop stress sigma θ And vertical ground stress p v Solving the ratio of K of the surrounding rock unit zp:
K＝σ θ /p v (13)
the calculation result is saved to the variable zone. Extra (zp).
Further, step S4 adopts a triangle progression method to synthesize an artificial earthquake motion time course, and is realized by using a Matlab platform:
s4.1: creating a main program for synthesizing artificial earthquake motion reaction time course acceleration, setting variables for saving reaction duration time T, reaction interval dT and frequency spectrum upper limit omega u Lower spectral limit omega l Frequency domain sampling interval Δω; generating an empty array variable acc for storing the reaction acceleration;
s4.2: according to the random vibration theory, inputting the maximum reaction distribution of a simple substance point system in a stable random process:
in which a is g (t) is a stationary random process, initial phaseIs a random phase angle, omega, uniformly distributed in (0, 2 pi) k A is the circular frequency of the kth frequency, A k (ω k ) For the amplitude of the kth frequency, it is calculated by:
s (omega) k ) As a function of power spectral density, Δω is the frequency sampling interval, calculated by:
Δω＝2π/T (16)
Wherein T is the duration of the reaction;
s4.3: creating an intensity envelope function subprogram, calling in a main program for synthesizing artificial earthquake motion reaction time course acceleration, inputting a time variable t into the subprogram, and returning an intensity envelope function value f (t); calculating the starting time t of acceleration plateau in the intensity envelope function according to the equivalent shock level M and the equivalent shock center distance D 1 Acceleration plateau end time t 2 And attenuation coefficient c:
the intensity envelope function is constructed as follows:
s4.4: creating a power density function subroutine, and calling in a main program for synthesizing artificial earthquake motion reaction time course acceleration; introducing a frequency variable omega into a subroutine k Return power spectral density function value S (omega) k ) The method comprises the steps of carrying out a first treatment on the surface of the Setting variables for saving the self-vibration period T and the field characteristic period T of the structure g Damping ratio ζ, acceleration response spectrum maximum S max The decay index gamma of the descending segment and the reaction overrun probability P; reaction spectrum S for synthesizing design acceleration time course a The method comprises the following steps:
the power spectral density function is calculated by:
s4.5: and (3) operating a main program for synthesizing the acceleration of the artificial earthquake motion reaction time course, calling the subprogram, synthesizing the artificial earthquake motion acceleration time course, and storing the artificial earthquake motion acceleration time course into an array variable acc:
s4.6: creating a low-pass filtering main program, and generating an array variable input_acc to store the artificial earthquake motion acceleration time course acc in the step S5.5; generating variable preserving filter order n and cut-off frequency omega n ；
Using the Butterworth low-pass filter function of Matlab platform, generating transfer function coefficients a and b of the filter, returning a row vector with length (n+1):
filtering through transfer function coefficients, and returning the acceleration reaction time interval filtered_acc after filtering;
s4.7: creating a base line correction main program, and generating an array variable iput_acc to store the acceleration response time interval filtered_acc filtered in the step S5.6; calculating running average acceleration using a sampling sliding window algorithm: setting the size k of a sliding window, and returning to an array of local k-point average values, wherein each average value is calculated on adjacent elements of the iput_acc through the sliding window with the length k; when k is even, the window is centered on the current and previous element; when there are not enough elements to fill the window, the window size will automatically truncate at the endpoint; when the window is truncated, only averaging the elements filling the window; returning to the running average acceleration smooth_a, the acceleration response time course after baseline correction is:
modified_a＝iput_acc-smooth ed_a (23)
wherein modified_a is the acceleration response time course; input_acc is a set of variables; the smoothed_a is the return running average acceleration;
s4.8: creating a time-course conversion subprogram which can be called in any main program; converting the acceleration reaction time course into a speed reaction time course and a displacement reaction time course; the acceleration response time course is transmitted to the subprogram and stored as an array variable iput_acc, and a speed variable vel and a displacement variable disp are initialized:
Calculating a velocity reaction time course and a displacement reaction time course by adopting numerical integration as return values:
s4.9: creating a Fourier transform subprogram which can be called in any main program; converting the time domain earthquake motion signal into a frequency domain signal by adopting a fast Fourier transform function of a matrix calculation numerical platform for analysis; and (4) transmitting an acceleration response time course to the subroutine, saving the acceleration response time course as an array variable iput_acc, and returning an acceleration amplitude frequency spectrum.
In a second aspect, the present invention provides an electronic terminal, including a processor and a storage medium; the storage medium is used for storing instructions; the processor being operative according to the instructions to perform the steps of the method according to any one of claims 1 to 6.
In a third aspect, the invention provides a computer-readable storage medium on which a computer program is stored, characterized in that the program, when being executed by a processor, carries out the steps of the method according to any one of claims 1 to 6.
Compared with the prior art, the invention has the beneficial effects that:
(1) The method is based on a finite difference numerical simulation platform, deformation parameters of a structural surface are calculated, fault contact surface model parameters are calibrated, a geomechanical model is built, a near-fault tunnel model is built on the basis of the calculation result of the geomechanical model, a secondary stress field of the near-fault tunnel model after the tunnel is simulated, and a stress concentration coefficient is calculated, so that structural stress boundary conditions formed by faults can be fully considered according to the occurrence and properties of faults in actual engineering, and stability control of surrounding rocks of the near-fault tunnel is facilitated.
(2) According to the method, the artificial earthquake motion acceleration time course can be synthetically designed according to the parameters of the standard response spectrum, and the dynamic boundary conditions of the geomechanical model are simulated through the signal processing modes such as filtering, baseline adjustment, fourier transformation and the like, so that the artificial earthquake motion acceleration time course meeting the design requirements is generated, and the efficiency of the earthquake-resistant design of the tunnel structure is improved.
(3) The invention fully plays the high efficiency of the numerical simulation platform and the adaptability to different projects, calculates the stress field and the dynamic response of the surrounding rock of the near-fault tunnel under the earthquake action according to the artificial earthquake wave generated by the design response spectrum as the dynamic boundary condition, and has the advantages of easy operation and realization of the process method, specific and visible calculation result and easy understanding and judgment.
Drawings
The accompanying drawings, which are incorporated in and constitute a part of this specification, illustrate embodiments of the invention and together with the description, serve to explain the principles of the invention.
The invention may be more clearly understood from the following detailed description taken in conjunction with the accompanying drawings, in which:
FIG. 1 is a flow chart of a method for numerical analysis of near-fault tunnel stress field and seismic dynamic response provided by an embodiment of the invention;
FIG. 2 is a graph of normal deformation of a rock sample for a near-fault tunnel stress field and a numerical analysis method of seismic dynamic response provided by an embodiment of the present invention;
FIG. 3 is a structural plane closure deformation graph of a near-fault tunnel stress field and seismic dynamic response numerical analysis method provided by an embodiment of the invention;
FIG. 4 is a schematic diagram of a near-fault tunnel stress field calculation for a near-fault tunnel stress field and a numerical analysis method of seismic dynamic response provided by an embodiment of the invention;
FIG. 5 is a graph of a near-fault tunnel numerical model for analyzing the stress field and dynamic response of the near-fault tunnel surrounding rock under the action of an earthquake, provided by the embodiment of the invention;
FIG. 6 is a graph of a near-fault tunnel surrounding rock stress field and a dynamic response free fault tunnel surrounding rock stress field under analysis of seismic action provided by an embodiment of the invention;
FIG. 7 is a graph of near-fault tunnel surrounding rock stress field for analyzing near-fault tunnel surrounding rock stress field and dynamic response under the action of an earthquake provided by an embodiment of the invention;
FIG. 8 is a designed horizontal acceleration response spectrum for analyzing the stress field and dynamic response of near-fault tunnel surrounding rock under the action of earthquake, provided by the embodiment of the invention;
FIG. 9 is a graph of artificially synthesized earthquake motion acceleration response time for analyzing stress fields and dynamic responses of surrounding rocks of a near-fault tunnel under the action of an earthquake, which is provided by the embodiment of the invention;
FIG. 10 is a graph of baseline correction of acceleration response time course of near-fault tunnel surrounding rock stress field and dynamic response under analysis of seismic action provided by an embodiment of the invention;
FIG. 11 is a corrected acceleration response time chart for analyzing the stress field and dynamic response of near-fault tunnel surrounding rock under the action of an earthquake provided by an embodiment of the invention;
FIG. 12 is a graph of artificially synthesized earthquake velocity response time for analyzing stress fields and dynamic responses of surrounding rocks of near-fault tunnels under the action of earthquake, provided by the embodiment of the invention;
FIG. 13 is a graph of artificially synthesized earthquake motion displacement response time for analyzing stress fields and dynamic responses of surrounding rocks of near-fault tunnels under the action of earthquake, provided by the embodiment of the invention;
FIG. 14 is a corrected velocity response time chart for analyzing the stress field and dynamic response of near-fault tunnel surrounding rock under the action of an earthquake provided by an embodiment of the invention;
FIG. 15 is a graph of corrected displacement response time for analyzing the stress field and dynamic response of near-fault tunnel surrounding rock under the action of an earthquake provided by an embodiment of the invention;
FIG. 16 is a Fourier spectrum of artificial synthetic earthquake motion for analyzing stress fields and dynamic responses of near-fault tunnel surrounding rock under the action of earthquake provided by the embodiment of the invention.
Detailed Description
The following detailed description of the technical solutions of the present invention is made by the accompanying drawings and specific embodiments, and it should be understood that the specific features of the embodiments and embodiments of the present application are detailed descriptions of the technical solutions of the present application, and not limiting the technical solutions of the present application, and the technical features of the embodiments and embodiments of the present application may be combined with each other without conflict. The following examples are only for more clearly illustrating the technical aspects of the present invention, and are not intended to limit the scope of the present invention.
The term "and/or" is herein merely an association relationship describing an associated object, meaning that there may be three relationships, e.g., a and/or B, may represent: a exists alone, A and B exist together, and B exists alone. In addition, the character "/" herein generally indicates that the front and rear associated objects are an "or" relationship.
Example 1
FIG. 1 is a flow chart of a method of numerical analysis of near-fault tunnel stress field and seismic dynamic response in accordance with example 1 of the present invention. The flow chart merely shows the logical sequence of the method according to the present embodiment, and the steps shown or described may be performed in a different order than shown in fig. 1 in other possible embodiments of the invention without mutual conflict.
This embodiment is an exemplary implementation of the present invention, and provides a method for analyzing the values of the stress field and the seismic dynamic response of a tunnel near fault, which can be applied to a terminal and can be executed by an electronic terminal, which can be implemented by software and/or hardware, and which can be integrated in the terminal, for example: any smart phone, tablet computer or computer device with communication function. As shown in fig. 1, the method of this embodiment specifically includes the following steps:
Step 1: based on the Flac3D6.0 platform, establishing a numerical model of a complete rock sample without a structural surface and a numerical model of a rock sample with a structural surface, wherein the numerical model of the rock sample with the structural surface comprises the structural surface, obtaining deformation parameters of the structural surface by comparing the numerical model of the complete rock sample without the structural surface with the numerical model of the rock sample with the structural surface, applying a speed boundary condition according to the deformation parameters of the structural surface, and applying a normal stiffness K of the structural surface n And shear stiffness K s Obtaining a fault contact surface model by using a Flac3D software model library or a finite difference method, and calibrating the normal rigidity K of the structural surface to the fault contact surface model n And shear stiffness K s ；
Step 2: according to the normal rigidity of the structural surface and the shear rigidity of the structural surface, a geomechanical model is constructed, and the geomechanical model comprises a plurality of units, wherein each unit comprises a plurality of nodes and a plurality of surface units; assuming that surrounding rock materials in the geomechanical model obey a line elastic strength criterion, and that the fault contact surface model obeys a Moire-Coulomb strength criterion; giving deformation parameters of a structural surface to the geomechanical model, applying normal constraint conditions and ground stress boundary conditions according to boundary coordinates of the geomechanical model, performing simulation calculation on the geomechanical model, setting displacement and speed of nodes of the geomechanical model to be zero, and reserving calculation results of an initial ground stress field;
Step 3: s2, constructing a near-fault tunnel model and excavating a tunnel based on the calculation result of the initial ground stress field in the step S; simulating a secondary stress field of the near-fault tunnel model after tunnel excavation and calculating a stress concentration coefficient;
step 4: based on a Matlab platform, according to the related requirements on tunnel earthquake-resistant design in highway tunnel earthquake-resistant design specification (JTG 2232-2019) and underground structure earthquake-resistant design standard (GB/T51336-2018), artificial earthquake waves meeting the requirements of the design specification are obtained by synthesizing and designing artificial earthquake acceleration time course according to an earthquake acceleration response spectrum and performing filtering and baseline correction;
step 5: and (3) applying a viscoelastic artificial boundary and a damping coefficient to the geomechanical model based on the calculation result of the stress concentration coefficient of the secondary stress field in the step (3), applying artificial seismic waves to the geomechanical model boundary as seismic power boundary conditions, recording the peak acceleration and the peak stress concentration coefficient of the unit in the reaction process, and analyzing the dynamic response of the geomechanical model to obtain a stress field diagram distributed by the peak stress concentration coefficient of the near-fault tunnel model and a horizontal acceleration reaction time course.
In step 1, a numerical model of a complete rock sample without a structural surface and a numerical model of a rock sample with a structural surface are respectively built based on a Flac3D6.0 platform, a coordinate axis z is taken as a rotation axis, a circular plane is positioned on an xoy coordinate plane, and a radius r=1m and a height h=4m are set. The numerical model of the complete rock sample without structural face uses the following material parameters: e=5gpa, v=0.20, c=14 MPa,the structural material adopts the following material parameters: e=2gpa, v=0.25, c=3 MPa, ++>
Generating an upper boundary and a lower boundary, applying a rate of 10 to the upper boundary and the lower boundary -5 Speed boundary conditions of m/s; the realization is realized by adopting a fish language:
generating a global mapping variable surface;
the numerical model of the complete rock sample without the structural surface and the numerical model of the rock sample with the structural surface comprise node bodies, the node bodies are traversed, and if a node body gp is positioned at the upper boundary, the node bodies are saved to a surface;
applying (1-10) x 10 to global map variable surface 5 m/s speed boundary condition, simulating the loading rate of a uniaxial compression test;
Traversing the node body, recording the unbalanced force in the z direction on the global mapping variable surface by using a function gp.force.unbal.z (gp) of the maximum unbalanced force, and storing the unbalanced force in the z direction to a variable zforce;
calculate normal stress AxialStress = zforce/pi r 2 ；
Setting artificial convergence conditions, and calculating convergence when the displacement of the numerical model of the complete rock sample without the structural surface and the numerical model of the rock sample with the structural surface reaches the termination strain.
When computing convergence, sigma n -ΔV r Sum sigma n -ΔV t The curve is shown in fig. 2. According to DeltaV j ＝ΔV t -ΔV r To obtain DeltaV j Further, the sigma is plotted n -ΔV j Curves, as shown in figure 3. In the present embodiment, V m =0.03 mm. From the following componentsAndcalculated to obtain K n =20 MPa/cm, by->Calculated to obtain K s ＝8MPa/cm。
Specifically, in the step 1, the model parameters of the fault contact surface are as follows: k (K) n ＝20MPa/cm，K s =8mpa/cm. In the present embodiment, it is assumed that the fault inclination angle β=60°, and the distance d=2r of the fault from the tunnel 0 Horizontal and vertical stress straight p h ＝p v =5 MPa, as shown in fig. 4.
In this embodiment, a structured grid numerical model of a straight wall semicircular arch tunnel is built, with radius r=1.2m, as shown in fig. 5. And (5) establishing a xoz plane model of the tunnel section according to the real size, and enabling the geometric center of the model to be located at the origin of coordinates. The thickness of the model was set to 1.0. Surrounding rocks with the range of 4-6 times of tunnel radius are arranged on the periphery of an established tunnel model, the upper disc and the lower disc of the surrounding rocks are respectively named in groups, surface units are created between the upper disc and the lower disc of a fault to generate a fault contact surface model, and K is given to the fault contact surface model n 、K s Mechanical parameters of the iso-fault. And monitoring points are respectively arranged on the vault, the arch shoulder, the straight wall, the arch foot and the bottom plate of the tunnel.
Further, deformation parameters of the structural surface are given to the geomechanical model, normal constraint conditions and ground stress boundary conditions are applied according to boundary coordinates of the geomechanical model, initial ground stress balance is calculated, displacement and speed are cleared, and a calculation result of an initial ground stress field is reserved.
Specifically, in step 3, a tunnel model is excavated on the basis of the calculation result of the initial ground stress field. And simulating a secondary stress field of the surrounding rock model after tunnel excavation, calculating stress concentration coefficients, and comparing and analyzing displacement stress fields of the non-fault tunnel and the near-fault tunnel, as shown in fig. 6 and 7.
Specifically, in step 4, the analysis object of the embodiment is a near-fault tunnel at a certain burial depth. The basic intensity of earthquake is VIII, and the type of earthquake fortification is B according to the highway tunnel earthquake proof design rule (JTG 2232-2019) and the underground structure earthquake proof design standard (GB/T51336-2018). The fortification target is that the E1 earthquake action (75 years of reproduction period) should not be damaged, and the probability of exceeding 50 years is 63%. Shock-resistant importance coefficient C i Ground vibration peak acceleration a=0.3 g, ground vibration peak acceleration adjustment coefficient C s =1.0, the earth's surface level is to design the seismic peak acceleration a sh ＝C s C i A=0.13 g. The peak acceleration of underground horizontal to design earthquake motion is 1/2 of the peak acceleration of earth surface horizontal to design earthquake motion, namely A uh ＝A sh ＝0.07g。
Further, byCalculated acceleration plateau start time t 1 =3s, acceleration plateau end time t 2 =17s, attenuation coefficient c=0.16. Structural damping ratio ζ=0.05, damping adjustment coefficient C d =1.0, then the horizontal acceleration response spectrum maximum S max ＝C d A uh =0.18 g. Characteristic period T g =0.4 s, according to ∈>The plot of the design response spectrum is shown in figure 8.
Further, a synthetic artificial seismic main routine is run, in which an intensity envelope function subroutine and a power density function subroutine are invoked. The acceleration of the reaction time course is shown in FIG. 9.
Further, the baseline correction main routine is run, the sampling sliding window algorithm calculates running average acceleration, and sets the sliding window size k=10, as shown in fig. 10. The time-course conversion subroutine is called to display the acceleration after baseline correction as shown in fig. 11, and the velocity and displacement reaction time course before and after correction as shown in fig. 12 to 15.
Further, a low-pass filter main routine is run, setting the cut-off frequency to 15Hz. The fourier transform subroutine is called to obtain the corrected pre-and post-fourier spectra as shown in fig. 16.
Specifically, in step 5, a viscoelastic artificial boundary and a damping coefficient are applied to the model, and a seismic dynamic boundary condition is applied to the model boundary according to the relative position of the seismic source and the tunnel engineering.
Preferably, in the present embodiment, the pulling is calculated and applied in each time step in the same manner as the boundary load is appliedAttraction t n And t s 。
Preferably, in this embodiment, local damping is used to simulate structural damping in the dynamic response. The mass of the mesh or node changes at a particular time in the vibration cycle, but the overall mass is conserved. When the direction of the velocity changes, the mass increases; when the direction of acceleration changes, the mass decreases. In each vibration cycle, the kinetic energy increment is lost at the speed extremum, the lost energy is proportional to the maximum instantaneous strain energy, the ratio is independent of speed and frequency, and only critical damping is involved. In this embodiment, the critical damping ζ=0.05 of the structure, then the model local damping α L ＝πζ＝0.157。
Further, the artificial seismic waves synthesized in the step 4 are applied to the model as dynamic boundary conditions, and the calculation is finished when the dynamic time step is 30 s. And respectively drawing peak stress concentration coefficient stress field diagrams of the non-fault tunnel and the near-fault tunnel. And drawing a horizontal acceleration response time course according to the record result of the monitoring point.
Example 2
The embodiment of the invention also provides an electronic terminal, which comprises a processor and a storage medium; the storage medium is used for storing instructions; the processor is configured to operate according to the instructions to perform the steps of the method of:
step 1: based on the Flac3D6.0 platform, establishing a numerical model of a complete rock sample without a structural surface and a numerical model of a rock sample with a structural surface, wherein the numerical model of the rock sample with the structural surface comprises the structural surface, obtaining deformation parameters of the structural surface by comparing the numerical model of the complete rock sample without the structural surface with the numerical model of the rock sample with the structural surface, applying a speed boundary condition according to the deformation parameters of the structural surface, and applying a normal stiffness K of the structural surface n And shear stiffness K s Obtaining a fault contact surface model by using a Flac3D software model library or a finite difference method, and calibrating the normal rigidity K of the structural surface to the fault contact surface model n And shear stiffness K s ；
Step 2: according to the normal rigidity of the structural surface and the shear rigidity of the structural surface, a geomechanical model is constructed, and the geomechanical model comprises a plurality of units, wherein each unit comprises a plurality of nodes and a plurality of surface units; assuming that surrounding rock materials in the geomechanical model obey a line elastic strength criterion, and that the fault contact surface model obeys a Moire-Coulomb strength criterion; giving deformation parameters of a structural surface to the geomechanical model, applying normal constraint conditions and ground stress boundary conditions according to boundary coordinates of the geomechanical model, performing simulation calculation on the geomechanical model, setting displacement and speed of nodes of the geomechanical model to be zero, and reserving calculation results of an initial ground stress field;
step 3: s2, constructing a near-fault tunnel model and excavating a tunnel based on the calculation result of the initial ground stress field in the step S; simulating a secondary stress field of the near-fault tunnel model after tunnel excavation and calculating a stress concentration coefficient;
step 4: based on a Matlab platform, according to the related requirements on tunnel earthquake-resistant design in highway tunnel earthquake-resistant design specification (JTG 2232-2019) and underground structure earthquake-resistant design standard (GB/T51336-2018), artificial earthquake waves meeting the requirements of the design specification are obtained by synthesizing and designing artificial earthquake acceleration time course according to an earthquake acceleration response spectrum and performing filtering and baseline correction;
Step 5: and (3) applying a viscoelastic artificial boundary and a damping coefficient to the geomechanical model based on the calculation result of the stress concentration coefficient of the secondary stress field in the step (3), applying artificial seismic waves to the geomechanical model boundary as seismic power boundary conditions, recording the peak acceleration and the peak stress concentration coefficient of the unit in the reaction process, and analyzing the dynamic response of the geomechanical model to obtain a stress field diagram distributed by the peak stress concentration coefficient of the near-fault tunnel model and a horizontal acceleration reaction time course.
The electronic terminal provided by the embodiment of the invention can execute the numerical analysis method of the stress field and the earthquake dynamic response of the near-fault tunnel provided by any embodiment of the invention, and has the corresponding functional modules and beneficial effects of the execution method.
Example 3
The embodiment of the invention also provides a computer readable storage medium having stored thereon a computer program which when executed by a processor performs the steps of:
step 1: based on the Flac3D6.0 platform, establishing a numerical model of a complete rock sample without a structural surface and a numerical model of a rock sample with a structural surface, wherein the numerical model of the rock sample with the structural surface comprises the structural surface, obtaining deformation parameters of the structural surface by comparing the numerical model of the complete rock sample without the structural surface with the numerical model of the rock sample with the structural surface, applying a speed boundary condition according to the deformation parameters of the structural surface, and applying a normal stiffness K of the structural surface n And shear stiffness K s Obtaining a fault contact surface model by using a Flac3D software model library or a finite difference method, and calibrating the normal rigidity K of the structural surface to the fault contact surface model n And shear stiffness K s ；
Step 2: according to the normal rigidity of the structural surface and the shear rigidity of the structural surface, a geomechanical model is constructed, and the geomechanical model comprises a plurality of units, wherein each unit comprises a plurality of nodes and a plurality of surface units; assuming that surrounding rock materials in the geomechanical model obey a line elastic strength criterion, and that the fault contact surface model obeys a Moire-Coulomb strength criterion; giving deformation parameters of a structural surface to the geomechanical model, applying normal constraint conditions and ground stress boundary conditions according to boundary coordinates of the geomechanical model, performing simulation calculation on the geomechanical model, setting displacement and speed of nodes of the geomechanical model to be zero, and reserving calculation results of an initial ground stress field;
step 3: s2, constructing a near-fault tunnel model and excavating a tunnel based on the calculation result of the initial ground stress field in the step S; simulating a secondary stress field of the near-fault tunnel model after tunnel excavation and calculating a stress concentration coefficient;
step 4: based on a Matlab platform, according to the related requirements on tunnel earthquake-resistant design in highway tunnel earthquake-resistant design specification (JTG 2232-2019) and underground structure earthquake-resistant design standard (GB/T51336-2018), artificial earthquake waves meeting the requirements of the design specification are obtained by synthesizing and designing artificial earthquake acceleration time course according to an earthquake acceleration response spectrum and performing filtering and baseline correction;
Step 5: and (3) applying a viscoelastic artificial boundary and a damping coefficient to the geomechanical model based on the calculation result of the stress concentration coefficient of the secondary stress field in the step (3), applying artificial seismic waves to the geomechanical model boundary as seismic power boundary conditions, recording the peak acceleration and the peak stress concentration coefficient of the unit in the reaction process, and analyzing the dynamic response of the geomechanical model to obtain a stress field diagram distributed by the peak stress concentration coefficient of the near-fault tunnel model and a horizontal acceleration reaction time course.
The computer readable storage medium provided by the embodiment of the invention is stored with a computer program, and can execute the numerical analysis method of the tunnel stress field and the earthquake dynamic response of the near fault provided by any embodiment of the invention, and the method has the corresponding functional modules and beneficial effects of the execution method.
In the description of the present invention, it should be understood that the terms "center", "longitudinal", "lateral", "upper", "lower", "front", "rear", "left", "right", "vertical", "horizontal", "top", "bottom", "inner", "outer", etc. indicate orientations or positional relationships based on the orientations or positional relationships shown in the drawings, are merely for convenience in describing the present invention and simplifying the description, and do not indicate or imply that the devices or elements referred to must have a specific orientation, be configured and operated in a specific orientation, and thus should not be construed as limiting the present invention. Furthermore, the terms "first," "second," and the like, are used for descriptive purposes only and are not to be construed as indicating or implying a relative importance or implicitly indicating the number of technical features indicated. Thus, a feature defining "a first", "a second", etc. may explicitly or implicitly include one or more such feature. In the description of the present invention, unless otherwise indicated, the meaning of "a plurality" is two or more.
In the description of the present invention, it should be noted that, unless explicitly specified and limited otherwise, the terms "mounted," "connected," and "connected" are to be construed broadly, and may be either fixedly connected, detachably connected, or integrally connected, for example; can be mechanically or electrically connected; can be directly connected or indirectly connected through an intermediate medium, and can be communication between two elements. The specific meaning of the above terms in the present invention can be understood by those of ordinary skill in the art in a specific case.
The foregoing is merely a preferred embodiment of the present invention, and it should be noted that modifications and variations could be made by those skilled in the art without departing from the technical principles of the present invention, and such modifications and variations should also be regarded as being within the scope of the invention.

Claims (7)

1. The numerical analysis method of the stress field of the tunnel and the earthquake dynamic response of the near fault is characterized by comprising the following steps:
s1: the method comprises the steps of pre-constructing a numerical model of a complete rock sample without a structural surface and a numerical model of a rock sample with the structural surface, wherein the numerical model of the rock sample with the structural surface comprises the structural surface, obtaining deformation parameters of the structural surface by comparing the numerical model of the complete rock sample without the structural surface with the numerical model of the rock sample with the structural surface, applying a speed boundary condition according to the deformation parameters of the structural surface, calculating normal stiffness and shear stiffness parameters of the structural surface, obtaining a fault contact surface model from a model library, and calibrating the normal stiffness of the structural surface and the shear stiffness of the structural surface for the fault contact surface model;
S2: constructing a geomechanical model according to the normal stiffness of the structural surface and the shear stiffness of the structural surface, wherein the geomechanical model comprises a plurality of units, and each unit comprises a plurality of nodes and a plurality of surface units; the deformation parameters of the structural surface are given to the geomechanical model, normal constraint conditions and ground stress boundary conditions are applied according to the boundary coordinates of the geomechanical model, simulation calculation is carried out on the geomechanical model, the displacement and the speed of the nodes of the geomechanical model are set to be zero, and an initial environment is provided for calculation of a next model;
s3: constructing a near-fault tunnel model based on the initial environment and excavating a tunnel; simulating a secondary stress field of the near-fault tunnel model after tunnel excavation and calculating a stress concentration coefficient;
s4: synthesizing and designing an artificial earthquake motion acceleration time course according to the earthquake motion acceleration response spectrum, and filtering and correcting a base line to obtain an artificial earthquake wave meeting the requirement of a design specification;
s5: and applying a viscoelastic artificial boundary and a damping coefficient to the geomechanical model based on the calculation result of the stress concentration coefficient of the secondary stress field, applying the artificial seismic wave to the geomechanical model boundary as a seismic power boundary condition, recording the peak acceleration and the peak stress concentration coefficient of the unit in the reaction process, and analyzing the dynamic response of the geomechanical model to obtain a peak stress concentration coefficient distribution stress field diagram and a horizontal acceleration reaction time course of the near-fault tunnel model.
2. The method for analyzing the tunnel stress field and the seismic dynamic response of the near fault according to claim 1, wherein the step S1 specifically comprises:
s1.1: generating a numerical model of a complete rock sample without a structural surface and a numerical model of a rock sample with a structural surface by adopting software, and setting the ratio of the height to the diameter of the numerical model of the complete rock sample without the structural surface and the numerical model of the rock sample with the structural surface to be 2:0;
s1.2: comparing the numerical model of the complete rock sample without the structural surface with the numerical model of the rock sample with the structural surface to obtain the deformation parameters of the structural surface, wherein the normal deformation of the numerical model of the complete rock sample without the structural surface is always smaller than that of the numerical model of the rock sample with the structural surface, and the stress-deformation curve sigma in the numerical model of the complete rock sample without the structural surface n -ΔV j The slope of the slope is larger than the numerical model of the rock sample containing the structural plane, and a normal constraint condition is applied to the boundary of the numerical model of the rock sample containing the structural plane according to the boundary coordinates of the numerical model of the rock sample containing the structural plane;
s1.3: simulationThe loading rate of the uniaxial compression test is that a speed boundary condition is applied to a numerical model of a complete rock sample without a structural surface and a numerical model of a rock sample with a structural surface, respectively, a round plane is positioned on an xoy coordinate plane by taking a coordinate axis z as a rotation axis, an upper boundary and a lower boundary are generated, and the application rate of the upper boundary and the lower boundary is 10 -5 Speed boundary conditions of m/s; the method adopts programming language to realize simulation uniaxial compression test, and comprises the following steps:
s1.3.1: generating a global mapping variable surface;
s1.3.2: the numerical model of the complete rock sample without the structural surface and the numerical model of the rock sample with the structural surface comprise node bodies, the node bodies are traversed, and if a node body gp is positioned at the upper boundary, the node bodies are saved to a surface;
s1.3.3: applying (1-10) x 10 to global map variable surface 5 m/s speed boundary condition, simulating the loading rate of a uniaxial compression test;
s1.3.4: traversing the node body, recording the unbalanced force in the z direction on the global mapping variable surface by using a function gp.force.unbal.z (gp) of the maximum unbalanced force, and storing the unbalanced force in the z direction to a variable zforce;
s1.3.5: calculate normal stress AxialStress = zforce/pi r 2 ；
Wherein AxialStress represents normal stress; zforce represents the z-direction imbalance force; r represents the radius of the numerical model of the complete rock sample without structural face and the numerical model of the rock sample with structural face; pi represents the circumference ratio;
s1.3.6: setting a manual convergence condition, namely a condition that the software considers that the calculation is finished when calculating, and finishing the calculation when the displacement of the numerical model of the complete rock sample without the structural surface and the numerical model of the rock sample with the structural surface reaches the termination strain;
Obtaining deformation DeltaV of a complete rock sample r Deformation DeltaV of rock sample containing structural surface t Calculating normal closing deformation DeltaV of structural surface j ：
ΔV j ＝ΔV t -ΔV r (1)
In DeltaV r Representing deformation of the complete rock sample; deltaV t Representing deformation of a rock sample containing a structural face; deltaV j Representing normal closing deformation of the structural face;
s1.4: fitting σ using the equation n -ΔV j Stress-deformation relationship curve:
v in m Representing the maximum normal closing deformation of the structural surface, sigma n -ΔV j An asymptote of the stress-deformation relationship curve; sigma (sigma) n Representing normal stress;representing the normal stiffness of the structural face fitted by the equation;
s1.5: solving the tangential slope K of the asymptote equation at the origin according to the equation fitted by equation (2) ni I.e. the initial normal stiffness of the structural face:
s1.6: calibrating the normal stress sigma according to the rock sample parameters containing the structural surface n Under the action, the normal rigidity of the structural surface and the shear rigidity of the structural surface of the fault contact surface model:
wherein L is the trend length of the fault; JRC is a fault structural planeRoughness coefficient of (c); JCS is the wall rock strength of the structural face; phi (phi) r Residual friction angle of the fault structural surface; k (K) n Is the normal rigidity of the structural surface; k (K) s Is the shear stiffness of the structural face.
3. The method for analyzing the tunnel stress field and the seismic dynamic response of the near fault according to claim 1, wherein the step S2 specifically comprises:
S2.1: if the section of the tunnel is simple in shape, such as a round shape and a semicircular arch shape, a structural grid is adopted to build a model; if the section shape of the tunnel is complex, such as a straight wall three-heart arch shape and a horseshoe shape, an unstructured grid is adopted to build a model; establishing a tunnel model according to the real size, enabling the geometric center to be located at the origin of coordinates, and setting the thickness to be 1.0;
s2.2: surrounding rocks with the range of 4-6 times of tunnel radius are arranged on the periphery of an established tunnel model, the surrounding rocks comprise an upper disc and a lower disc, the upper disc and the lower disc are respectively named in a grouping way, and a surrounding rock model is constructed;
s2.3: constructing a fault contact surface model based on a finite difference numerical platform; imparting a structural plane normal stiffness K to a fault interface model n Shear stiffness of structural face K s 。
S2.4: the tunnel model, the surrounding rock model and the fault contact surface model jointly form a geomechanical model, a normal constraint condition and a ground stress boundary condition are applied according to boundary coordinates of the geomechanical model, the normal speed of a surface unit of the geomechanical model boundary is set to be 0, and horizontal stress and vertical stress of 0-50 MPa are applied to the geomechanical model.
4. The method for analyzing the stress field and the seismic dynamic response of the near-fault tunnel according to claim 1, wherein the step S3 adopts a programming language to simulate the secondary stress field of the near-fault tunnel model after the tunnel is excavated, and specifically comprises the following steps:
S3.1: traversing the near-fault tunnel model by using a programming language, and extracting the stress state of each unit zp:
middle sigma x The horizontal stress of the unit body in the numerical calculation model is calculated; sigma (sigma) z The vertical stress of the unit body in the numerical calculation model is calculated; τ xz The shear stress of the unit body in the numerical calculation model is calculated; zp is the pointer variable of each unit; xx is stress in the x direction of the unit body, yy is stress in the y direction of the unit body, and xy is shear stress of the unit body; zone is a function of the unit cell; stress is the stress function of the unit cell;
s3.2: calculating the distance R from the center of each unit body to the geometric center point of the tunnel:
wherein pos is a position coordinate vector of the unit body;
s3.3: taking a straight wall semicircular arch as an example without losing generality, according to the section shape of the tunnel, the sine and cosine of the connecting line from the center of each unit body to the origin of coordinates are calculated by using a programming language:
s3.3.1: the unit body is positioned at the arc-shaped part, namely zone.pos.z (zp) > 0.0
S3.3.2: the unit body is positioned at the straight wall part, namely-1.2 < zone.pos.z (zp) < 0.0
S3.3.3: the unit body is located at the corner site, i.e. zone.pos.z (zp) < -R, and math.abs (zone.pos.x (zp)) > R
S3.3.4: the unit body is located at the bottom plate position, i.e. zone.pos.z (zp) < -R, and math.abs (zone.pos.x (zp)) < R
S3.4: calculating radial stress sigma of surrounding rocks around the tunnel according to a rotating shaft formula r And hoop stress sigma θ ：
S3.5: defining stress concentration coefficient K as hoop stress sigma θ And vertical ground stress p v Solving the K of the unit zp of the surrounding rock:
K＝σ θ /p v (13)
the calculation result is saved to the variable zone. Extra (zp).
5. The method for analyzing the numerical value of the stress field and the earthquake dynamic response of the tunnel near fault as claimed in claim 1, wherein the artificial earthquake motion acceleration time course is synthesized by adopting a trigonometric series method and is realized by using a matrix calculation numerical platform, and the method specifically comprises the following steps:
s4.1: creating an acceleration main program for synthesizing an artificial earthquake motion reaction time course, and setting variables for saving the reaction duration T, the reaction interval dT and the upper limit omega of a frequency spectrum u Lower spectral limit omega l Frequency domain sampling interval Δω; generating an empty array variable acc to store the acceleration of the artificial earthquake motion reaction time course;
s4.2: according to the random vibration theory, inputting the maximum reaction distribution of a simple substance point system in a stable random process:
in which a is g (t) is a smooth random process; initial phaseIs a random phase angle uniformly distributed in (0, 2 pi); omega k A circular frequency that is the kth frequency; length is a length metric function; a is that k (ω k ) For the amplitude of the kth frequency, it is calculated by:
s (omega) k ) As a function of power spectral density, Δω is the frequency sampling interval, calculated by:
Δω＝2π/T (16)
wherein T is the duration of the reaction;
s4.3: creating an intensity envelope function subprogram 1, calling the intensity envelope function subprogram in an acceleration main program for synthesizing an artificial earthquake motion reaction time course, transmitting a time variable t to the subprogram, and returning an intensity envelope function value f (t); calculating the starting time t of acceleration plateau in the intensity envelope function according to the equivalent shock level M and the equivalent shock center distance D 1 Acceleration plateau end time t 2 And attenuation coefficient c:
the intensity envelope function is constructed as follows:
s4.4: creating a power density function subroutine in the synthesis of artificial earthquakesCalling in a main program of dynamic response time acceleration; introducing a frequency variable omega into a subroutine k Return power spectral density function value S (omega) k ) The method comprises the steps of carrying out a first treatment on the surface of the Setting variables for saving the self-vibration period T and the field characteristic period T of the structure g Damping ratio ζ, acceleration response spectrum maximum S max The decay index gamma of the descending segment and the reaction overrun probability P; reaction spectrum S for synthesizing design acceleration time course a The method comprises the following steps:
the power spectral density function is calculated by:
S4.5: and (3) operating a main program for synthesizing the acceleration of the artificial earthquake motion reaction time course, calling the subprogram, synthesizing the artificial earthquake motion acceleration time course, and storing the artificial earthquake motion acceleration time course into an array variable acc:
s4.6: creating a low-pass filtering main program, and generating an array variable input_acc to store an artificial earthquake motion acceleration time course acc; generating variable preserving filter order n and cut-off frequency omega n ；
Using the Butterworth low-pass filter function of the matrix calculation numerical platform to generate transfer function coefficients a and b of the filter, and returning a row vector H (z) with the length of (n+1):
filtering through transfer function coefficients, and returning the acceleration reaction time interval filtered_acc after filtering;
s4.7: creating a base line correction main program, and generating an array variable iput_acc to store a filtered acceleration response time interval filtered_acc; setting the size k of a sliding window, and returning to an array of local k-point average values, wherein each average value is calculated on adjacent elements of the iput_acc through the sliding window with the length k; when k is even, the window is centered on the current and neighboring elements; when there are not enough elements to fill the window, the window size will automatically truncate at the endpoint; when the window is truncated, only averaging the elements filling the window; returning to the running average acceleration smoothened_a, the baseline corrected acceleration response time course modified_a may be calculated by:
modified_a＝iput_acc-smoothed_a (23)
S4.8: creating a time-course conversion subprogram which can be called in any main program; converting the acceleration reaction time course into a speed reaction time course and a displacement reaction time course, and completely describing the reaction time course for synthesizing the artificial earthquake motion; the acceleration response time course is transmitted to the subprogram and stored as an array variable iput_acc, and a speed variable vel and a displacement variable disp are initialized:
wherein zeros is an initializing matrix function;
calculating a velocity reaction time course and a displacement reaction time course by adopting numerical integration as return values:
wherein vel (i) is a velocity response time course; disp (i) is the shift reaction schedule;
s4.9: creating a Fourier transform subprogram which can be called in any main program; converting the time domain earthquake motion signal into a frequency domain signal by adopting a fast Fourier transform function of a matrix calculation numerical platform for analysis; and (4) transmitting an acceleration response time course to the subroutine, saving the acceleration response time course as an array variable iput_acc, and returning an acceleration amplitude frequency spectrum.
6. An electronic terminal is characterized by comprising a processor and a storage medium;
the storage medium is used for storing instructions;
the processor being operative according to the instructions to perform the steps of the method according to any one of claims 1 to 5.
7. A computer readable storage medium, on which a computer program is stored, characterized in that the program, when being executed by a processor, implements the steps of the method according to any one of claims 1-6.
CN202311813634.4A 2023-12-27 2023-12-27 Numerical analysis method, terminal and medium for stress field and earthquake dynamic response of near-fault tunnel Active CN117763914B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202311813634.4A CN117763914B (en) 2023-12-27 2023-12-27 Numerical analysis method, terminal and medium for stress field and earthquake dynamic response of near-fault tunnel

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202311813634.4A CN117763914B (en) 2023-12-27 2023-12-27 Numerical analysis method, terminal and medium for stress field and earthquake dynamic response of near-fault tunnel

Publications (2)

Publication Number Publication Date
CN117763914A true CN117763914A (en) 2024-03-26
CN117763914B CN117763914B (en) 2024-06-21

Family

ID=90312128

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202311813634.4A Active CN117763914B (en) 2023-12-27 2023-12-27 Numerical analysis method, terminal and medium for stress field and earthquake dynamic response of near-fault tunnel

Country Status (1)

Country Link
CN (1) CN117763914B (en)

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106198208A (en) * 2016-06-21 2016-12-07 中电建路桥集团有限公司 A kind of interbedding of soft and hard rocks surrounding rock tunnel method for analog construction considering layer reason effect
US20170284911A1 (en) * 2016-03-31 2017-10-05 Powerchina Huadong Engineering Corporation Limited Integrated style shear apparatus for rock structural plane and a shear experimental method for rock structural plane
US20210263003A1 (en) * 2019-07-15 2021-08-26 China University Of Mining And Technology Discrete element method for modelling a fracture evolution of a roadway surrounding rock
CN113435087A (en) * 2021-06-24 2021-09-24 中铁二院工程集团有限责任公司 Method for analyzing local stability of cave surrounding rock
CN113919196A (en) * 2021-09-26 2022-01-11 中国石油大学(华东) Reservoir three-dimensional stress field simulation method, simulation system, terminal and storage medium
CN114755310A (en) * 2022-04-26 2022-07-15 中国地质大学（武汉） Method for predicting evolution rule of rock mechanics layer of fractured reservoir
CN115659741A (en) * 2022-10-25 2023-01-31 温州大学 SV wave oblique incidence lower mountain tunnel portal section structure seismic response analysis method

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20170284911A1 (en) * 2016-03-31 2017-10-05 Powerchina Huadong Engineering Corporation Limited Integrated style shear apparatus for rock structural plane and a shear experimental method for rock structural plane
CN106198208A (en) * 2016-06-21 2016-12-07 中电建路桥集团有限公司 A kind of interbedding of soft and hard rocks surrounding rock tunnel method for analog construction considering layer reason effect
US20210263003A1 (en) * 2019-07-15 2021-08-26 China University Of Mining And Technology Discrete element method for modelling a fracture evolution of a roadway surrounding rock
CN113435087A (en) * 2021-06-24 2021-09-24 中铁二院工程集团有限责任公司 Method for analyzing local stability of cave surrounding rock
CN113919196A (en) * 2021-09-26 2022-01-11 中国石油大学(华东) Reservoir three-dimensional stress field simulation method, simulation system, terminal and storage medium
CN114755310A (en) * 2022-04-26 2022-07-15 中国地质大学（武汉） Method for predicting evolution rule of rock mechanics layer of fractured reservoir
CN115659741A (en) * 2022-10-25 2023-01-31 温州大学 SV wave oblique incidence lower mountain tunnel portal section structure seismic response analysis method

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
张志国;肖明;陈俊涛;: "大型地下洞室地震灾变过程三维动力有限元模拟", 岩石力学与工程学报, no. 03, 15 March 2011 (2011-03-15) *

Also Published As

Publication number Publication date
CN117763914B (en) 2024-06-21

Similar Documents

Publication Publication Date Title
Chen et al. An efficient nonlinear octree SBFEM and its application to complicated geotechnical structures
CN108710153B (en) Wave number domain method for magnetic full tensor gradient inversion underground three-dimensional magnetic distribution
Ford et al. A nonhydrostatic finite-element model for three-dimensional stratified oceanic flows. Part I: Model formulation
Casarotti et al. CUBIT and seismic wave propagation based upon the spectral-element method: An advanced unstructured mesher for complex 3D geological media
CN105138731B (en) A kind of decomposition of hydrate causes submarine slope unstability evaluation system and method
Mejia et al. Dynamic analysis of earth dams in three dimensions
CN112596106B (en) Method for carrying out gravity and vibration joint inversion on density interface distribution under spherical coordinate system
CN111553107A (en) Liquefiable site pile foundation random earthquake response analysis and safety evaluation method
CN108983285A (en) A kind of a variety of source wavefield analogy methods and device based on moment tensor
CN105717547A (en) Anisotropy medium magnetotelluric meshless value simulation method
CN114139335B (en) Viscous sound wave simulation method based on single relaxation time lattice Boltzmann model
CN104597496B (en) A kind of three dimensions method for homing of 2-d seismic data medium velocity
CN117057110A (en) Construction method of P-wave induced elliptic tunnel surrounding rock dynamic stress concentration coefficient calculation model
CN115795974A (en) Air-ground integrated monitoring method for evolution of underground water seepage field of deep-buried tunnel
CN117725354B (en) Rapid forward and backward modeling method and system combining large data volume gravity and gravity gradient
CN114167511B (en) Bit field data rapid inversion method based on continuous expansion downward continuation
CN117454725B (en) Offshore wind power foundation seismic load simulation method and equipment based on superunit condensation
CN117763914B (en) Numerical analysis method, terminal and medium for stress field and earthquake dynamic response of near-fault tunnel
CN109444973A (en) Gravity forward modeling accelerated method under a kind of spherical coordinate system
Cullen On the accuracy of the semi‐geostrophic approximation
Zhang et al. Rotation errors in numerical manifold method and a correction based on large deformation theory
Khoei et al. SUT-DAM: An integrated software environment for multi-disciplinary geotechnical engineering
CN115392032A (en) GIS-MPM seamless integrated dynamic three-dimensional geological model construction method
Shah et al. Soil structure interaction analysis methods-State of art-Review
Mandal et al. 2d finite element analysis of rectangular water tank with separator wall using direct coupling

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant