A rigid wall servoing method based on discrete element particle motion

Abstract

The application discloses a rigid wall servo method based on discrete element particle motion, and first establishes a discrete element biaxial compression test model composed of a rigid wall and particles, and makes the wall contact stress tend to a target stress by controlling the displacement speed of the rigid wall; the displacement speed of the rigid wall is composed of a displacement speed initial value and a displacement speed correction value, wherein the displacement speed initial value is calculated based on the difference between the wall contact stress of a current time step and the target stress, and the displacement speed correction value is calculated based on the difference between the wall contact stress of the current time step and the wall contact stress of a previous time step; when a set calculation termination condition is met, the discrete element calculation is terminated. Compared with a traditional servo method, the servo precision of the method provided by the application is significantly improved and has little influence on the calculation efficiency.

Description

Technical Field

[0001] This invention relates to the field of discrete element method (DEM) technology in computational geotechnical engineering, and more particularly to a rigid wall servo method based on discrete element particle motion. Background Technology

[0002] Discrete Element Method (DEM) has been widely applied in geotechnical engineering, materials design, agriculture, and other fields. Its basic principle is to solve for the resultant force on each particle based on the interaction forces between particles, and to calculate the particle velocity, displacement, and other information at the next time step based on Newton's laws of motion. The motion of the particles reflects the macroscopic motion and deformation of the material. Because the objects simulated by the DEM do not need to satisfy the continuity assumption, it can naturally simulate the large deformation and failure processes of materials such as soil and rock.

[0003] In Discrete Element Methods (DEM), two methods are used to apply loads to materials: flexible membrane boundaries and rigid wall boundaries. Flexible membrane boundaries are formed by a series of cemented particles. External forces act on the membrane particles and the load is directly applied through the contact between the membrane particles and the internal particles. Because the membrane particles can rotate, the applied load distribution is more uniform. However, this method is more difficult to use and has a certain impact on computational efficiency. Rigid wall boundaries consist of one or more massless rigid walls. One rigid wall comes into contact with multiple particles. The load applied by the rigid wall is defined as the ratio of the sum of the normal contact forces perpendicular to the rigid wall to the area of ​​the rigid wall. This type of boundary uses a servo method to control the rigid wall to squeeze the particles at a specific speed, achieving indirect load application through the contact force between the rigid wall and the particles. This method is convenient to use and has a relatively low impact on computational efficiency, and is widely used in various simulations.

[0004] According to the force-displacement law of contact, the normal force acting on the wall depends on the amount of overlap between the wall and the particles. Therefore, the core of the rigid wall servo method lies in controlling the wall displacement velocity so that the overlap in the next time step is just enough to generate the required load. Due to the large number of particles in contact with the wall and the high computational cost of traversal, traditional servo methods, to ensure computational efficiency, assume that the particles do not move during the simulation and calculate the wall displacement velocity based on the difference between the stress acting on the wall at the current time step and the target stress. However, this method produces a large error when the particle movement speed is high. Therefore, it is necessary to design a new method that can consider particle motion while ensuring computational efficiency to improve servo accuracy. Summary of the Invention

[0005] Purpose of the invention: To address the shortcomings of existing technologies, this invention proposes a rigid wall servo method based on discrete element particle motion. By controlling the displacement velocity of the rigid wall, load can be applied in any situation. This solves the problem that traditional servo methods do not consider the influence of particle motion when calculating the rigid wall motion, and significantly improves servo accuracy while ensuring computational efficiency.

[0006] Technical solution: The rigid wall servoing method based on discrete element particle motion of the present invention includes the following steps:

[0007] Step (1) establishes a discrete element biaxial compression test model consisting of a rigid wall and particles; the process is as follows:

[0008] Step (1.1): Use mutually perpendicular rigid walls to form a region, create boundary conditions, and determine the width and height of the region enclosed by the rigid walls;

[0009] Step (1.2) involves generating overlapping particles within the rigid wall component area;

[0010] Step (1.3) sets a small contact stiffness value and a large damping value for the particle-particle contact;

[0011] Step (1.4) involves iterative calculations to cause the overlapping particles to bounce off under the action of contact force, while the discrete element biaxial compression test model is balanced under the action of damping. The calculation is terminated when the unbalanced force of the discrete element biaxial compression test model is less than the set value, thus obtaining the initial discrete element biaxial compression test model.

[0012] Step (2) sets initial parameters, including target load, contact parameters between rigid wall and particles, and calculation termination conditions; the process is as follows:

[0013] Step (2.1): Set the target loads for the rigid walls on the left and right sides during the loading process of the discrete element biaxial compression test model. The movement of the rigid walls on the left and right sides is controlled by a servo method. The upper wall moves downward and the lower wall moves upward. The movement speed is fixed and meets the quasi-static loading conditions.

[0014] In step (2.2), set the particle-particle contact parameters and cancel the damping set in step (1.3); the calculation termination condition is that the axial strain reaches the set value, and the axial strain calculation method is as follows:

[0015]

[0016] In the formula, h is the current height of the discrete element biaxial compression test model; h0 is the initial height of the discrete element biaxial compression test model;

[0017] Step (3) involves performing discrete element method (DEM) calculations, calling the servo method at each time step during the calculation process; causing the rigid wall to compress the particles at a set displacement velocity to apply the load, the process being as follows:

[0018] Step (3.1): Assuming the current time step is the nth time step and represented by a superscript in a variable, call the built-in function of the discrete element method to obtain the following parameters for the nth time step: the number of contacts N between the rigid wall and the particles. nContact stress of rigid walls Time step Δt n Define N n-1 The number of contacts between the rigid wall and the particles at time step n-1. Let Δt be the contact stress of the rigid wall at time step n-1. n-1 For a time step of n-1; first, for variable N... n-1 , Δt n-1 An initialization process is performed as follows:

[0019]

[0020] Step (3.2): Calculate the initial value of the rigid wall displacement velocity according to the servo method.

[0021]

[0022] In the formula, A is the rigid wall area of ​​the discrete element biaxial compression test model in step (1); k is the contact stiffness in step (2); σ req The target load in step (2) is α; α is the relaxation factor, which is empirically taken as α = 0.8; N n N represents the number of contacts between the rigid wall and the particles; when its value is 0, take N as the contact quantity. n =1;

[0023] Step (3.3) is based on the difference in contact stress of the rigid wall between the current time step and the previous time step. Calculate the average displacement velocity of particles in contact with the rigid wall. And take this average displacement velocity as Correction value, to obtain the wall displacement velocity

[0024]

[0025] In step (3.4), to prevent excessive displacement velocity of the rigid wall from causing particle penetration, a velocity threshold is set based on the minimum particle radius. When the wall displacement velocity in step (3.3) is... When the speed threshold is exceeded, the following formula is used for correction:

[0026]

[0027] In the formula, u max R is the speed threshold. min The minimum particle radius.

[0028] Step (3.5), update the n-1 time step variable in step (3.1):

[0029]

[0030] Step (3.6): Update the rigid wall position vector based on the displacement velocities of the left and right rigid walls.

[0031]

[0032] In the formula, This is the position vector at time step (n+1). Let n be the displacement vector at time step n. This is the unit normal vector of the plane containing the rigid wall;

[0033] Step (4): When the calculation termination condition set in step (2) is met, the discrete element calculation terminates.

[0034] In step (1.4), calculations are iterated based on Newton's second law to make the overlapping particles bounce off under the action of contact force.

[0035] To ensure calculation accuracy, the contact parameters between the rigid wall and particles in step (2) are set to the parameters corresponding to the linear elastic model, and no restrictions are placed on the contact between particles.

[0036] In step (2), the contact parameters between the rigid wall and the particles include contact stiffness.

[0037] In step (2), the contact parameters between the rigid wall and the particles also include the coefficient of friction and the bond strength.

[0038] In step (4), the discrete element calculation terminates when the axial strain reaches the set value.

[0039] Working Principle: This invention utilizes a discrete element method for rigid wall servoing. It controls the displacement velocity of the rigid wall to bring the wall contact stress closer to the target stress. The rigid wall displacement velocity consists of an initial value and a correction value. The initial value is calculated based on the difference between the wall contact stress and the target stress at the current time step, while the correction value is calculated based on the difference between the wall contact stress at the current and previous time steps. Compared to traditional servo methods, the proposed method significantly improves servo accuracy with minimal impact on computational efficiency.

[0040] Beneficial effects: Compared with the prior art, the present invention has the following advantages:

[0041] (1) Traditional servo algorithms assume that the particles do not move, resulting in low accuracy. The servo method based on discrete elements in this invention further considers the influence of particle motion when calculating the displacement velocity of the wall, which significantly improves the accuracy of servo compared with traditional servo methods.

[0042] (2) Discrete element models usually have a large number of particles. The servo method proposed in this invention indirectly considers the influence of particle motion by calculating the difference between the contact stress of the rigid wall in the current time step and the previous time step, thereby avoiding the particle traversal problem and ensuring computational efficiency. Attached Figure Description

[0043] Figure 1 This is a flowchart of the rigid wall servo method according to an embodiment of the present invention;

[0044] Figure 2 This is a schematic diagram of the discrete element biaxial compression test used in this invention;

[0045] Figure 3 The stress time history curves of the walls on both sides of the biaxial compression test using the method of this invention and the conventional method are shown.

[0046] Figure 4 This is a time history curve showing the relative error between the method of this invention and the conventional method. Detailed Implementation

[0047] like Figure 1 As shown in the embodiment, the rigid wall servoing method based on discrete element particle motion of the present invention includes the following steps:

[0048] Step (1) establishes a discrete element biaxial compression test model consisting of a rigid wall and particles. The process is as follows:

[0049] Step (1.1) uses mutually perpendicular rigid walls to form a region, forming boundary conditions, and determines the width and height of the region enclosed by the rigid walls; in this embodiment, four mutually perpendicular rigid walls are established to constrain particle displacement.

[0050] Step (1.2) generates overlapping particles within the area enclosed by the rigid wall. In this embodiment, a series of overlapping particles are generated within an area with a width of 750 mm and a height of 1800 mm.

[0051] Step (1.3) sets a small contact stiffness value and a large damping value for the particle-particle contact to facilitate the subsequent generation of the initial model;

[0052] Step (1.4) involves iterative calculations based on Newton's second law to automatically separate overlapping particles under contact force. Simultaneously, the discrete element biaxial compression test model is balanced under damping. The run terminates when the unbalanced force of the model is less than the set default value. The resulting initial discrete element biaxial compression test model is shown below. Figure 2 As shown.

[0053] Step (2) involves setting the target load, contact parameters between the rigid wall and particles, and calculation termination conditions in this embodiment. The process is as follows:

[0054] In step (2.1), the target load of the rigid walls on the left and right sides during the loading process of the discrete element biaxial compression test model is set to 100 kPa. The movement of the rigid walls on the left and right sides is controlled by the servo method. The upper wall moves downward and the lower wall moves upward. The movement speed is fixed and meets the quasi-static loading conditions.

[0055] In step (2.2), set the particle-particle contact parameters and cancel the damping set in step (1.3); the calculation termination condition is that the axial strain reaches 15%, and the axial strain calculation method is as follows:

[0056]

[0057] In the formula, h is the current height of the discrete element biaxial compression test model; h0 is the initial height of the discrete element biaxial compression test model.

[0058] Step (3) involves performing discrete element method (DEM) calculations. At each time step during the calculation process, a servo method is called to cause the rigid wall to squeeze the particles at a set displacement velocity, thereby applying the load. The specific steps are as follows:

[0059] Step (3.1): Assuming the current time step is the nth time step and represented by a superscript in a variable, call the built-in function of the discrete element method to obtain the following parameters for the nth time step: the number of contacts N between the rigid wall and the particles. n Contact stress of rigid walls Time step Δt n Define N n-1 Let n be the number of contacts between the rigid wall and the particles at time step n-1. Contact stress of rigid wall at time step n-1 Δt n-1 For time steps of n-1, and for variable N n-1 , Δt n-1 An initialization process is performed as follows:

[0060]

[0061] Step (3.2): Calculate the initial value of the rigid wall displacement velocity according to the servo method.

[0062]

[0063] In the formula, A is the rigid wall area of ​​the discrete element biaxial compression test model in step (1); k is the contact stiffness in step (2); σ req The target load size in step (2) is α; α is the relaxation factor, which is empirically taken as α = 0.8; when N is obtained in step (3.1) n When N is 0, take N.n =1;

[0064] Step (3.3) is based on the difference in contact stress of the rigid wall between the current time step and the previous time step. Calculate the average displacement velocity of particles in contact with the rigid wall. And take this average displacement velocity as Correction value, to obtain the wall displacement velocity

[0065]

[0066] In step (3.4), to prevent excessive displacement velocity of the rigid wall from causing particle penetration, a velocity threshold is set based on the minimum particle radius. When the wall displacement velocity in step (3.3) is... When the speed threshold is exceeded, the following formula is used for correction:

[0067]

[0068]

[0069] In the formula, u max R is the speed threshold. min The minimum particle radius;

[0070] Step (3.5), update the n-1 time step variable in step (3.1):

[0071]

[0072] Step (3.6): Update the rigid wall position vector based on the displacement velocities of the left and right rigid walls.

[0073]

[0074] In the formula, This is the position vector at time step (n+1). Let n be the displacement vector at time step n. This is the unit normal vector of the plane containing the rigid wall.

[0075] Step (4): When the axial strain reaches 15%, the discrete element calculation terminates.

[0076] The above process was simulated using both the method of this invention and the traditional servo method. During the loading process, the time history curves of the contact stress between the left and right walls using the two servo methods were obtained from... Figure 3 It can be seen that the wall confining pressure of the traditional servo method is always greater than the target confining pressure, while the confining pressure of the servo method disclosed in this invention fluctuates around the target confining pressure; from Figure 4The given relative error time history curve shows that the method disclosed in this invention reduces the relative error by two orders of magnitude.

[0077] The above demonstrates that the method proposed in this invention has higher computational accuracy than traditional servo methods. To compare the computational efficiency of the two servo methods, simulations of this embodiment were performed on the same computer while keeping other conditions consistent. The computational efficiency of the two methods is shown in Table 1. The method disclosed in this invention has the same number of running steps per unit time as the traditional method, indicating that the method of calculating the wall displacement velocity correction value based on the wall contact stress difference between the current time step and the previous time step significantly improves servo accuracy while maintaining computational efficiency.

[0078] Table 1. Computational efficiency of different methods in the embodiments of the present invention.

[0079]

Claims

1. A rigid wall servoing method based on discrete element particle motion, characterized in that: Including the following steps: Step (1) establishes a discrete element biaxial compression test model consisting of a rigid wall and particles; the process is as follows: Step (1.1): Use mutually perpendicular rigid walls to form a region, create boundary conditions, and determine the width and height of the region enclosed by the rigid walls; Step (1.2) involves generating overlapping particles within the rigid wall component area; Step (1.3) sets a small contact stiffness value and a large damping value for the particle-particle contact; Step (1.4) involves iterative calculations to cause the overlapping particles to bounce off under the action of contact force, while the discrete element biaxial compression test model is balanced under the action of damping. The calculation is terminated when the unbalanced force of the discrete element biaxial compression test model is less than the set value, thus obtaining the initial discrete element biaxial compression test model. Step (2) involves setting the target load, contact parameters between the rigid wall and the particles, and the calculation termination conditions; the process is as follows: Step (2.1): Set the target loads for the rigid walls on the left and right sides during the loading process of the discrete element biaxial compression test model. The movement of the rigid walls on the left and right sides is controlled by a servo method. The upper wall moves downward and the lower wall moves upward. The movement speed is fixed and meets the quasi-static loading conditions. In step (2.2), set the particle-particle contact parameters and cancel the damping set in step (1.3); the calculation termination condition is that the axial strain reaches the set value, and the axial strain calculation method is as follows: ； In the formula, This represents the current height of the discrete element biaxial compression test model. The initial height of the discrete element biaxial compression test model; Step (3) involves performing discrete element method (DEM) calculations, calling the servo method at each time step during the calculation process; causing the rigid wall to compress the particles at a set displacement velocity to apply the load, the process being as follows: Step (3.1): Assuming the current time step is the nth time step and represented by a superscript in a variable, call the built-in function of the discrete element method to obtain the following parameters for the nth time step: the number of contacts between the rigid wall and the particles. Contact stress of rigid walls Time step ;definition Let n be the number of contacts between the rigid wall and the particles at time step n-1. The contact stress of the rigid wall at time step n-1 is... The time step is n-1; Step (3.2): Calculate the initial value of the rigid wall displacement velocity according to the servo method. : ； In the formula, The area of ​​the rigid wall in the discrete element biaxial compression test model in step (1); The magnitude of the contact stiffness in step (2), The target load size in step (2); It is a relaxation factor; When it is 0, take ; Step (3.3) is based on the difference in contact stress of the rigid wall between the current time step and the previous time step. Calculate the average displacement velocity of the particles in contact with the rigid wall. And take the average displacement velocity as Correction value, to obtain the wall displacement velocity : ； ； Step (3.4): Set a velocity threshold based on the minimum particle radius. When the wall displacement velocity in step (3.3) When the speed threshold is exceeded, the following formula is used for correction: ； ； In the formula, For speed threshold, The minimum particle radius; Step (3.5), update the n-1 time step variable in step (3.1): ； Step (3.6): Update the rigid wall position vector based on the displacement velocities of the left and right rigid walls. ； In the formula, This is the position vector at time step (n+1). Let n be the position vector at time step n. This is the unit normal vector of the plane containing the rigid wall; Step (4): When the calculation termination condition set in step (2) is met, the discrete element calculation terminates.

2. The rigid wall servoing method based on discrete element particle motion according to claim 1, characterized in that: In step (1.4), calculations are iterated based on Newton's second law to make the overlapping particles bounce off under the action of contact force.

3. The rigid wall servoing method based on discrete element particle motion as described in claim 1, characterized in that: In step (2), the contact parameters between the rigid wall and the particles include contact stiffness.

4. The rigid wall servoing method based on discrete element particle motion according to claim 1, characterized in that: In step (2), the contact parameters between the rigid wall and the particles also include the coefficient of friction and the bonding strength.

5. The rigid wall servoing method based on discrete element particle motion according to claim 1, characterized in that: In step (2), the contact parameters between the rigid wall and the particles are the parameters corresponding to the linear elastic model.

6. The rigid wall servoing method based on discrete element particle motion according to claim 1, characterized in that: In step (2.1), the upper wall moves downward and the lower wall moves upward, with a fixed speed and satisfying the quasi-static loading condition.

7. The rigid wall servoing method based on discrete element particle motion according to claim 1, characterized in that: In step (3.1), the number of contacts between the rigid wall and the particles is... Contact stress of rigid walls Time step Perform an initialization: 。 8. The rigid wall servoing method based on discrete element particle motion according to claim 1, characterized in that: In step (3.1), the discrete element correlation fish function is called to obtain the following parameters for time step n: the number of contacts between the rigid wall and the particles. Contact stress of rigid walls Time step .

9. The rigid wall servoing method based on discrete element particle motion according to claim 1, characterized in that: In step (3.4), a velocity threshold is set based on the minimum particle radius to prevent the rigid wall displacement velocity from being too large, which could cause the particles to penetrate.

10. The rigid wall servoing method based on discrete element particle motion according to claim 1, characterized in that: In step (4), the discrete element calculation terminates when the axial strain reaches the set value.

Images