Contact detection algorithm for arbitrarily shaped blocks based on contact theory in 2D-DDA

Contact theory proposed by Shi expresses the contact problem mathematically by establishing an entrance block and reference point, which greatly reduces the complexity of the contact problem. In this paper, a contact detection algorithm based on contact theory for blocks of arbitrary shapes is proposed and integrated into two-dimensional discontinuous deformation analysis (2D-DDA). First, a new vertex-to-edge contact cover detection algorithm that considers contact transfer information is proposed. Second, a filtering algorithm is designed to preserve contact covers based on contact transfer information. Third, an efficient algorithm is proposed to combine contact patterns based on contact covers. The proposed algorithm is verified by using classic examples, and the results are compared with the results of classic 2D-DDA. The findings confirm that the proposed algorithm is feasible and efficient.


Introduction
Discontinuous deformation analysis (DDA) [1], an implicit discrete element method that can simulate the static and dynamic mechanical behaviors of discrete block systems, has been successfully applied to rock engineering, particularly landslides [2][3][4][5], rock fracturing failure [6][7][8], and underground engineering [9][10][11]. Similar to the distinct element method [12] and other discontinuum-based numerical methods, DDA must solve the complex contact problem between discrete blocks, which greatly differs from continuum-based methods. Furthermore, DDA builds a global equilibrium equation based on the principal of minimum potential energy and guarantees the correct handling of contact between blocks by open/close iteration. Therefore, DDA achieves better results compared with other solutions in the large-displacement and deformation problems of block systems. While twodimensional DDA (2D-DDA) has been successfully applied and valuable work [13][14][15][16] has been done to solve contact problems in three-dimensional DDA (3D-DDA), poor performance is still obtained because the contact between 3D polyhedral blocks is much more complex than that in 2D cases. Contact theory [17] has been proposed to overcome the contact problem; this theory provides a unified framework that could be used to describe the contact between blocks in any dimension and remarkably simplifies the complexity of the contact problem [18]. However, research on contact theory remains in the early stages.
Contact detection in DDA consists of three steps: rough contact detection, precise contact detection, and contact transfer. Rough contact detection searches for all block pairs that may come into contact using one of four methods, namely, binary search [19], double-ended spatial sorting (DESS) [20], multi-shell cover [21], and multi-cover searching [22]. The former two methods use one bounding box or sphere to approximate one block, while the latter two methods use more shape covers for subpartitions of each block to obtain more information for precise contact detection. All block pairs found in this step are then forwarded to the precise contact detection step.
The precise contact detection step determines the exact contact positions of block pairs to be used when building the global stiffness equation. Two contact types, namely, vertex-to-vertex and vertexto-edge, have been studied. The parallel edge-to-edge contact type can be converted into two vertexto-edge contacts. Classic 2D-DDA uses distance and angle criteria to determine contact types. Some improvement studies [23][24][25][26] based on the classic contact detection algorithm have been reported.
The contact transfer step transfers the contact state of the same block pairs in the last time step to current step. For example, if a vertex-to-edge contact appears in adjacent time steps, the contact state in the last step should be used as the initial contact state of current step. Open, friction, lock, and tension describe the four contact states; of these, the last three contact states are collectively called the closed state.
In this paper, a 2D contact detection algorithm for blocks of arbitrary shapes is implemented. This paper is organized as follows. Section 2 gives a brief introduction of contact theory and the implementation of the algorithm. Section 3 verifies the proposed algorithm by using some classic examples. Finally, Section 4 gives the conclusions.

Introduction to contact theory
In contact theory [13], contact is represented and proven by mathematics and geometry. The complex contact relationships of two blocks can be represented by the contact between a reference point (i.e., any point in space) and an entrance block E(A, B). The entrance block is defined by contact theory as follows. For any two blocks A and B, select any point a 0 in space (usually the center of block A), and then fix block B. Next, translate the reference point a0 and block A synchronously and record the positions of the reference point when the two blocks intersect; the entrance block will be composed of the recorded positions of the reference point. As shown in Figure 1, the center of block A is used as the reference point, and the generated entrance block of blocks A and B is a hexagon. If the reference point is in or on the boundary of the entrance block, two triangular blocks come into contact. In general, only the boundary of the entrance block, i.e., the dashed boundary line segments of the hexagon in Figure 1, is used in precise contact detection. Each dashed line segment is a contact cover that represents a vertex-to-edge contact.

Calculation of contact cover
As described above, 2D contact covers are line segments. According to contact theory, the contact relationships of two blocks can be determined by the relationship between the reference point and the 3 entrance block. Thus, if the vertex is used as a reference point for each vertex-edge pair, all contact covers can be obtained for every pair of blocks. In this paper, contact covers are retained if they pass following steps. First, no part of the half space of the edge and the angle space of the vertex must be shared, as shown in Figure 2, i.e., the dot products of the two edge vectors of the angle and the inner normal of the edge must satisfy Eq. (1), where ℎ 1 and ℎ 2 are the edge vectors of the angle and is the inner normal of the edge. Second, the projection point of the vertex is on the edge and the penetration should not be not larger than the distance threshold. Third, the contact is retained if the penetration is smaller than the negative threshold but used to build the stiffness equation in the last step. Fourth, the contact is retained if the projection point is close to the end points of the edge. Fifth, the contact is retained if it is in a close state in last step or its penetration is larger than the last step, as shown in Figure 3.

Identification of contact patterns
Two types of contact patterns may be defined for 2D contact: vertex-to-vertex contact and vertex-toedge contact. Classic 2D-DDA iterates all possible contacts, including (1)

all pairs of vertices of blocks A and B, (2) all pairs of vertices of block A and edges of block B, and (3) all pairs of edges of block A and vertices of block B.
Theoretically, every vertex-to-edge contact may be checked three times. Taking this possibility into account, this paper proposes a contact algorithm based on contact theory that searches all contact covers and then combines contact patterns by contact cover. The process includes five steps: 1. Find all vertex-to-edge contact covers whose projection point of the vertex is close to the endpoints of the edge. 2. Tag every vertex-to-edge contact cover in the last step by the index of the endpoint that the projection point is close to for each block. 3. Search for the vertices of each block in which more than one vertex-to-edge contact cover exists and a vertex-to-vertex contact pattern is present. 4. Find all pairs of vertex-to-edge contact covers that belong to one vertex-to-vertex contact pattern and whose vertices are in two blocks. 5. Mark all remaining vertex-to-edge contact covers as a vertex-to-edge contact. Figure 4 shows three blocks, namely, A, B, and C. Two vertex-to-edge contact covers exist for blocks A and B, i.e., C(v 4, e1) and C(v4,e3), where C refers to the contact cover, vi refers to the i-th vertex, and ej refers to the j-th edge. These two contact covers are related to vertex v4; thus, a vertexto-vertex contact pattern exists for vertices v1 and v4. The vertex-to-edge contact covers of blocks B and C are C(v7, e6) and C(v6,e9), and each of these covers is related to the other's vertex; thus, a vertexto-vertex contact pattern is generated for vertices v6 and v7.

Verification of the analytical solution by the sliding block model
We first test the proposed algorithm with the sliding model and then compare the results with those of analytical solutions. The model has two blocks: a slope block with inclination of 45° and a small triangular block on the slope block. The physical parameters of the blocks include a density of 2800 kg/m 3 , elastic modulus of 2 GPa, and Poisson's ratio of 0.25. The numerical control parameters include a step time of 0.02 s, maximum step displacement ratio of 0.1, and contact spring stiffness of 4 GPa. This example is simulated for 40 steps with different friction angles of 10°, 20°, and 30°. The displacement results of the proposed algorithm are in good agreement with those of the analytical solutions, as shown in Figure 5.

Comparison with current precise contact algorithms
The engineering application case (i.e., example 23) is simulated with the classic contact algorithm proposed by Shi, the classic algorithm with vertex-to-edge contact cover [26], and the proposed algorithm. Ref. [26] replaces the rough contact algorithm with DESS [20] and then uses the same precise contact algorithm workflow employed in the classic algorithm but includes the vertex-to-edge contact cover based on contact theory. This precise contact algorithm is implemented in this paper to compare results. In classic example 23, the physical parameters of the blocks include a density of 2548 kg/m 3 , elastic modulus of 327 MPa, and Poisson's ratio of 0.24. The physical parameters of the joints include a friction angle of 20°, cohesion of 0.05 MPa, and tensile strength of 0.03 MPa. The numerical control parameters include an automatic step time, maximum step displacement ratio of 0.001, and contact spring stiffness of 1 GPa.
The proposed algorithm and the classic algorithm with vertex-to-edge contact cover [26] are integrated into the classic DDA program presented by Shi to ensure consistent conditions. The time consumptions of the proposed algorithm, the classic algorithm with vertex-to-edge contact cover, and the classic algorithm are 261, 341, and 291 s, respectively. In addition, the average vertex-to-edge contact cover numbers of these three algorithms over all steps are 3756, 4996, and 4792, respectively. The proposed algorithm requires the least computation time and generates the smallest number of vertex-to-edge contact covers. Compared with the two other precise contact detection algorithms that may check every vertex-to-edge contact twice for vertex-to-vertex or vertex-to-edge contact patterns, the proposed algorithm guarantees that every vertex-to-edge contact is checked only once, which reduces the computation cost. Moreover, contact cover filtering allows the proposed algorithm to generate as few vertex-to-edge contact covers as possible, which further reduces the time cost for the subsequent formation of linear equations. As shown in Figure 6, the simulation results of the three algorithms are generally consistent, but slight differences may be observed after 1000 steps. Differences among the simulation results may be attributed to occasional deviations of the contact detection results due to logic differences among the contact algorithms. While the deviations of the simulation results of different algorithms are inevitable, especially for complex cases, this topic deserves further investigation in future research.

Conclusion
In this paper, we propose a precise contact detection algorithm based on contact theory integrated in 2D-DDA. Compared with the classic precise contact detection algorithm in 2D-DDA, the proposed algorithm has two advantages. First, it reduces the complexity of searching contact covers by ensuring that every vertex-to-edge contact is checked only once. Second, the algorithm adopts contact patterns generated by searching contact covers. The proposed algorithm is verified to be feasible and efficient by using classic 2D-DDA examples.