Discontinuous mechanical analysis of manifold element strain of rock slope based on open source Gempy

Three-dimensional geological modeling is an important means to transparently study the three-dimensional geological spatial structure and change trend of major complex rock mass projects, and plays a very important role in the decision-making and mining of major projects. Today's modeling methods rely too much on the required data, the process is cumbersome and the time is relatively long. Therefore, based on the above discussion, this article uses a complete open source code geological modeling method here, based on the implicit potential field interpolation method of the open source GemPy package, using this interpolation algorithm can be relatively simple to build a complex full three-dimensional geological model, Including fault network, fault surface interaction, unconformity and dome structure. The article describes in detail how to construct the faults and strata of the 3D geological model, and how to add topographic maps, and use the generated model as input data, use MOOSE to voxelize it, and export the required data, use NMM Perform stability analysis. This method has the advantages of simple operation, fast modeling speed, and visual interactive operation. The establishment of a three-dimensional geological model of the fractured rock mass was very effective, laying a solid foundation for the subsequent stability


Introduction
In recent years, major national rock mass projects represented by the Sichuan-Tibet Railway and the West-East Power Transmission Hydropower Project have been carried out successively. The safety of the construction and service process of the project plays an important role in the country's economic development and national security, so it is complicated The engineering geological modeling of the company plays a key role in the construction and service of rock mass engineering. In the process of rock mining, the establishment of a three-dimensional geological model is of great significance to rock stability analysis and rock mining planning. In different time periods of the rock mass, the information of the stratum will change greatly, and faults and joints will also appear. Therefore, we usually capture the knowledge of underground geological features in the form of geological models, that is, the background of geometric structures. Three-dimensional representation. With the continuous development of science and technology, there are many ways to generate 3D geological models. However, the main problem is that the current commercial software charges too high. On the basis of the geological model, considering the impact of joints and faults in the rock mass on the safety of the rock mass engineering, efficient continuous-discontinuous mechanical analysis of the rock mass engineering is to serve the safety design and construction of the rock mass engineering The key link. In recent years, 3D geological modeling has continuously received attention and research from experts and scholars. Many scholars at home and abroad have carried out a lot of research and practice in this area. In terms of large-scale commercial software, such as GoCAD, Vulan, Gemcom, etc., these softwares have matured in modeling technology. At present, the widely used geological modeling software also includes international software brands ArcGis, Revite, Surface, and domestic software brands MapGis. The rock engineering mechanics analysis software includes UDEC, 3DEC, PFC, and HOSS based on finite discrete element method (FDEM).
In China, experts and scholars have made long-term progress in the technical research of 3D geological modeling.
Based on the source of the data and the analysis of the data, the modeling method is summarized. Among them, He Huaijian et al. [1] determined the stratigraphic interface according to the order of the strata boundary points in the geological drilling data; He Yuhao, Chen Yong et al. [2] focused on the 3D modeling data source and data interpolation, panel modeling, volume Six aspects of meta-modeling, hybrid integrated modeling, automatic modeling and dynamic model update, summarize the research progress of orebody 3D modeling methods.
Based on the existing theories, a more streamlined modeling mechanism is proposed to better generate 3D geological models. Among them, Pang Qinggang et al. [3] proposed a three-dimensional stratum modeling technology that can deal with non-uniform data based on multi-scale CSRBFs interpolation method, which can carry out high-precision and rapid reconstruction of stratum surface models for non-uniform data; Zhu Liangfeng et al. [4] Using borehole data as the source data for stratum modeling, a simple, practical, fast and robust 3D modeling method was realized; Yang Wenhuan et al. [5] A 3D stratum model construction method based on a generalized triangular prism model of a tetrahedral grid, Carried out very effective modeling and calculation for mine geological analysis; Guo Jiateng et al. [6] introduced radial basis implicit function and proposed an automatic ore body modeling method based on implicit surface. This method has a high degree of automation, a smooth model and There is no topology error, which provides a new solution for 3D software development. Zhu Huafa et al. [7] combined the geological radar detection image and borehole data for the three-dimensional modeling of the stratum, and constructed a grid interpolation algorithm based on the surface smoothing condition of the processed geological radar data. And carried out three-dimensional stratum modeling for a certain area. Qu Honggang, Wang Yong et al. [8] proposed a method to generate a three-dimensional geological model based on boundary expression by introducing the topological relationship (adjacent, association, and inclusion) between spatial elements (polygon-arc-node) in the section . And based on this method, a 3D geological modeling and visualization system was developed and successfully applied.
Create complex 3D geological models based on simple theories. Among them, Hou Sanitation et al. [9] put forward the overall construction of fault planes based on the wire frame model based on the analysis of the characteristics of the plane geological map data, and realized the three-dimensional construction of complex fault planes; M. Jing et al. [10] used the network The shape profile diagram realizes the creation and visualization of large-scale complex three-dimensional engineering geological bodies; Stability analysis of rock and soil slopes based on 3D geological modeling. Among them, Linda et al. [11] used NURBS technology to establish a three-dimensional geological model of the study area, and realized the three-dimensional slope stability analysis. Through the above analysis, it can be found that the 3D geological model has been widely used in various fields of earth science.
Based on the GemPy open source geological modeling program, this paper develops a new model data management and export platform, and combines it with the manifold element based on the voxel crack model to realize the whole process analysis from geological modeling to mechanical calculation. The new process has two major advantages. First, the geological modeling program is an open source program, which has better secondary development characteristics and data security characteristics; second, the mechanical calculation uses the manifold element method based on the voxel crack model to avoid The traditional finite element method is difficult to generate computational grids for complex fractured rock masses. The numerical example in Section 3 shows that although the method developed in this paper is relatively preliminary, it has shown good applicability for rock mass engineering.

GemPy geological modeling 2.1 Introduction to GemPy modeling
GemPy is a Python-based, open source geological modeling library. It can build complex 3D geological models based on a potentially powerful implicit method, including various features such as fold structures, fault networks, and unconformities. 3D model created by GemPy is shown in Figure 1. GemPy can perform geological modeling of three-dimensional scenes such as deposited strata, discontinuous faults, and folds. And the generated geological model is the same as the software modeling. The reason why it can generate very realistic and complex geological structures is based on its powerful interpolation algorithm. This interpolation algorithm was proposed by Lajaunie et al. [12] (1997) and by Calcagno et al. [13] (2008) Extended universal cokriging method. Two data types need to be input when modeling, to define the three-dimensional coordinates of the boundary point between the ground plane and the fault plane, and the direction of the pole perpendicular to the surface tilt at any point in the space. GemPy can also define the topology between the layers Relations, make the model better displayed, and use matplotlib to achieve two-dimensional model display and use VTK for interactive display of three-dimensional models. GemPy's stochastic geological modeling method is one of its most advanced features. This is it Different from other commercial software, GemPy uses the above powerful performance to perfectly present the 3D geological model, laying a solid foundation for future stability analysis.

Stratum and fault formation
In GemPy, the method for constructing strata and faults is the same. Using implicit interpolation algorithm, you can construct strata and faults only by giving the required data. The most important core is to define at least two coordinates for each level. Point and a direction, otherwise, the definition will be wrong. When defining strata and faults, you can directly define the names of the strata and faults in the Python program, as well as the points and directions of the required layers, by calling the function to put the points belonging to the same layer on one layer, and then use the visualization tool Display the defined stratigraphic faults. As shown in Figure 2 and Figure 3, the defined formations and faults are displayed.

Construction of topographic maps
When using the open source GemPy library for 3D geological modeling, the construction of topographic maps is different from the above two layer sequence construction methods. There are two types of topographic map construction, namely the construction of real topographic maps and the construction of random topographic maps. First of all, for the construction of random topographic maps, you can call the geo_model.set_topography function in the open source GemPy library to generate random topographic maps based on the fractal grid, and modify the fractal dimension fd, terrain height difference d_z, range, and corresponding resolution Any coordinated and beautiful random topographic map can be generated. As shown in Figure 4. Secondly, for the real topographic map, it needs to be loaded from the corresponding real topographic raster file. The first step is to download the raster file of the real terrain that you need, and find the required terrain map on some websites that can download the terrain map data, such as geospatial data cloud, Google Earth, OpenTopography, etc. This article chooses to download the required topographic map from OpenTopography, first enter the website, check the area and submit the required data to get the required address data, because the topographic map needs to be displayed in the program with a file ending in .tif, So you can convert the downloaded data, such as https://mygeodata.cloud/converter/kml -to-tiff, After the conversion is completed, because the range of the topographic map is different from the range of the area built in this article, it needs to be modified The scope of the topographic map is adapted to the scope of the construction area. Import the converted data into the terrain processing software QGIS, use the registration tool to process the terrain, modify the lower left corner and the lower right corner of the graph to be the scope of your own modeling. And save it, as shown in Figure 5. In this way, the required real topographic map files are introduced into the program, and the topography is visualized, as shown in Figure 6.

Construction of a complete 3D geological model
The 3D geological model constructed this time has four strata and two faults, and the construction of real terrain and random terrain has been completed. The range of the model built is (0, 5000, 0, 5000, 0, 5000), and the resolution is (80, 80, 80), the dip angles of the faults are 75 degrees and 45 degrees, and the strata from top to bottom are sandstone, shale, schist and gneiss. Among them, for the real terrain, the terrain data is downloaded from OpenTopography. The downloaded area is located near Linzhi, Tibet, my country, and the QGIS is imported to modify the range to 5000, so as to visualize in the program. For random terrain, set fd to 1.9 and height difference d_z to 3000-5000, so as to call the geo_model.set_topography function to visualize the topography. As shown in Figure 7.

Model data export
The data of the established 3D geological model can be exported through the MOOSE program, and can be stored as a file in code form for subsequent use.

Voxelization conversion of model data
On the basis of the above-mentioned CSV file derived from the model, it needs to be further transformed into a voxel model required for manifold element calculation. Intuitively, voxels are three-dimensional pixels. And on the basis of the voxel model, through the boundary recognition algorithm, the voxels at the boundary of different formations are determined and marked as 0, which facilitates the subsequent application of the penalty function spring to ensure the continuity of the displacement field. The conversion from CSV file to voxel model, voxel model slicing and boundary search can be realized by the following algorithm. % Input：voxel_resolution, lattic_size, voxel_blockID, xyz_coor % Output：J voxel_model = zeros(voxel_resolution); il = length(voxel_blockID); for i = 1:il idx = int32(ceil(xyz_coor(i,1:3)./lattic_size)); idx_linear = sub2ind(voxel_resolution, idx(i,1), idx(i,2), idx(i,3)); voxel_model ( voxel_resolution is the resolution of the voxel mode, lattice_size is the size of the voxels, voxel_blockID is the material block index of the voxels, and xyz_coor are the coordinates of the centers of the voxels. voxel_model is the voxel model transformed from the CSV file, as shown in following figure. slice_voxel is one of the slices from the voxel model as shown in Figure. boundary_slce is the slice_voxel incorporated with boundaries as shown in following.

Mechanical analysis
At the transition interface of different rock formations, the displacement field is continuous, but because of the different parameters such as elastic modulus and Poisson's ratio of different lithologies, the strain field is discontinuous. In order to accurately describe this problem of continuous displacement and discontinuous strain, we use the manifold method to calculate. First, establish a regular grid as shown in the figure, which can cover the entire calculation object and establish the mathematical coverage required for manifold element analysis. Based on the above-mentioned digital image, the non-continuous line cutting mathematical coverage of the material interface is used to obtain the physical coverage.
Penalty function springs are set on both sides of the rock formation interface to ensure the continuity of the displacement field on both sides of the rock formation interface. At the same time, it is satisfied that the strain field is discontinuous because of the different physical coverings used in different manifold units. By solving the manifold element linear equations, the global displacement field distribution can be obtained (as shown in the figure), and the strain distribution can be obtained by derivation, as shown in the figure. It can be seen that the calculated displacement field is continuous, while the obtained strain field is discontinuous.
(a) vertical displacement fieldand (b) vertical strain field Figure 12. It can be seen that the displacement field is smooth and continuous,but because of the different lithology in the model, the strain at the rock layer interface is discontinuous.

Conclusion
This paper develops a new process from geological modeling to rock engineering mechanics analysis. Geological modeling is completed by the open source Gempy program. 3D geological modeling is a complex task. The processing of faults and strata is particularly important. Existing domestic and foreign geological modeling methods are not completely open source and have certain limitations. This article Innovative use of GemPy for 3D geological modeling of fractured rock masses, faster and better display of 3D geological models of fractured rock masses, and easier to add faults, joints, folds, etc. to the 3D geological model. It has better secondary development potential and guarantees data security in the process of geological modeling. The mechanical analysis is completed by the manifold element method, and regular grids can be used to uniformly deal with the continuous displacement and discontinuous strain of rock masses of different lithologies at the rock layer interface. This article illustrates the applicability of the geological modeling-mechanical analysis introduced in this article through a rock mass slope with different lithology. At present, mechanical calculations are still limited to two-dimensional problems. The analysis and calculation of three-dimensional rock mass engineering and the simulation of dynamic evolution of faults under load will be the next work.