Method for identifying dominant flow channel in three-dimensional fracture network based on topology network

Abstract

The invention relates to the technical field of fracture network fracture water channel identification, and discloses an identification method for analyzing a dominant flow channel in a three-dimensional fracture network based on a topological network. According to the method, probability distribution parameters such as geometric occurrence, gap width and density of a rock mass multi-scale fracture system in a target area are obtained through field surveying and mapping and three-dimensional scanning, and a spatial topological structure of the three-dimensional fracture system is reconstructed by adopting Monte Carlo method simulation and discrete fracture network modeling technology; secondly, abstracting the fracture network into a three-dimensional topological graph model based on a graph theory principle, defining fracture center points and cross points as graph nodes, and converting fracture sections into weighted edges; the method breaks through the continuous medium hypothesis limitation of traditional seepage analysis, and has important engineering application value in the fields of deep geological energy storage reservoir seepage risk assessment, shale gas fracture network optimization design, rock slope stability analysis and the like; compared with a traditional numerical simulation method, the method has the technical advantages of being high in calculation efficiency, good in prediction precision, high in multi-scale applicability and the like.

Description

Technical Field

[0001] The present invention relates to the technical field of fracture network fracture water channel identification, and in particular to a method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network. Background Art

[0002] In the field of fracture network water channel identification, accurately determining the dominant flow paths within fractures is crucial for fluid resource development, groundwater disaster prevention and control, and pollution control. The formation of dominant flow channels within three-dimensional fracture networks is influenced by the complex coupling of fracture geometry and topological structure. Currently, existing research focuses primarily on analyzing the hydraulic conductivity, fracture size, and distribution patterns of fracture networks.

[0003] However, in a complex three-dimensional fracture network environment, how to accurately identify and quantify these factors that have a key impact on flow remains an unsolved problem.

[0004] On the one hand, the intrinsic connection between the distribution of flow resistance and fracture geometry lacks both sufficient theoretical explanation and practical verification. This makes it difficult to accurately grasp the flow patterns of fluids in fracture networks in practical engineering applications, and thus cannot provide a reliable basis for relevant decision-making.

[0005] On the other hand, traditional simulation methods have obvious limitations in this field. When faced with large-scale, highly complex fracture networks, their computational efficiency is extremely low. Due to the need to process massive amounts of data and complex calculation processes, it not only consumes a lot of time and computing resources, but also easily accumulates errors during the calculation process, affecting the accuracy of the final results. In addition, traditional simulation methods find it difficult to capture the microscopic mechanism of the flow channelization effect and are unable to deeply analyze the fine flow behavior of fluids in fracture networks. This has created a major obstacle to a deep understanding of the movement laws of fracture water and the optimization of engineering design.

[0006] To address the above issues, the present invention proposes a method for identifying dominant flow channels in three-dimensional fracture networks based on a topological network. This method effectively improves the accuracy of identifying seepage paths by combining the statistical distribution characteristics of the geometric parameters of the three-dimensional fracture network to construct an accurate model; converts the complex network into a streamlined topological graph model, significantly improving computational efficiency; uses a comprehensive model to calculate fracture flow resistance, more accurately quantifying the resistance of each fracture segment; and uses a specific algorithm combination to search and identify multiple minimum resistance paths, providing rich data references for hydraulic property analysis. At the same time, the method reduces computational complexity while ensuring accuracy, and has wide applicability and scalability to meet the needs of different engineering scenarios. Summary of the Invention

[0007] In order to solve the technical problems raised in the background technology, the present invention provides a method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network.

[0008] The present invention is implemented by the following technical solution: a method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network, comprising the following steps:

[0009] Step 1: Obtain the statistical distribution characteristics of geometric parameters of the rock mass fracture system in the target area, and construct a three-dimensional discrete fracture network model based on the statistical distribution characteristics;

[0010] Step 2: Remove isolated fractures and disconnected fractures in the three-dimensional fracture network to generate a backbone fracture network model;

[0011] Step 3: converting the three-dimensional fracture network model into a weighted topological graph model, wherein fracture centers and fracture intersections are graph nodes, and fracture segments are weighted edges;

[0012] Step 4: Use the modified cubic law and Barton model to calculate the three-dimensional resistance coefficient of the fracture, quantify the influence of fracture opening and roughness on flow resistance, and define the edge weight of the weighted edge as the flow resistance value;

[0013] Step 5: Use the A star heuristic algorithm to search for the path with the minimum flow resistance from the starting point to the end point;

[0014] Step 6: After obtaining the initial minimum resistance path through the KSP path bifurcation strategy, perform deviation search on the nodes on the path one by one, obtain the first K minimum resistance paths in turn, and establish a candidate set of dominant percolation paths.

[0015] Furthermore, the geometric parameter statistical distribution characteristics of the rock mass fracture system in step 1 include: fracture size distribution characteristics, fracture aperture distribution characteristics, fracture roughness distribution characteristics, and fracture occurrence distribution characteristics;

[0016] In step 1, the distribution parameters of the constructed three-dimensional discrete fracture network model are as follows: the fracture trace length and the aperture obey the log-normal distribution, and the probability density function of the log-normal distribution is expressed as:

[0017]

[0018] μ1 is the mean of the lognormal distribution, σ1 is the standard deviation of the normal distribution, and σ2 is the standard deviation of the lognormal distribution.

[0019] Furthermore, in step 2, removing isolated cracks and non-connected cracks includes:

[0020] An isolated fissure is a fissure that does not intersect with any other fissures;

[0021] A non-connected fracture is a fracture that does not belong to any connected path connecting the inlet and outlet boundaries.

[0022] Furthermore, the steps of using the adjacency matrix to remove isolated cracks and disconnected cracks are as follows:

[0023] Establish an adjacency matrix that reflects the positional relationship of the fractures in the fracture network; if fracture i and fracture j directly intersect, then a(i, j) is stored as 1; if the two fractures do not intersect, then it is stored as 0. By traversing all fractures, the adjacency matrix A of the fracture network can be obtained;

[0024] A=[a ij ] n×n

[0025] All elements of the rows and columns corresponding to the isolated crack adjacency matrix are 0, so they can be removed by judging whether each element of the row or column is 0;

[0026] Simplify each crack into a node, starting from the node on the boundary of x = 0 and ending at the node on x = L (L is the length of the boundary domain). Determine a starting node and search for all the adjacent nodes that have not been visited.

[0027] Then, starting from all the nodes adjacent to the starting node, we search for all the nodes adjacent to it that have not been visited, that is, we traverse from top to bottom in layers, and enter the next layer after completing the visit of one layer, until all connected nodes are traversed;

[0028] Finally, all nodes (cracks) between the starting point and the end point constitute the backbone of the three-dimensional discrete crack network.

[0029] Furthermore, in step 3, the three-dimensional fracture network model is converted into a topology model by the following steps:

[0030] The fracture network is abstracted into a weighted graph G = (V, E), where the weights of the weighted edges reflect the resistance of the fracture segment to fluid flow. For each fracture in the three-dimensional discrete fracture network, its geometric center point (i.e., the centroid or geometric center of the fracture in three-dimensional space) is obtained.

[0031] The center points of all cracks are recorded as a type of node in the topological graph, denoted as V c ={v c1 , v c2 ,…,v cn}, where n is the number of cracks in the backbone network;

[0032] Perform intersection detection on any two intersecting fractures in the three-dimensional fracture network. If the two fractures intersect, calculate the coordinates of the midpoint of the intersection line in space and record the midpoint as another type of node in the topological graph, denoted as V l ={v l1 , v l2 ,…,v lm}, where m is the number of midpoints of all intersecting lines;

[0033] Considering the flow direction of water in the fracture network, the inlet plane X = 0 is regarded as the source node S, and the outlet plane X = L is regarded as the target node T. When the fracture intersects the inlet and outlet boundaries, the midpoint of the intersection line is recorded as V b ;

[0034] Merge all the above nodes into a node set of the topology graph: V = V c ∪V l ∪V b ;

[0035] For the center point of the same crack and the end point of its intersection line: if the midpoint of the intersection line is indeed the intersection of the crack and other cracks, then an edge can be established between the center point and this midpoint, recorded as E c ;

[0036] For the center point and boundary node of the intersection line of the crack and the boundary: If the crack does have an intersection line or intersection point with the domain boundary, and its midpoint is regarded as a boundary node, a line can be established between the two to describe the connectivity between the crack and the boundary, which is recorded as E b ;

[0037] Merge all the above edges into a node set of the topological graph: E = E c ∪E b At this point, the three-dimensional fracture network has been transformed into a topological network G = (V, E).

[0038] Furthermore, in step 4, the modified cubic law is combined with the Barton model to calculate the fracture flow resistance, and the steps are as follows:

[0039] Under the assumption that a single fracture layer in a three-dimensional fracture network is simplified to a parallel plate model, the flow rate in the fracture is solved based on the modified cubic law, and the expression is as follows:

[0040]

[0041] Where: q is the flow rate per unit width (m 2 / s), e is the crack opening (m), g is the acceleration of gravity (m / s 2 ), v is the viscosity coefficient of water flow (m 2 / s), l is the crack length (m), J is the hydraulic gradient, ΔH is the head difference (m);

[0042] Considering the influence of crack surface roughness, the Barton model is used to correct the effective crack opening. The expression of the Barton model is:

[0043] e h =JRC 2.5 / (e m / e h ) 2

[0044] Where: JRC is the joint roughness coefficient, ranging from 0-20, e m is the average crack opening (μm), and eh is the hydraulic opening (μm).

[0045] Substituting the corrected hydraulic opening into the modified cubic law, the flow resistance coefficient of each fracture segment in the three-dimensional fracture network is obtained, and the coefficient is used as the weight parameter of the topological graph model for the subsequent minimum resistance path search. The expression is as follows:

[0046]

[0047] Furthermore, the identification algorithm of the least resistance path from the source node to the target node based on the A star heuristic search in step 5 includes the following sub-steps:

[0048] For a three-dimensional fracture network, the basic application form of the A star algorithm is as follows: the total resistance f(n) of each edge is composed of the edge weight g(n) and the heuristic function h(n):

[0049] f(n)=g(n)+h(n)

[0050] Assume that the nth edge is on the crack disk numbered m, and the edge weight g(n) is defined as the three-dimensional flow resistance coefficient Gm of the m crack. The heuristic function is defined as follows:

[0051]

[0052] Where, L n is the length of the side, is the sum of the lengths of the i edges on the fracture disk where the edge is located. Therefore, in the three-dimensional fracture network, the resistance f(n) of the nth edge in the topological network established based on the A star algorithm can be expressed as:

[0053]

[0054] Starting from the source node, with f(n)=g(n)+h(n) as the search criterion, the adjacent nodes are expanded in sequence. When the target node is searched, the minimum flow resistance path from the starting point to the end point in the percolation network can be obtained.

[0055] Furthermore, in step 6, identifying the first K paths with the least resistance based on the KSP path bifurcation strategy includes:

[0056] Based on the initial minimum resistance path, deviation search is performed on each node on the path to generate several feasible alternative paths; the resistance values of the alternative paths are sorted in ascending order;

[0057] The first K paths with the smallest flow resistance are selected from the above sorting results to construct a candidate set of dominant seepage channels.

[0058] The present invention proposes an electronic device, comprising a processor, a memory, and a computer program stored in the memory and configured to be executed by the processor. When the processor executes the computer program, the method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network is implemented.

[0059] The present invention proposes a computer-readable storage medium, which includes a stored computer program, wherein when the computer program is running, the device where the computer-readable storage medium is located is controlled to execute the above-mentioned method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network.

[0060] Compared with the prior art, the present invention has the following beneficial effects:

[0061] This paper constructs a precise discrete fracture network model by combining the statistical distribution characteristics of the geometric parameters of a three-dimensional fracture network. This model fully reflects fracture characteristics such as size, aperture, roughness, and fracture occurrence, providing a scientific basis for subsequent flow resistance calculations and effectively improving the accuracy of seepage path identification.

[0062] The construction of the backbone network fracture network model proposed in this invention transforms the originally complex three-dimensional fracture network into a more streamlined topology model, providing a clear and concise network structure for subsequent path search and greatly improving computational efficiency.

[0063] This paper combines the modified cubic law with the Barton model to calculate fracture flow resistance, fully accounting for the effects of fracture aperture, fracture size, and surface roughness on flow resistance. This comprehensive calculation method more accurately quantifies the resistance of each fracture segment in a fracture network, providing more precise weighting data for subsequent flow path optimization.

[0064] This paper uses the A star heuristic algorithm to search for the path of minimum flow resistance in a three-dimensional fracture network, and combines it with the KSP path bifurcation strategy to further identify the top K minimum resistance paths. This combined approach not only efficiently identifies the mainstream channel but also captures suboptimal paths, generating multiple candidate seepage channels. By comparing and analyzing multiple paths, the concentration and distribution characteristics of seepage in the fracture network can be effectively reflected, providing more reference data for hydraulic property analysis.

[0065] Furthermore, by employing a topological graph model and heuristic algorithms, particularly the combination of the A-star algorithm and the KSP path bifurcation strategy, the present invention significantly reduces computational complexity while maintaining accurate seepage identification. Compared to traditional full-grid numerical simulation methods, the present invention offers greater computational efficiency and can rapidly identify the primary flow channels and multiple alternative paths within a fracture network, significantly improving the efficiency of seepage analysis.

[0066] This method is not only suitable for identifying seepage paths in complex three-dimensional fracture networks, but can also be applied to groundwater flow prediction, rock mass stability analysis, underground engineering design, and other fields. The method is highly versatile and scalable, allowing for adjustment and optimization based on different rock mass structures and engineering requirements. BRIEF DESCRIPTION OF THE DRAWINGS

[0067] Figure 1 A flowchart of a method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network provided in an embodiment of the present application;

[0068] Figure 2 Part (a) is a schematic diagram of the fracture network flow resistance method in two-dimensional space in an embodiment of the present application. Figure 2 Part (b) is a schematic diagram of the crack topology network in three-dimensional space;

[0069] Figure 3 Part (a) is the target three-dimensional discrete crack network diagram in the example of the present invention, Figure 3 (b) is a schematic diagram of the resistance of the crack network topology. Figure 3 Part (c) is a schematic diagram of minimum flow resistance. Figure 3 Part (d) is a schematic diagram of the flow distribution in the three-dimensional fracture network;

[0070] Figure 4 The upper part is a schematic diagram of the frequency of each crack in the crack network in the first 100 paths in the example of the present invention. Figure 4 The lower part is a comparison chart of the total resistance value and the time trend of reaching the boundary of the first ten minimum resistance paths in the fracture network. DETAILED DESCRIPTION

[0071] The present invention will be further described below in conjunction with the accompanying drawings and specific implementation methods. It should be noted that, under the premise of no conflict, the various embodiments or technical features described below can be arbitrarily combined to form new embodiments.

[0072] Example 1:

[0073] This case study takes the left bank slope rock mass of Jinping I Hydropower Station in China as the main research object.

[0074] like Figure 1 As shown, an example of the present invention provides a method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network, comprising the following steps:

[0075] Step 1: Obtain the statistical distribution characteristics of geometric parameters of the rock mass fracture system in the target area, and construct a three-dimensional discrete fracture network model based on the statistical distribution characteristics;

[0076] The statistical distribution characteristics of the geometric parameters of the rock mass fracture system include: fracture size distribution characteristics, fracture aperture distribution characteristics, fracture roughness distribution characteristics, and fracture occurrence distribution characteristics. The distribution parameters of the constructed three-dimensional discrete fracture network model are as follows:

[0077] The crack trace length and opening obey the lognormal distribution. The probability density function of the lognormal distribution is expressed as

[0078]

[0079] μ1 is the mean of the lognormal distribution, σ1 is the standard deviation of the normal distribution, and σ2 is the standard deviation of the lognormal distribution.

[0080] Table 1 Geometric parameters of rock mass fracture system in target area

[0081]

[0082] The target area is constructed with a size of 100m using the Monte Carlo method. 3 The three-dimensional discrete fracture network generates four groups of fractures. The results are as follows Figure 3 Part (a) of .

[0083] Step 2: Remove isolated fractures and disconnected fractures in the 3D fracture network to generate a backbone fracture network model;

[0084] An adjacency matrix is constructed to reflect the positional relationships of the fractures in the fracture network. If fractures i and j directly intersect, a(i, j) is stored as 1; if the two fractures do not intersect, it is stored as 0. By traversing all fractures, the adjacency matrix A of the fracture network can be obtained.

[0085] A=[a ij ] n×n

[0086] All elements of the rows and columns corresponding to the isolated crack adjacency matrix are 0, so they can be removed by judging whether each element of the row or column is 0.

[0087] Each fracture is simplified into a node, starting with the node on the boundary X = 0 and ending with the node on X = L (where L is the boundary domain length). Once a starting node is identified, all adjacent, unvisited nodes are sequentially searched. Starting from these adjacent nodes, all adjacent, unvisited nodes are sequentially searched, layer by layer, from top to bottom. After completing a layer, the next layer is entered, and all connected nodes are traversed until all connected nodes have been traversed. Finally, all nodes (fractures) between the starting and ending points constitute the backbone of the three-dimensional discrete fracture network.

[0088] Step 3: Convert the three-dimensional fracture network model into a weighted topological graph model, where fracture centers and fracture intersections are graph nodes, and fracture segments are weighted edges.

[0089] The fracture network is abstracted into a weighted graph G = (V, E), where the weight of the weighted edge reflects the resistance of the fracture segment to the fluid flow. For each fracture in the three-dimensional discrete fracture network, its geometric center point (i.e., the centroid or geometric center of the fracture in three-dimensional space) is obtained. The center points of all fractures are recorded as a type of node in the topological graph, denoted as V c ={v c1 , v c2 ,…,v cn}, where n is the number of cracks in the backbone network, such as Figure 1 As shown in the schematic diagram.

[0090] Perform intersection detection on any two intersecting fractures in the three-dimensional fracture network. If the two fractures intersect, calculate the coordinates of the midpoint of the intersection line in space and record the midpoint as another type of node in the topological graph, denoted as V l ={v l1 , v l2 ,…,v lm}, where m is the number of midpoints of all intersecting lines, such as Figure 1 The black dashed line and blue nodes in part (a).

[0091] Considering the flow direction of water in the fracture network, the inlet plane X = 0 is regarded as the source node S, and the outlet plane X = L is regarded as the target node T. When the fracture intersects the inlet and outlet boundaries, the midpoint of the intersection line is recorded as V b Merge all the above nodes into a node set of the topology graph: V = V c ∪V l ∪V b

[0092] For the center point of the same crack and the end point of its intersection line: if the midpoint of the intersection line is indeed the intersection of the crack and other cracks, then an edge can be established between the center point and this midpoint, recorded as E c .

[0093] For the center point and boundary node of the intersection line of the crack and the boundary: If the crack does have an intersection line or intersection point with the domain boundary, and its midpoint is regarded as a boundary node, a line can be established between the two to describe the connectivity between the crack and the boundary, which is recorded as E b .

[0094] Merge all the above edges into a node set of the topological graph: E = E c ∪E b At this point, the three-dimensional fracture network has been transformed into a topological network G(V, E). The result of constructing the topological graph in this example is as follows: Figure 3 As shown in part (b) of .

[0095] Step 4: Use the modified cubic law and Barton model to calculate the three-dimensional resistance coefficient of the fracture, quantify the effect of fracture aperture and roughness on flow resistance, and define the edge weight of the weighted edge as the flow resistance value;

[0096] For this example three-dimensional fracture network, the hydraulic head difference is set to 0.5 m. The flow rate in the fracture is solved based on the modified cubic law. The specific expression is as follows:

[0097]

[0098] Where: q is the flow rate per unit width (m 2 / s), e is the crack opening (m), g is the acceleration of gravity (m / s 2 ), v is the viscosity coefficient of water flow (m 2 / s), l is the crack length (m), J is the hydraulic gradient, and △H is the head difference (m). Figure 3 Part (d) is a schematic diagram of the solution results.

[0099] Considering the influence of crack surface roughness, for this example, the crack roughness obeys the log-normal distribution, the mean value of the joint roughness (JRC) is 10, and the Barton model is used to correct the effective crack opening. The expression of the Barton model is:

[0100] e h =JRC 2.5 / (e m / e h ) 2

[0101] Where JRC is the joint roughness coefficient, ranging from 0 to 20, e mis the average crack opening (μm), e h is the hydraulic opening (μm).

[0102] Substituting the above-mentioned corrected hydraulic opening into the modified cubic law, the flow resistance coefficient of each fracture segment in the three-dimensional fracture network is obtained. This coefficient is then used as the weight parameter of the topological graph model for the subsequent search for the minimum resistance path. The expression is as follows:

[0103]

[0104] Step 5: Use the Astar heuristic algorithm to search for the path with the least flow resistance from the starting point to the end point;

[0105] For a three-dimensional fracture network, the basic application form of the Astar algorithm is as follows: the total resistance f(n) of each edge is composed of the edge weight g(n) and the heuristic function h(n):

[0106] f(n)=g(n)+h(n)

[0107] The resistance f(n) of the nth edge in the topological network established based on the Astar algorithm can be expressed as:

[0108]

[0109] Starting from the source node, with f(n) = g(n) + h(n) as the search criterion, the adjacent nodes are expanded in sequence. When the target node is found, the minimum flow resistance path from the starting point to the end point in the percolation network can be obtained. The following steps are used to Figure 2 Take part (a) as an example:

[0110] Add the starting point S to the open list, which stores the nodes to be expanded, as shown in the figure, which are nodes 1 and 2. Since the starting point S is fictitious, the first edge does not participate in the calculation. The selection of the entry boundary node depends on the size of nodes 1 to 3 and nodes 2 to 4: f(3) = 5, f(4) = 7;

[0111] Expand node 3 (with the least total resistance) and move node 3 to the visited list;

[0112] Expand node 7 (the open list has only one node), move node 7 to the visited list, and the resistance from node 7 to node 8 is f(8) = 5; the resistance from node 7 to node 10 is f(10) = 5;

[0113] Expand node 10 (total resistance f is minimum) and repeat the expansion until the end point T is moved to the visited list.

[0114] After stopping the expansion, the crack center point node is extracted from the closed list to obtain the minimum resistance path: S-3-7-12-18-20-T.

[0115] By applying the above steps to actual cases, we can get the following results: Figure 3 Part (c) shows the path of least flow resistance.

[0116] Step 6: After obtaining the initial minimum resistance path through the KSP path bifurcation strategy, perform deviation search on the nodes on the path one by one, obtain the first K minimum resistance paths in turn, and establish a candidate set of dominant percolation paths. Figure 1 Taking the fracture network shown in the figure as an example, the specific process is as follows:

[0117] The A* (A star) algorithm is used to identify the first minimum resistance path L1: S-3-7-12-18-20-T, and its total resistance value R1 is recorded.

[0118] The node blocking mechanism is used to offset the key nodes of the path (nodes with a node count greater than 1 in the open list). For example, the downstream branch node 12 of node 7 is blocked to perform local topology correction, forcing the algorithm to deviate from the current path to generate a suboptimal candidate path L2. The blocked nodes are then restored to the "visited" list to avoid repeated searches.

[0119] The above deviation-search-recovery operation is iterated until there are no expandable nodes in the open list or the top K paths are successfully obtained.

[0120] By globally sorting the total resistance values of each path {R1, R2, …, RK}, K paths with the least resistance are screened out, and their topological connection relationships and the corresponding crack geometric parameters (opening, trace length, JRC) are stored.

[0121] By applying this step to the actual case 1, the result is as follows Figure 4 The path of least flow resistance shown is the dominant seepage channel in the fracture.

[0122] Overall, the A-star algorithm leverages heuristic functions to effectively guide search directions, demonstrating improved search efficiency and computational performance in complex three-dimensional discrete fracture networks. These results demonstrate the potential and advantages of the A-star algorithm for identifying dominant seepage paths in high-dimensional and highly heterogeneous fracture networks, providing a viable technical path for subsequent large-scale three-dimensional seepage simulations and engineering applications.

[0123] The above embodiments are only preferred embodiments of the present invention and cannot be used to limit the scope of protection of the present invention. Any non-substantial changes and replacements made by technicians in this field on the basis of the present invention fall within the scope of protection required by the present invention.

Claims

1. A method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network, characterized in that: The steps include: Step 1: Obtain the statistical distribution characteristics of geometric parameters of the rock mass fracture system in the target area, and construct a three-dimensional discrete fracture network model based on the statistical distribution characteristics; Step 2: Remove isolated fractures and disconnected fractures in the three-dimensional fracture network to generate a backbone fracture network model; Step 3: converting the three-dimensional fracture network model into a weighted topological graph model, wherein fracture centers and fracture intersections are graph nodes, and fracture segments are weighted edges; Step 4: Use the modified cubic law and Barton model to calculate the three-dimensional resistance coefficient of the fracture, quantify the influence of fracture opening and roughness on flow resistance, and define the edge weight of the weighted edge as the flow resistance value; Step 5: Use the A star heuristic algorithm to search for the path with the minimum flow resistance from the starting point to the end point; Step 6: After obtaining the initial minimum resistance path through the KSP path bifurcation strategy, perform deviation search on the nodes on the path one by one, obtain the first K minimum resistance paths in turn, and establish a candidate set of dominant percolation paths.

2. The method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network according to claim 1, characterized in that: The geometric parameter statistical distribution characteristics of the rock mass fracture system in step 1 include: fracture size distribution characteristics, fracture aperture distribution characteristics, fracture roughness distribution characteristics and fracture occurrence distribution characteristics; In step 1, the distribution parameters of the constructed three-dimensional discrete fracture network model are as follows: the fracture trace length and the aperture obey the log-normal distribution, and the probability density function of the log-normal distribution is expressed as: μ1 is the mean of the lognormal distribution, σ1 is the standard deviation of the normal distribution, and σ2 is the standard deviation of the lognormal distribution.

3. The method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network according to claim 1, wherein: In step 2: The isolated fissure is a fissure that does not intersect with any other fissures; The non-connected crack is a crack that does not belong to any connected path connecting the two boundaries of the inlet and the outlet.

4. The method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network according to claim 1, wherein: Use the adjacency matrix to remove isolated cracks and disconnected cracks. The steps are as follows: Establish an adjacency matrix that reflects the positional relationship of the fractures in the fracture network; if fracture i and fracture j directly intersect, then a(i, j) is stored as 1; if the two fractures do not intersect, then it is stored as 0. By traversing all fractures, the adjacency matrix A of the fracture network can be obtained; A=[a ij ] n×n All elements of the rows and columns corresponding to the isolated crack adjacency matrix are 0, and then they are removed by judging whether each element of the row or column is 0; Simplify each crack into a node, starting from the node on the x=0 boundary and ending at the node on x=L. Determine a starting node and search for all the adjacent nodes that have not been visited. Then, starting from all the nodes adjacent to the starting node, we search for all the nodes adjacent to it that have not been visited, that is, we traverse from top to bottom in layers, and enter the next layer after completing the visit of one layer, until all connected nodes are traversed; Finally, all the nodes between the starting point and the end point constitute the backbone of the three-dimensional discrete fracture network.

5. The method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network according to claim 1, wherein: In step 3, the three-dimensional fracture network model is converted into a topology model by the following steps: The fracture network is abstracted into a weighted graph G = (V, E), where the weight of the weighted edge reflects the resistance of the fracture segment to fluid flow. For each fracture in the three-dimensional discrete fracture network, its geometric center point is obtained. The center point is the centroid or geometric center of the fracture in three-dimensional space. The center points of all cracks are recorded as a type of node in the topological graph, denoted as V c ={v c1 , v c2 ,…,v cn }, where n is the number of cracks in the backbone network; Perform intersection detection on any two intersecting fractures in the three-dimensional fracture network. If the two fractures intersect, calculate the coordinates of the midpoint of the intersection line in space and record the midpoint as another type of node in the topological graph, denoted as V l ={v l1 , v l2 ,…,v lm }, where m is the number of midpoints of all intersecting lines; Considering the flow direction of water in the fracture network, the inlet plane X = 0 is regarded as the source node S, and the outlet plane X = L is regarded as the target node T. When the fracture intersects the inlet and outlet boundaries, the midpoint of the intersection line is recorded as V b ; Merge all the above nodes into a node set of the topology graph: V = V c ∪V l ∪V b ; For the center point of the same crack and the end point of its intersection line: if the midpoint of the intersection line is indeed the intersection of the crack and other cracks, then an edge can be established between the center point and this midpoint, recorded as E c ; For the center point and boundary node of the intersection line of the crack and the boundary: If the crack does have an intersection line or intersection point with the domain boundary, and its midpoint is regarded as a boundary node, a line can be established between the two to describe the connectivity between the crack and the boundary, which is recorded as E b ; Merge all the above edges into a node set of the topological graph: E = E c ∪E b At this point, the three-dimensional fracture network is transformed into a topological network G = (V, E).

6. The method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network according to claim 1, wherein: In step 4, the modified cubic law is combined with the Barton model to calculate the fracture flow resistance, and the steps are as follows: Under the assumption that a single fracture layer in a three-dimensional fracture network is simplified to a parallel plate model, the flow rate in the fracture is solved based on the modified cubic law, and the expression is as follows: Where: q is the flow rate per unit width (m 2 / s), e is the crack opening (m), g is the acceleration of gravity (m / s 2 ), v is the viscosity coefficient of water flow (m 2 / s), l is the crack length (m), J is the hydraulic gradient, ΔH is the head difference (m); Considering the influence of crack surface roughness, the Barton model is used to correct the effective crack opening. The expression of the Barton model is: and h =JRC 2.5 / (and m / and h ) 2 Where: JRC is the joint roughness coefficient, ranging from 0-20, e m is the average crack opening (μm), e h is the hydraulic opening (μm). Substituting the corrected hydraulic opening into the modified cubic law, the flow resistance coefficient of each fracture segment in the three-dimensional fracture network is obtained, and the coefficient is used as the weight parameter of the topological graph model for the subsequent minimum resistance path search. The expression is as follows:

7. The method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network according to claim 1, wherein: The identification algorithm of the least resistance path from the source node to the target node based on the A star heuristic search in step 5 includes the following sub-steps: For a three-dimensional fracture network, the basic application form of the Astar algorithm is as follows: the total resistance f(n) of each edge is composed of the edge weight g(n) and the heuristic function h(n): f(n)=g(n)+h(n) Assuming that the nth edge is on the crack disk numbered m, the edge weight g(n) is defined as the three-dimensional flow resistance coefficient G of the m crack m , the heuristic function is defined as follows Where, L n is the length of the side, is the sum of the lengths of the i edges on the fracture disk where the edge is located. Therefore, in the three-dimensional fracture network, the resistance f(n) of the nth edge in the topological network established based on the A star algorithm can be expressed as: Starting from the source node, with f(n)=g(n)+h(n) as the search criterion, the adjacent nodes are expanded in sequence. When the target node is searched, the minimum flow resistance path from the starting point to the end point in the percolation network can be obtained.

8. The method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network according to claim 1, wherein: In step 6, identifying the first K least-resistance paths based on the KSP path bifurcation strategy includes: Based on the initial minimum resistance path, deviation search is performed on each node on the path to generate several feasible alternative paths; the resistance values of the alternative paths are sorted in ascending order; The first K paths with the smallest flow resistance are selected from the above sorting results to construct a candidate set of dominant seepage channels.

9. An electronic device, characterized in that: The method comprises a processor, a memory, and a computer program stored in the memory and configured to be executed by the processor, wherein when the processor executes the computer program, the method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network as described in any one of claims 1 to 8 is implemented.

10. A computer-readable storage medium, characterized in that The computer-readable storage medium includes a stored computer program, wherein, when the computer program is executed, the device where the computer-readable storage medium is located is controlled to execute the method for identifying dominant flow channels in a three-dimensional fracture network based on a topological network as described in any one of claims 1 to 8.

Images