Numerical calculation method for sand dolomite seepage failure process

Abstract

The invention discloses a numerical calculation method for a sand dolomite seepage failure process. The method comprises the following steps: generating a sand dolomite model sample; generating a numerical simulation model based on fluid-solid coupling; the movement direction of the fluid is set, and the seepage failure process of the sanded dolomite under different water flow conditions is achieved; generating model diagrams of different time nodes in the test process, and connecting the images together to reproduce the migration process of the filling particles; writing a function to count the mass loss condition of the sandy dolomite model sample; recording the porosity change in the particle migration process, and drawing a porosity change curve chart; and acquiring an average value of the fluid grid permeability coefficients by using the stored model file, and drawing a permeability coefficient change curve. According to the method, numerical simulation and a fluid-solid coupling model are combined, the seepage failure process of the sanded dolomite under the action of different water flows can be accurately analyzed, comprehensive analysis of the seepage failure process of the sanded dolomite is achieved, and the coupling relation between the water flows, particle migration and the like is revealed.

Description

Technical Field

[0001] The invention belongs to the technical field of engineering geology and geotechnical engineering, and particularly relates to a numerical calculation method for the seepage failure process of sandy dolomite. Background Art

[0002] The seepage damage process of sandy dolomite is a complex water-rock interaction phenomenon, and the sudden water (mud) problem caused by it is one of the major geological disasters in tunnel construction in karst areas. Existing studies on the seepage damage of sandy dolomite are mostly conducted in indoor model tests. Because the water pressure is relatively high, steel is usually used to make permeameters to ensure the safety of the test. The indoor model test method cannot directly observe the movement process of the filling fine particles. Moreover, existing studies mostly focus on the analysis of macroscopic geological characteristics or empirical formula predictions, and lack a comprehensive study of the microscopic dissolution mechanism and seepage damage characteristics of sandy dolomite. With the rapid development of numerical simulation technology, numerical models based on fluid-solid coupling can provide assistance for studying the seepage damage process of sandy dolomite. Summary of the invention

[0003] The present invention aims to analyze the seepage damage process of sandy dolomite by constructing a fluid-solid coupling model based on particle flow simulation, which can be used to predict and control the seepage damage risk of sandy dolomite in tunnel construction.

[0004] The technical solution of the present invention is as follows: A numerical calculation method for the seepage failure process of sandy dolomite comprises the following steps: Step 1, generate a structural model: program in the PFC software to establish a rock container model, fill the container model with the particle model generated by the PFC software, and obtain a sandy dolomite model sample; Step 2, generate a fluid-solid coupling model and divide the grid: Use the CFD module of the PFC software to simulate the fluid-solid coupling between particles and fluids, and use the Python language for secondary development to divide the grid and add the grid to the rock container model to form a numerical simulation model based on fluid-solid coupling; Step 3, simulating the seepage destruction process: setting the movement direction of the fluid, and controlling the magnitude of the hydraulic gradient value to achieve the seepage destruction process of sandy dolomite under different water flow conditions; Step 4, set the interval time to automatically save the model file, export the model file, generate model images at different time nodes during the test, connect these images together, and reproduce the migration process of filling particles when penetration damage occurs; Step 5, write a function to calculate the mass loss of the sandy dolomite model sample; Step 6, using the porosity change function in the PFC software, record the porosity change during the particle migration process; export the porosity change recorded in the history, and draw a porosity change curve; Step 7: Use the saved model file to obtain the average value of the fluid grid permeability coefficient and draw the permeability coefficient change curve.

[0005] Furthermore, in step 1, the FISH programming language is used to program and establish a rock container model, and the wall of the rock container is set to be a smooth rigid body and impermeable, and the rigidity of the wall is greater than the rigidity of the particles.

[0006] The sandy dolomite particles are regarded as complete spheres, and the PFC software is used to generate two particle models with different particle sizes, which are used as skeleton particles and filling particles respectively.

[0007] Furthermore, in step 2, after the particles are generated in the PFC software, they are given motion and mechanical properties, and then fluid calculations are performed in the CFD module.

[0008] The motion and mechanical properties given to the sandy dolomite rock particles include particle mass, inertia, contact mechanics model (such as contact stiffness, friction coefficient, bonding strength, etc.), external forces (such as gravity, confining pressure, fluid force), etc., to ensure that the particles can be subjected to force and move in the seepage field. Each particle is regarded as an independent rigid sphere, governed by Newton's second law, and has the ability to move freely and contact each other. The fluid calculation steps in the CFD module belong to the existing technology and will not be repeated.

[0009] Furthermore, in step 3, an initial water pressure is given to the inlet boundary at 0S, and the outlet boundary water pressure is set to a fixed value of 0. Thereafter, the inlet boundary water pressure is increased at regular intervals to realize the permeability destruction process of sandy dolomite under different water flow conditions with different hydraulic gradient values.

[0010] Furthermore, in step 5, a function is written to measure the total mass of the current model in real time, and the current total mass loss is obtained by comparing the total mass of the current model with the total mass of the initial model; the mass loss data is imported into an Excel spreadsheet, and a curve of the current mass loss is drawn based on the data.

[0011] Furthermore, the statistical data of the filling particle loss at different time points are derived, the mass loss rate of the filling particles is calculated, and the mass loss rate variation curve is drawn.

[0012] Furthermore, in step 6, the porosity change of the historical records is exported in an Excel table, and a porosity change curve is drawn according to the data.

[0013] Beneficial effects: The present invention combines numerical simulation and fluid-solid coupling model to accurately analyze the seepage damage process of sandy dolomite under different water flows. It realizes a comprehensive analysis of the seepage damage process of sandy dolomite, reveals the coupling relationship between water flow and particle migration, etc., can help predict and control the seepage damage risk of sandy dolomite in tunnel construction, and provide a scientific risk assessment and optimization design basis for tunnel construction. BRIEF DESCRIPTION OF THE DRAWINGS

[0014] Figure 1 The structure of the rock mass container model is shown; Figure 2 The particle model structure generated by the PFC software is shown; Figure 3 The structure of the sandy dolomite model is shown; Figure 4 The migration process of some particles during the penetration damage process is shown in Figure 1; Figure 4 The model diagram when a is 0S, Figure 4 The model diagram when b is 0.1S, Figure 4 The model diagram when c is 0.2S, Figure 4 The model diagram when d is 0.3S, Figure 4 The model diagram when e is 0.4S, Figure 4 Model diagram when f is 0.5; Figure 5 The second part shows the migration process of some particles during the penetration damage process; Figure 5 The model diagram when a is 0.6S, Figure 5 The model diagram when b is 0.7S, Figure 5 The model diagram when c is 0.8S, Figure 5 Model diagram when d is 1.0S; Figure 6 The migration process of a single packing particle is shown; Figure 7 The model mass loss variation curve is shown; Figure 8 The curve of the model mass loss rate variation is shown; Fig. 9 The curve of the model permeability variation is shown; Fig.10 The structure of the measuring sphere is shown; Fig.11 The model porosity variation curve is shown. DETAILED DESCRIPTION

[0015] The present invention will be further described below in conjunction with the accompanying drawings and specific embodiments.

[0016] Example

[0017] A numerical calculation method for the seepage failure process of sandy dolomite comprises the following steps: Step 1: Use the FISH programming language to create a cylindrical rock container model, such as Figure 1 As shown. In this embodiment, after repeated calculations, the size of the rock container model is 100mm×400mm and the diameter is 200mm. In order to simulate the actual situation, the model is expanded by 50mm at both ends in the y-axis direction. The confining pressure is set to 500kpa to simulate the burial depth of 50m underground, and the wall on the right side of the rock container model is removed to simulate the loss of particles. The wall of the rock container model is designed to be a smooth and rigid impermeable material with a stiffness slightly higher than that of the particles.

[0018] like Figure 1 As shown in the figure, the particle model is used to act as the soil particles in reality, and the rock container model is like a steel cylinder, which restricts the soil particles from scattering. In the numerical model, the wall on the left side of the rock container model can be set to allow water to enter. A high-pressure water (or a certain head) is applied abstractly at the left end of the rock container model to allow the water to flow through the internal soil particles. Because the wall on the right side of the rock container model is removed, the particles can be carried away with the water flow.

[0019] In this example, the particles of sandy dolomite are regarded as perfect spheres, and a random uniform distribution model containing two particles of different particle sizes is constructed using PFC software, wherein the diameters of the skeleton particles and the filling particles in the rock container model are set to 40 mm and 20 mm, respectively. The mass ratio of the two particles is set to 1:1, and the initial porosity is set to 0.35. Under these parameter conditions, when the unbalanced force between the particles is close to 0, the model generates a total of 8303 particles. The particle model generated by PFC software is as follows: Figure 2 Its main parameters are shown in Table 1.

[0020] Table 1 Main parameters of particle model

[0021] Step 2, using the CFD module of the PFC software to simulate the fluid-solid coupling between particles and fluid. After the particle unit is generated in the PFC software, its motion and mechanical properties are given, and then the fluid calculation is performed in the CFD. The specific calculation process is the existing technology. The calculation setting parameters for the simulated particles of sandy dolomite in this embodiment are shown in Tables 1 and 2.

[0022] Use Python language for secondary development to divide the grid. Since the particle size of the skeleton particles is 40 mm, the size of a single grid is set to a 40 mm cube, and a total of 300 grids are generated, such as Figure 3As shown in the figure, the grid is added to the rock container model, and finally a numerical simulation model based on fluid-solid coupling is formed, namely the sandy dolomite model.

[0023] Step 3, set the fluid movement direction to be from left to right along the y-axis, and control the hydraulic gradient value (pressure difference) during simulation to achieve fluid movement. After repeated tests, the initial water pressure at the inlet boundary is given at 0S, and the hydraulic gradient value (i.e. Figure 4 , 5 i) is set to 0.1, and the hydraulic gradient value increases by 0.1 every 0.1S, while the outlet boundary water pressure is set to a fixed value of 0. After the model is initialized, the CFD time step is set to 1.0e-5S, and the DEM time step is set to 1.0e-6S. The fluid model and calculation time step parameter settings are shown in Table 2.

[0024] Table 2 Fluid model and calculation time step parameter settings

[0025] DEM calculation time step: The discrete element (particle) updates the motion equation in each iteration step, and a smaller step size is required to ensure the stability of the calculation. In each tiny time step, the program will calculate the external forces and torques acting on the particles based on the updated motion equations, and then update the particle's velocity, position and angular velocity, repeating the cycle to simulate the dynamic evolution of the particles.

[0026] CFD calculation time step: The time increment for iterative update of the fluid force field, which is usually larger than the DEM step size, is used to calculate the drag force of the fluid on the particles, pressure distribution, etc. When the fluid calculation is completed, it is coupled with the DEM to exchange information to obtain the new position and porosity of the particles, etc., forming a cyclic iteration of fluid-solid coupling.

[0027] When the penetration damage occurs and the loss of filling particles is basically stable, the pressurization is stopped and the final model running time is 1S. Figure 5 shown.

[0028] Whether the loss of filling particles is stable can be judged by the mass loss rate. When the mass loss rate curve remains stable, it can be judged that it is basically stable.

[0029] Step 4: During the penetration damage process, save the model file every 0.01S, and then generate the model diagrams at different time points during the test according to the exported file. By connecting these images together, the migration process of the filling particles when the penetration damage occurs can be reproduced, such as Figure 4 a, b, c, d, e, f, and Figure 5 In this simulation, the total simulation time is 1S, with a total of 1 million steps.

[0030] From the migration process of the filling particles when the seepage damage occurs, it can be seen that at the beginning of the test, the filling particles gradually migrated with the water flow toward the right outlet, and the filling particles at the right outlet flowed out earlier than the filling particles on the left. Some of the filling particles will be blocked by the skeleton particles and stuck in the pores and no longer move, while the other part of the particles will completely flow out with the water flow. From a macroscopic perspective, the water flows straight from the inlet (left) boundary to the outlet (right) boundary. But from a microscopic perspective, the water flow cannot pass through the sandy dolomite rock mass in a straight line like that, but flows back and forth in the pores of the skeleton particles, flowing along the high water pressure to the part with lower water pressure, and carrying out fine particles.

[0031] The movement trajectory of a single filling particle, such as Figure 6 As shown in the figure, it is found that the movement of a single filling particle in the sandy dolomite rock mass under seepage conditions is highly random. In this experiment, the particles are affected by four forces: mutual adhesion between particles, gravity, surrounding rock pressure, and water flow drag force. The combined force of these forces determines the randomness of particle movement. Because the porosity of the sample is uneven, the water flow velocity and drag force are different at different locations, and the interaction force between each particle and the surrounding rock pressure are also different, resulting in the sum of the forces on a single particle not being exactly the same. In general, the movement trend of a single particle is to migrate from the inlet boundary of the model along the water flow direction to the outlet, but its specific migration trajectory cannot be predicted.

[0032] Step 5, write a function to count the mass loss of the sanded dolomite model specimen. When modeling, the wall and the calculation area have been set to the same length, so that when the filling particles flow out of the wall with the water, they will be automatically deleted. By comparing the total mass of the current model with the total mass of the initial model, the current total mass loss can be obtained. After importing these data into an Excel spreadsheet and converting the units, the current mass loss change curve can be drawn. In this way, the loss of particles during the seepage process can be recorded in real time and the rules can be analyzed. The mass loss change curve of the sanded dolomite model is shown below. Figure 7 Table 3 shows the total loss mass and percentage at each stage.

[0033] To write a function in the PFC software using the FISH programming language to track mass loss in real time, follow these steps: 1. Define a variable to record the initial quality When the model is first generated, all the filling particles are traversed and the mass is added to get initial_mass.

[0034] 2. Write a FISH function named CalcMassLoss This function is called automatically after each time step: Initialize a temporary variable current_mass = 0; Traverse all existing filling particles in the model and add their individual masses to current_mass; Calculate mass_loss = initial_mass - current_mass; Store or export mass_loss to log, txt file, Excel, or save in PFC history.

[0035] 3. Determine which filler particles are taken out of the model When the filling particles exceed the right boundary, they are automatically deleted, which means that the filling particles are lost from the outlet; the total mass after deleting the filling particles will be reduced accordingly.

[0036] 4. Through the above function, the mass of the remaining filling particles in the model at the current moment can be obtained, thereby obtaining the mass loss amount.

[0037] Table 3 Mass loss and percentage at each stage

[0038] From Table 3 and Figure 7 According to the data, the migration of filling particles in the sandy dolomite mass is instantaneous and nonlinear. Some filling particles are migrated out along the pores of the skeleton particles under the seepage of water, while the other part is blocked by the skeleton particles. This process can be roughly divided into three stages: slow seepage stage (0S-0.2S), fast seepage stage (0.2S-0.8S) and stable seepage stage (0.8S-1.0S).

[0039] (1) Slow seepage stage (0S-0.2S): When the initial water pressure is applied to the sample at 0S (i=0 at this time), the model is under the action of the mutual adhesion of particles, gravity, surrounding rock pressure and water flow drag force, and the sum of other forces is less than the water flow drag force. Therefore, the filling particles at the outlet on the right side of the model are taken out of the seepage channel, resulting in a rapid rise of the mass loss curve at the beginning. With the continuous migration of the filling particles at the outlet of the model, the filling particles inside also gradually begin to migrate slowly and flow out. At this time, the curve is relatively flat. The total mass loss in this stage is 156.2g, accounting for about 5% of the mass loss in the entire seepage stage.

[0040] (2) Rapid seepage stage (0.2S-0.8S): As the filling particles continue to migrate and flow out during the slow seepage stage, the gaps in the model become larger and the permeability is further enhanced. When t=0.2S (i=0.3), the critical hydraulic gradient is reached, the mass loss of the particles suddenly increases, and the seepage changes suddenly. As the hydraulic gradient continues to increase, the filling particles are quickly and continuously carried out by the water flow, forming a seepage channel. The total mass loss in this stage is 2618.5g, accounting for about 84.2% of the mass loss in the entire seepage stage.

[0041] (3) Stable seepage stage (0.8S-1.0S): When t = 0.8S (i = 0.9), the filling particles basically no longer migrate, and the porosity and permeability coefficient basically no longer change. The total mass loss in this stage is 339.2g, accounting for about 10.8% of the mass loss in the entire seepage stage.

[0042] According to the saved model file, the statistical data of the filling particle loss at different time nodes can be exported, and the mass loss rate of the model particles can be calculated based on this. As shown in Table 4, it includes the total mass of the particles M i , total cumulative loss ∑△m i 、Mass loss per 0.1S△m i And the mass loss rate △m i / △t. Finally, the Figure 8 The mass loss rate variation curve is shown.

[0043] The particle loss refers only to the loss of 20mm filling particles, because the loss of skeleton particles is very small. It is mainly due to the appearance of seepage channels caused by the loss of filling particles, which further increases the permeability coefficient and leads to permeability damage. Therefore, the loss of filling particles is equivalent to the loss of the model, and the loss of skeleton particles is not considered.

[0044] Table 4 Model mass loss rate

[0045] Depend on Figure 8 It can be seen that the mass loss rate of the model shows an overall trend of increasing first and then decreasing during the entire seepage process, with the mass loss rate being the largest at 0.4S, with a maximum value of about 5029g / s. At 0.7S, a peak value of about 4672g / s also appears.

[0046] Step 6, recording the porosity change during the particle migration process according to the porosity change function in the PFC software; The porosity change function records the particle migration process, which is achieved by the following steps: 1. Define the measurement area Three measuring balls of equal size are constructed along the y-axis direction by using the FISH function. They are named measuring ball No. 1, measuring ball No. 2 and measuring ball No. 3. The radius of each ball is 80 mm and the center of each ball is located on the center axis of the y-axis of the model. Fig.10 As shown, measuring balls No. 1, No. 2 and No. 3 measure the porosity of the left, middle and right sides of the model respectively.

[0047] 2. Set the porosity of the measurement area The pore volume of the area can be obtained by subtracting the total volume of particles in the area where the measuring ball is located from the volume of the measuring ball, and the porosity can be obtained by dividing the volume of the measuring ball.

[0048] 3. History Record Function Set the porosity of the measurement area where the measuring ball is located as a history item, and set the PFC software to automatically save the porosity values ​​of the area at intervals. Then export these data at the end or midway of the simulation to obtain discrete point data of porosity changes over time.

[0049] In order to clearly show the change process of porosity, the porosity change records were exported in Excel and plotted as follows: Fig.11 The green curve represents the monitoring result of measuring ball No. 1 in the left area of ​​the model, the red curve represents the monitoring result of measuring ball No. 2 in the middle area of ​​the model, and the blue curve represents the monitoring result of measuring ball No. 3 in the right area of ​​the model.

[0050] The porosity change data of the exported historical records needs to be processed to uniform decimal places before drawing the curve.

[0051] In general, the porosity represented by the three curves first increases slowly, then increases sharply, and finally tends to be flat.

[0052] The measuring ball No. 3 represented by the blue curve is located on the right side of the model, and the filling particles are lost first. Therefore, the porosity change represented by the blue curve rises sharply and reaches the highest amplitude in the rapid seepage stage.

[0053] The measuring ball No. 2 represented by the red curve is located in the middle of the model. The red curve drops several times in the rapid seepage stage. This is because the filling particles on the left side of the model flow to the right with the water flow, and are blocked by the skeleton particles in the middle and have not moved yet, causing the porosity curve to drop briefly. As the hydraulic gradient increases, the filling particles flow out from the gaps in the skeleton particles to the right outlet with the water flow, and the porosity gradually increases.

[0054] The measuring ball No. 1 represented by the green curve is located on the left side of the model. The change trend of the green curve is similar to that of the red curve, but in the rapid seepage stage (t=0.2S-0.8S), as the hydraulic gradient increases, the filling particles on the left side migrate and flow out rapidly, causing the porosity to rise sharply.

[0055] from Fig.11 It can be seen that with the passage of time, the final values ​​of the three porosity change curves will approach the same stable value.

[0056] Step 7: Use the saved model file to obtain the average value of the fluid grid permeability coefficient and draw the permeability coefficient change curve.

[0057] The model file saves the physical quantities in the fluid motion equations such as velocity field, pressure field, porosity, etc., which can then be calculated according to Darcy's formula:

[0058] Where v is the local velocity, μ is the fluid viscosity, is the pressure gradient and k is the permeability coefficient.

[0059] Calculate the k in each fluid grid at each moment, take the average value to get the overall k at that moment, and then plot the k at different moments into a curve to get the permeability coefficient change curve.

[0060] Under the action of pressure difference, the filling particles are continuously migrated out, resulting in a gradual increase in the porosity of the model. At the same time, this process will cause more filling particles to be taken out of the seepage channel, and the permeability of the model will be further increased. Using the saved model file, obtain the average value of the fluid grid permeability coefficient every 0.01S and plot it as a continuous curve, see Fig. 9 .

[0061] according to Fig. 9 As shown in the figure, the development process of permeability coefficient can be roughly divided into three stages: slow seepage stage, fast seepage stage and stable seepage stage, which experienced a trend from slow to fast and then to stable, which is roughly the same as the change of porosity. During the entire seepage process, the permeability coefficient increased from the initial 9.9e-06m / s to the final 193.2e-06m / s.

Claims

1. A numerical calculation method for the seepage failure process of sandy dolomite, characterized in that: The following steps are involved: Step 1, generate a structural model: program in the PFC software to establish a rock container model, fill the container model with the particle model generated by the PFC software, and obtain a sandy dolomite model sample; Step 2, generate a fluid-solid coupling model and divide the grid: Use the CFD module of the PFC software to simulate the fluid-solid coupling between particles and fluids, and use the Python language for secondary development to divide the grid and add the grid to the rock container model to form a numerical simulation model based on fluid-solid coupling; Step 3, simulating the seepage destruction process: setting the movement direction of the fluid, and controlling the magnitude of the hydraulic gradient value to achieve the seepage destruction process of sandy dolomite under different water flow conditions; Step 4, set the interval time to automatically save the model file, export the model file, generate model images at different time nodes during the test, connect these images together, and reproduce the migration process of filling particles when penetration damage occurs; Step 5, write a function to calculate the mass loss of the sandy dolomite model sample; Step 6, using the porosity change function in the PFC software, record the porosity change during the particle migration process; export the porosity change recorded in the history, and draw a porosity change curve; Step 7: Use the saved model file to obtain the average value of the fluid grid permeability coefficient and draw the permeability coefficient change curve.

2. The numerical calculation method for the seepage failure process of sandy dolomite according to claim 1 is characterized in that: In step 1, the FISH programming language is used to program and establish a rock container model, and the wall of the rock container is set to be a smooth rigid body and impermeable, and the rigidity of the wall is greater than that of the particles; The sandy dolomite particles are regarded as complete spheres, and the PFC software is used to generate two particle models with different particle sizes, which are used as skeleton particles and filling particles respectively.

3. The numerical calculation method for the seepage failure process of sandy dolomite according to claim 2 is characterized in that: In step 2, after the particles are generated in the PFC software, they are given motion and mechanical properties, and then fluid calculations are performed in the CFD module.

4. The numerical calculation method for the seepage failure process of sandy dolomite according to claim 2 is characterized in that: In step 3, the initial water pressure at the inlet boundary is given at 0S, and the water pressure at the outlet boundary is set to a fixed value of 0. Thereafter, the water pressure at the inlet boundary is increased at regular intervals to realize the penetration destruction process of sandy dolomite under different water flow conditions with different hydraulic gradient values.

5. The numerical calculation method for the seepage failure process of sandy dolomite according to claim 2 is characterized in that: In step 5, a function is written to measure the total mass of the current model in real time. By comparing the total mass of the current model with the total mass of the initial model, the current total mass loss is obtained. The mass loss data is imported into an Excel spreadsheet, and a curve of the current mass loss is drawn based on the data.

6. The numerical calculation method for the seepage failure process of sandy dolomite according to claim 5 is characterized in that: In step 5, the statistical data of the filling particle loss at different time points are exported, the mass loss rate of the filling particles is calculated, and a mass loss rate variation curve is drawn.

7. The numerical calculation method for the seepage failure process of sandy dolomite according to claim 2 is characterized in that: In step 6, the porosity change of the historical records is exported in an Excel table, and a porosity change curve is drawn based on the data.

Images