CN115906704A - Single-joint stress-seepage coupling numerical calculation method - Google Patents
Single-joint stress-seepage coupling numerical calculation method Download PDFInfo
- Publication number
- CN115906704A CN115906704A CN202211584079.8A CN202211584079A CN115906704A CN 115906704 A CN115906704 A CN 115906704A CN 202211584079 A CN202211584079 A CN 202211584079A CN 115906704 A CN115906704 A CN 115906704A
- Authority
- CN
- China
- Prior art keywords
- contact
- unit
- equation
- pressure
- stress
- 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.)
- Pending
Links
- 238000010168 coupling process Methods 0.000 title claims abstract description 25
- 230000008878 coupling Effects 0.000 title claims abstract description 24
- 238000005859 coupling reaction Methods 0.000 title claims abstract description 24
- 238000004364 calculation method Methods 0.000 title claims abstract description 23
- 238000000034 method Methods 0.000 claims abstract description 24
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 claims description 8
- 230000002457 bidirectional effect Effects 0.000 claims description 3
- 230000009471 action Effects 0.000 claims description 2
- 239000007787 solid Substances 0.000 abstract description 2
- 238000010586 diagram Methods 0.000 description 8
- 208000010392 Bone Fractures Diseases 0.000 description 7
- 239000011435 rock Substances 0.000 description 7
- 230000008859 change Effects 0.000 description 5
- 230000008569 process Effects 0.000 description 4
- 238000004088 simulation Methods 0.000 description 4
- 238000006073 displacement reaction Methods 0.000 description 3
- 230000000694 effects Effects 0.000 description 3
- 239000011159 matrix material Substances 0.000 description 2
- 230000007246 mechanism Effects 0.000 description 2
- 230000035699 permeability Effects 0.000 description 2
- 208000003044 Closed Fractures Diseases 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000002939 conjugate gradient method Methods 0.000 description 1
- 230000001808 coupling effect Effects 0.000 description 1
- 239000012530 fluid Substances 0.000 description 1
- 239000010438 granite Substances 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 239000011148 porous material Substances 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
Images
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
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
The invention discloses a method for calculating a single-joint stress-seepage coupling numerical value, which belongs to the field of fluid-solid coupling numerical value calculation and comprises the following steps: acquiring input parameters of single-joint stress-seepage coupling, and performing discretization processing on the input parameters to obtain a discretization unit, wherein the discretization unit comprises: a contact unit and a non-contact unit; establishing a contact equation based on the contact unit, and obtaining contact stress and deformation based on the contact equation; updating the opening degree of the untouched unit based on the deformation to obtain an updated result, performing contact judgment on the updated result, establishing a Reynolds equation of the untouched unit based on the judged result, obtaining the pressure and the flow velocity of the untouched unit based on the Reynolds equation, judging the convergence based on the pressure, and if the untouched unit does not converge, continuously establishing a contact equation until the contact equation converges. The invention can solve the technical problems of high requirement on grids and long solving time, has high calculation speed and simple grid division, and realizes the automatic processing of a computer.
Description
Technical Field
The invention belongs to the field of fluid-solid coupling numerical calculation, and particularly relates to a single-joint stress-seepage coupling numerical calculation method.
Background
The coupling effect exists between the rock mass seepage field and the stress field, and the displacement field and the stress field are influenced by seepage load and influence the distribution of seepage water heads in the rock mass, because numerous joints are distributed in the rock mass, the stress effect causes the change of the geometric characteristics of the joints, such as the opening or closing of the joint opening, the change of the joint roughness and the like, thereby influencing the water permeability of the joints to cause the great change of the permeability of the whole rock mass. The change of the seepage field can change the seepage force, thereby influencing the joint stress field, and the two influence each other. Joint seepage and stress coupling are important thermal issues of rock hydraulics.
The joint stress-seepage coupling mechanism is very complicated due to the rough and uneven surface of the rock joint. Many experts and scholars at home and abroad generalize the joint roughness into one or two variables to establish a stress-seepage coupling model. These models neglect the effect of the joints indicating a fine asperity distribution and opening distribution on the stress-vadose coupling mechanism, leading to unavoidable errors. In order to accurately capture the influence of joint roughness, many scholars establish a joint fine model by using numerical techniques such as CFD (computational fluid dynamics) and discrete unit method and simulate the joint stress-seepage coupling process. However, the numerical simulation method needs a lot of time to establish the joint discrete grid model, and the numerical solution has high requirements on the grid, has long solution time, and often encounters the problem of non-convergence. Therefore, CFD simulation is cumbersome to apply.
Disclosure of Invention
The invention aims to provide a single-joint stress-seepage coupling numerical calculation method to solve the problems of high grid requirement and long solving time in the prior art.
In order to achieve the purpose, the invention provides a single-joint stress-seepage coupling numerical calculation method, which comprises the following steps of:
acquiring input parameters of single-joint stress-seepage coupling, carrying out discretization processing on the input parameters to obtain a discretization unit, wherein the discretization unit comprises: a contact unit and a non-contact unit;
establishing a contact equation based on the contact unit, and obtaining contact stress and deformation based on the contact equation;
updating the opening degree of the non-contact unit based on the deformation to obtain an updated result, performing contact judgment on the updated result, establishing a Reynolds equation of the non-contact unit based on the judged result, obtaining the pressure and the flow rate of the non-contact unit based on the Reynolds equation, judging the convergence based on the pressure, and if the pressure is not converged, continuing to establish a contact equation until the convergence.
Preferably, the discretizing the input parameter includes:
acquiring input parameters of single joint stress-seepage coupling, wherein the input parameters comprise: and discretizing the upper joint and the lower joint by adopting a rectangular unit to obtain a discretization unit.
Preferably, the contact equation comprises: an action force balance equation, a deformation coordination equation and an elasticity equation.
Preferably, before establishing the contact equation, the method further comprises:
carrying out initialization judgment on the discretization unit, judging whether the discretization unit is in contact or not, and if the discretization unit is in contact, determining that the discretization unit is a contact unit; if not, it is a non-contact unit.
Preferably, the process of obtaining the contact stress and the deformation amount comprises:
and solving the acting force balance equation by adopting a Gaussian iteration method based on the acting force balance equation to obtain the contact stress and the deformation of the convex body.
Preferably, the process of performing contact judgment on the update result includes:
if the opening degree of the updating result is smaller than 0, the touch unit is determined; and if the opening degree of the updating result is larger than 0, determining that the touch unit is not touched.
Preferably, the process of obtaining the pressure and flow rate of the non-contacting cell comprises:
based on the Reynolds equation and the cubic law, an equivalent weak integral form of the Reynolds equation is established by adopting a Galerkin method, a shape function is established by adopting rectangular units with bidirectional node interpolation, based on the shape function and the equivalent weak integral form, a discrete form of each rectangular unit and an equation set of water pressure of all unknown nodes are obtained, and based on the equation set, the pressure and the flow velocity of a non-contact unit are obtained.
Preferably, the process of judging the convergence based on the pressure includes:
if the opening and the pressure of the untouched unit meet the convergence relational expression, the convergence is determined; if the opening degree and the pressure of the non-contact means do not satisfy the convergence relational expression, convergence is not achieved.
The invention has the technical effects that:
the invention provides a single joint stress-seepage coupling numerical calculation method, which comprises the steps of obtaining input parameters of single joint stress-seepage coupling, and carrying out discretization processing on the input parameters to obtain a discretization unit, wherein the discretization unit comprises: a contact unit and a non-contact unit; establishing a contact equation based on the contact unit, and obtaining contact stress and deformation based on the contact equation; updating the opening degree of the non-contact unit based on the deformation to obtain an updated result, performing contact judgment on the updated result, establishing a Reynolds equation of the non-contact unit based on the judged result, obtaining the pressure and the flow rate of the non-contact unit based on the Reynolds equation, judging convergence based on the pressure, and if the pressure and the flow rate are not converged, continuing to establish a contact equation until the convergence.
The method can solve the technical problems of high grid requirement, long solving time and non-convergence in the prior art, and performs stress-seepage coupling simulation on a certain crack to obtain an opening distribution diagram, a flow velocity diagram and a pressure diagram under different normal stresses. The method has the advantages of high calculation speed, simple grid division and high division speed, realizes computer automatic processing, and considers the influence of fracture contact on seepage so as to ensure that the calculation is more reasonable.
Drawings
The accompanying drawings, which are incorporated in and constitute a part of this application, illustrate embodiments of the application and, together with the description, serve to explain the application and are not intended to limit the application. In the drawings:
FIG. 1 is a flow chart of a calculation method in an embodiment of the invention;
FIG. 2 is a schematic diagram of contact areas and non-contact areas of discrete slit units in an embodiment of the invention;
FIG. 3 is a schematic diagram of a split rock block fracture in an embodiment of the invention;
FIG. 4 is a schematic view of fracture morphology, elevation and opening measurements in an embodiment of the invention;
FIG. 5 is a schematic view of a relaxation profile without stress in an embodiment of the present invention;
FIG. 6 is a graph of flow velocity under various stresses in an embodiment of the present invention;
FIG. 7 is a graph of pressure distribution under various stresses in an embodiment of the present invention.
Detailed Description
It should be noted that the embodiments and features of the embodiments in the present application may be combined with each other without conflict. The present application will be described in detail below with reference to the embodiments with reference to the attached drawings.
It should be noted that the steps illustrated in the flowcharts of the figures may be performed in a computer system such as a set of computer-executable instructions and that, although a logical order is illustrated in the flowcharts, in some cases, the steps illustrated or described may be performed in an order different than presented herein.
Example one
As shown in fig. 1, the present embodiment provides a method for calculating a single-joint stress-seepage coupling value, including:
s1, inputting parameters: respectively inputting upper and lower joint coordinates (x, y) and height h (x, y)
S2, discretization treatment: the upper and lower joints are discretized with rectangular cells, each cell being denoted as (i, j), as shown in fig. 2.
S3, initialization judgment: judging whether each unit is contacted, and recording the contacted units as (i, j) m And the untouched cell is marked as (i, j) n
S4. Establishing contact area (i, j) m The acting force balance equation, the deformation coordination equation and the elastic equation are solved to obtain the force and the displacement of the contact unit.
The equilibrium equation of forces is as follows:
wherein F (i, j) m For stress distribution, P (i, j) n Is pore water pressure distribution, R 1 (i, j) is the upper radius of the pitch, R 2 (i, j) is the lower joint radius, see 2,R for the equivalent radius of the protrusion, R = (R) 1 (i,j)+R 2 (i, j))/2, Δ u is the deformation, v is the flow rate.
Solving the equilibrium equation set of forces by adopting a Gaussian iteration method to obtain the stress distribution F (i, j) of the convex body m And a deformation amount Δ u.
S5, according to the deformation value delta u obtained in the previous step, the opening b (i, j) of the non-contact area is calculated by using an opening calculation formula n And (6) updating.
The opening calculation formula is as follows:
b(i,j) n ＝b(i,j) n -Δu k (2)
s6, judging whether each unit is contacted according to the calculated result: if b (i, j) n < 0, then the cell is in contact, noted (i, j) m The opening degree is 0; if b (i, j) n > 0, the cell is not in contact,the contactless cell is denoted by (i, j) n
S7, establishing a non-contact unit (i, j) n Solving the discretized integral conduction equation of the Reynolds equation to obtain the pressure and the flow velocity of the non-contact area.
Considering that the water flow satisfies the Reynolds equation and the cubic law:
the equivalent weak integral form of the Reynolds equation is established by adopting a Galerkin method:
the rectangular unit adopting bidirectional node interpolation is adopted to establish a shape function:
in formula (5), the expression (6) brings about a discrete form of each rectangular unit,
K e P e ＝F e (7)
adding the coefficients of the corresponding numbers in the unit matrix to form a total conduction matrix in the formula (5), and obtaining an equation system for solving the water pressure of all unknown nodes:
KP＝F (11)
wherein [ K ]]Of (2) element(s)
Term constant->
Solving the formula (11) by adopting a pretreatment conjugate gradient method to obtain the water pressure on the untouched area unit
S8, solving the pressure value P (i, j) of the non-contact unit m
S9, judging whether convergence occurs: if b k -b k+1 ||<1e-4，||P k -P k+1 ||<1e-4, then convergence and calculation is completed, wherein b k And b k+1 Respectively representing the opening degree of k and k +1 iteration steps of the untouched unit, P k And P k+1 Respectively representing the pressure of k and k +1 iteration steps of the untouched unit;
if not, F (i, j) of the contact unit and the non-contact unit is updated m And P (i, j) n And establishing an acting force balance equation, a deformation coordination equation and an elastic equation of the contact unit again and solving to obtain the force and the deformation magnitude of the contact unit, and repeating S5-S8 until convergence.
And S10, finishing the calculation.
The specific implementation case is as follows:
(1) firstly, carrying out Brazilian splitting on a granite test piece with the thickness of 150mm multiplied by 150mm to obtain a fresh and complete crack, as shown in figure 3;
(2) measuring fracture appearance data, and measuring initial opening and contact radius. Firstly, sticking mark points (as shown in figure 4) on four side surfaces of an upper crack test piece and a lower crack test piece, and respectively scanning the rough surfaces of the upper crack and the lower crack and the contours under the closed crack by using a three-dimensional optical scanning system to obtain point cloud data; then, the crack is led into Geomaic studio software, and the upper and lower crack mark points are matched by using the splicing function of the software to simulate the crack closed state (figure 4); and finally, taking out the upper and lower fracture surfaces, calculating to obtain upper and lower profile elevations, subtracting the upper and lower joint elevations to obtain initial opening distribution data of the closed fracture, and if the difference between the upper joint and the lower joint at a certain point is less than 0, indicating that the contact is carried out, determining that the opening at the certain point is 0. And (4) above or below the contact radius, reducing the minimum elevation of each point to obtain a contact radius value.
(3) Input parameters and discretization processing: respectively inputting upper and lower joint coordinates (x, y), an elevation h (x, y), an initial opening b (i, j) and a contact radius R 1 (x, y) and R 2 (x, y). The upper and lower joints are discretized with rectangular cells, each cell being denoted as (i, j), as shown in fig. 2.
(4) Initialization judgment: judging whether each unit is contacted, and recording the contacted units as (i, j) m And the untouched cell is marked as (i, j) n
(5) Establishing contact areas (i, j) m And solving the acting force balance equation, the deformation coordination equation and the elastic equation to obtain the force and the displacement of the contact unit. And meanwhile, updating a new opening value and judging whether the contact is made.
(6) Establishing a contactless cell (i, j) n Solving the integral conduction equation of the discretization of the Reynolds equation to obtain the pressure and the flow velocity of the non-contact area. And solve for the pressure value of the non-contacted cell.
(7) Judging whether convergence occurs: if b k -b k+1 ||<1e-4，||P k -P k+1 ||<1e-4, converging and finishing calculation; if not, F (i, j) of the contact unit and the non-contact unit is updated m And P (i, j) n And (5) establishing an acting force balance equation, a deformation coordination equation and an elastic equation of the contact unit again and solving to obtain the force and the deformation of the contact unit, and repeating the steps (5) to (7) until convergence.
(8) And finishing the calculation.
(9) Different pressures were set to recalculate deformation and seepage. The normal stresses used for this calculation are respectively: 0. 1, 2, 3, 4 and 5MPa. Fig. 5, 6 and 7 show the opening, flow rate and pressure distribution, respectively, for different normal stresses.
The beneficial effects of the embodiment are that:
the embodiment can solve the technical problems of high grid requirement, long solving time and non-convergence in the prior art, and performs stress-seepage coupling simulation on a certain fracture to obtain an open distribution diagram, a flow velocity diagram and a pressure diagram under different normal stresses. The embodiment has the advantages of high calculation speed, simple grid division and high division speed, realizes computer automated processing, and considers the influence of fracture contact on seepage so as to make calculation more reasonable.
The above description is only for the preferred embodiment of the present application, but the scope of the present application is not limited thereto, and any changes or substitutions that can be easily conceived by those skilled in the art within the technical scope of the present application should be covered within the scope of the present application. Therefore, the protection scope of the present application shall be subject to the protection scope of the claims.
Claims (8)
1. A single-joint stress-seepage coupling numerical calculation method is characterized by comprising the following steps:
acquiring input parameters of single-joint stress-seepage coupling, and discretizing the input parameters to obtain a discretization unit, wherein the discretization unit comprises: a contact unit and a non-contact unit;
establishing a contact equation based on the contact unit, and obtaining contact stress and deformation based on the contact equation;
updating the opening degree of the non-contact unit based on the deformation to obtain an updated result, performing contact judgment on the updated result, establishing a Reynolds equation of the non-contact unit based on the judged result, obtaining the pressure and the flow rate of the non-contact unit based on the Reynolds equation, judging the convergence based on the pressure, and if the pressure is not converged, continuing to establish a contact equation until the convergence.
2. The single-joint stress-seepage coupling numerical calculation method of claim 1, wherein discretizing the input parameters comprises:
obtaining input parameters of single joint stress-seepage coupling, wherein the input parameters comprise: and discretizing the upper joint and the lower joint by adopting a rectangular unit to obtain a discretization unit.
3. The method of claim 1, wherein the contact equation comprises: an action force balance equation, a deformation coordination equation and an elasticity equation.
4. The method of claim 1, wherein the establishing the contact equation further comprises:
carrying out initialization judgment on the discretization unit, judging whether the discretization unit is in contact or not, and if the discretization unit is in contact, determining that the discretization unit is a contact unit; if not, it is a non-contact unit.
5. The method of claim 3, wherein the obtaining of the contact stress and the deformation comprises:
and solving the acting force balance equation by adopting a Gaussian iteration method based on the acting force balance equation to obtain the contact stress and the deformation of the convex body.
6. The method of claim 1, wherein the step of determining whether the updated result is in contact with the target object comprises:
if the opening degree of the updating result is smaller than 0, the touch unit is determined; and if the opening degree of the updating result is larger than 0, determining that the touch unit is not touched.
7. The method of claim 1, wherein obtaining the pressure and flow rate of the non-contacting element comprises:
based on the Reynolds equation and the cubic law, an equivalent weak integral form of the Reynolds equation is established by adopting a Galerkin method, a shape function is established by adopting rectangular units with bidirectional node interpolation, an equation set of a discrete form of each rectangular unit and water pressure of all unknown nodes is obtained based on the shape function and the equivalent weak integral form, and the pressure and the flow rate of the units which are not contacted are obtained based on the equation set.
8. The method of claim 1, wherein the step of determining convergence based on the pressure comprises:
if the opening degree and the pressure of the non-contact unit meet the convergence relational expression, convergence is achieved; if the opening and pressure of the non-contact means do not satisfy the convergence relational expression, the non-contact means does not converge.
Priority Applications (2)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202211584079.8A CN115906704A (en) | 2022-12-09 | 2022-12-09 | Single-joint stress-seepage coupling numerical calculation method |
| ZA2023/00724A ZA202300724B (en) | 2022-12-09 | 2023-01-17 | Single joint stress-seepage coupling numerical calculation method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202211584079.8A CN115906704A (en) | 2022-12-09 | 2022-12-09 | Single-joint stress-seepage coupling numerical calculation method |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN115906704A true CN115906704A (en) | 2023-04-04 |
Family
ID=85768826
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202211584079.8A Pending CN115906704A (en) | 2022-12-09 | 2022-12-09 | Single-joint stress-seepage coupling numerical calculation method |
Country Status (2)
| Country | Link |
|---|---|
| CN (1) | CN115906704A (en) |
| ZA (1) | ZA202300724B (en) |
-
2022
- 2022-12-09 CN CN202211584079.8A patent/CN115906704A/en active Pending
-
2023
- 2023-01-17 ZA ZA2023/00724A patent/ZA202300724B/en unknown
Also Published As
| Publication number | Publication date |
|---|---|
| ZA202300724B (en) | 2023-03-29 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Ooi et al. | A hybrid finite element-scaled boundary finite element method for crack propagation modelling | |
| CN106529038B (en) | A method of it is bolted from micro--grand yardstick model identification in conjunction with the tangential damping characteristic in portion | |
| CN104077440B (en) | A kind of method of determination faying face contact area and rigidity based on surface fitting | |
| CN107341316B (en) | Structural shape-topology combined optimization method under design related pressure load effect | |
| CN103335814B (en) | Correction method for inclination angle measurement error data of experimental model in wind tunnel | |
| CN112417573A (en) | Multi-objective optimization method for shield tunneling underneath existing tunnel construction based on GA-LSSVM and NSGA-II | |
| CN103956109B (en) | A kind of tunnel model test method that liner structure joint stiffness is variable | |
| CN110765695B (en) | Simulation calculation method for obtaining crack propagation path of concrete gravity dam based on high-order finite element method | |
| CN113158315B (en) | Rock-soil body parameter three-dimensional non-stationary conditional random field modeling method based on static sounding data | |
| CN109992885A (en) | A kind of casting limited fatigue life member prediction technique and system | |
| CN115146565B (en) | Hyperbolic arch dam foundation seepage evaluation method and system based on cohesion algorithm | |
| CN104750892A (en) | Three-dimensional modeling method for thickness-variable curved-surface part inner shape surface | |
| CN108345745A (en) | A kind of liquid hydrogen storage tank low temperature prestressing force Wetted modes analysis method based on fluid structurecoupling | |
| Zhang et al. | Voronoi cell finite element model to simulate crack propagation in porous materials | |
| CN115906704A (en) | Single-joint stress-seepage coupling numerical calculation method | |
| CN105045958A (en) | Implementation system and method of GPS (Global Positioning System) elevation fitting on the basis of BP (Back Propagation) neural network | |
| CN104133925A (en) | Nonlinear analysis method applicable to cockpit canopy structure static intensity | |
| CN114692441B (en) | Loess landslide stability prediction method, electronic equipment and storage medium | |
| Wowk et al. | Influence of p-method finite element parameters on predictions of crack front geometry | |
| CN108563613A (en) | A kind of Mathematical Modeling Methods using computer hyperspace | |
| CN108197397B (en) | Optimization design method for dynamic performance of fastening joint surface of aircraft engine | |
| Subcommittee | Committee report: Defining model calibration | |
| CN110532669A (en) | A method of it is modeled for Machine Joint Surfaces contact stiffness | |
| CN109992853A (en) | A kind of residual stress field numerical value method for reconstructing of surface peening metal parts | |
| CN117421939B (en) | Shale oil fracture system simulation agent method based on track piecewise linearization |
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 |