GPU-accelerated implicit geological modeling methods, systems, and devices

Abstract

Embodiments of the present disclosure disclose a GPU-accelerated implicit geological modeling method and an electronic device. A specific implementation of the method includes: generating a grid set of a calculation region according to geological modeling parameters and geological source data; constructing a linear equation set Ax=b; constructing a multi-grid based on the grid set and determining matrix coefficients of each grid layer; classifying grid points in the multi-grid into different color class sets; for each color class set, starting a same number of GPU threads as the number of grid points in the color class set, wherein each GPU thread is used for processing a Gauss-Seidel iteration calculation of a coefficient matrix row corresponding to a grid point; when all threads in a color class set complete a specified number of iterations or meet a convergence condition, switching to a next color class set until the entire multi-grid converges; and outputting a modeling result according to the converged multi-grid. Thus, a new implicit geological modeling method is provided.

Description

Technical Field

[0001] The present disclosure relates to the field of computer technology, and in particular to an implicit geological modeling method, device, and electronic device based on GPU acceleration. Background Art

[0002] With the continuous advancement of computer technology and the widespread availability of high-performance computing resources, computer-based geological modeling is playing an increasingly important role in resource exploration, environmental protection, disaster prediction, and other fields. Geological modeling, the process of describing, interpreting, and predicting geological entities, is a crucial component of geology, resource exploration, and environmental science. The primary task of geological modeling is to construct subsurface structure and property models based on geological data (such as borehole data, geophysical data, and geological images). This helps understand geological processes, guide resource exploration and management, and assess geological hazard risks.

[0003] Traditional geological modeling methods mainly include manual modeling and rule-based model fitting. These methods have problems such as slow calculation speed, strong subjectivity, and reliance on professional experience when processing large-scale geological data, and are difficult to meet the needs of modern geological engineering. Summary of the Invention

[0004] This disclosure section is provided to briefly introduce concepts that will be described in detail in the detailed description section below. This disclosure section is not intended to identify key features or essential features of the claimed technical solution, nor is it intended to limit the scope of the claimed technical solution.

[0005] In a first aspect, an embodiment of the present disclosure provides an implicit geological modeling method based on GPU acceleration, the method comprising: generating a grid set of a calculation area according to geological modeling parameters and geological source data; constructing a linear equation set Ax=b, wherein A is the coefficient matrix of the grid set, and the matrix vector b is used to characterize point constraints and direction constraints; constructing a multi-grid based on the grid set, and determining the matrix coefficients of each grid layer; classifying the grid points in the multi-grid into color class sets of different colors, wherein there is no mutual dependence between the grid points in each color class set; for each color class set, starting a GPU thread with the same number of grid points as the color class set, wherein each GPU thread is used to process the Gauss-Seidel iterative calculation of a coefficient matrix row corresponding to a grid point; when all threads in a color class set complete the specified number of iterations or meet the convergence conditions, switching to the next color class set until the entire multi-grid converges; and outputting the modeling results according to the multi-grid in the converged state.

[0006] In some embodiments, multiple GPU threads corresponding to a color class set are processed in parallel.

[0007] In some embodiments, classifying all grid points in the multi-grid into color class sets of different colors includes: classifying all grid points in the multi-grid based on vertex coloring and standard greedy coloring, wherein adjacent grid points belong to different color class sets.

[0008] In some embodiments, the method further includes: dividing the error solution of the multiple grids into a high-frequency part and a low-frequency part according to preset judgment rules; for the high-frequency part, using Jacobi iteration or Gauss-Seidel iteration for processing; for the low-frequency part, using iterative V loop for processing.

[0009] In some embodiments, the iterative V-loop processing includes: solving the linear equations on the finest grid using the Gauss-Seidel iterative algorithm; passing the residual of the current grid to the next layer of coarse grid through the interpolation operator, and repeating this process until the residual is restricted to the coarsest grid; accurately solving the linear equations in the coarsest grid to obtain the calculation result; returning the solution of the coarsest grid to the previous layer of fine grid through interpolation; using the new solution to perform Gauss-Seidel iteration to obtain a better solution, and repeating this process until the solution is returned to the finest grid; judging whether the termination condition is met, if the current result meets the preset termination error, the iteration ends; if the termination condition is not met, the Gauss-Seidel iterative calculation is continued with the current result as the initial value.

[0010] In some embodiments, constructing multiple grids based on the grid set and determining the matrix coefficients of each grid layer include: determining the coarse grid point set of the coefficient matrix A according to the coefficient matrix A, and determining the corresponding interpolation matrix P; determining the matrix coefficient Ac of the coarse grid according to the coefficient matrix A and the interpolation matrix P of the fine grid.

[0011] In some embodiments, the method of determining the coarse grid point set of the coefficient matrix A and determining the corresponding interpolation matrix P based on the coefficient matrix A includes: determining the strongly connected points based on the algebraic properties of the diagonal elements and non-zero elements of each row of the fine grid coefficient matrix; dividing the grid points into independent point sets based on the strongly connected points; dividing the rows in all fine grids into respective coarse grid point sets, and determining the intermediate transfer operator based on the clustering results; processing the intermediate transfer operator based on the Jacobi smoothing operator to obtain the interpolation matrix P.

[0012] In some embodiments, determining the matrix coefficients Ac of the coarse grid based on the coefficient matrix A of the fine grid and the interpolation matrix P includes: obtaining the coefficient matrix corresponding to each grid layer by calculating the following formula: c =P T AP

[0013] Where Ac is the coefficient matrix corresponding to the coarse grid, and A is the coefficient matrix corresponding to the fine grid of the previous layer.

[0014] In some embodiments, determining the strong connection points based on the algebraic properties of the diagonal elements and non-zero elements in each row of the fine grid coefficient matrix includes: using the following formula to determine the strong connection points corresponding to each row of the matrix, and the non-zero elements that meet the following formula are marked as the strong connection points of this row:

[0015]

[0016] Where, ε = 0.08 × 0.25 l , represents the current grid layer, a ij Represents the element in the i-th row and j-th column of the coefficient matrix.

[0017] In a second aspect, an embodiment of the present disclosure provides an electronic device comprising: one or more processors; a storage device for storing one or more programs, wherein when the one or more programs are executed by the one or more processors, the one or more processors implement the method described in the first aspect. BRIEF DESCRIPTION OF THE DRAWINGS

[0018] The above and other features, advantages, and aspects of the various embodiments of the present disclosure will become more apparent with reference to the following detailed description in conjunction with the accompanying drawings. Throughout the drawings, the same or similar reference numerals represent the same or similar elements. It should be understood that the drawings are schematic and that the originals and elements are not necessarily drawn to scale.

[0019] Figure 1 is a flow chart of an embodiment of the implicit geological modeling method based on GPU acceleration according to the present disclosure;

[0020] Figure 2 It is a schematic diagram of the coarsening of the algebraic multigrid algorithm;

[0021] Figure 3 It is a schematic diagram of the iterative process of the algebraic multigrid algorithm;

[0022] Figure 4 This is a diagram of vertex coloring using the Gauss-Seidel iterative method;

[0023] Figure 5 This is a schematic diagram of memory access for a GPU-based parallel algebraic multigrid solver;

[0024] Figure 6 is an exemplary overall flow chart;

[0025] Figure 7 This is a diagram of the CUDA multi-threaded architecture;

[0026] Figure 8It is a GPU parallel computing flow chart;

[0027] Figure 9 It is a schematic diagram of the basic structure of an electronic device provided according to an embodiment of the present disclosure. DETAILED DESCRIPTION

[0028] The following describes embodiments of the present disclosure in more detail with reference to the accompanying drawings. Although certain embodiments of the present disclosure are shown in the accompanying drawings, it should be understood that the present disclosure can be implemented in various forms and should not be construed as limited to the embodiments described herein. Rather, these embodiments are provided to provide a more thorough and complete understanding of the present disclosure. It should be understood that the drawings and embodiments of the present disclosure are for illustrative purposes only and are not intended to limit the scope of protection of the present disclosure.

[0029] It should be understood that the various steps described in the method embodiments of the present disclosure may be performed in different orders and / or in parallel. In addition, the method embodiments may include additional steps and / or omit the steps shown. The scope of the present disclosure is not limited in this respect.

[0030] As used herein, the term "including" and its variations are open-ended, i.e., "including but not limited to." The term "based on" means "based, at least in part, on." The term "one embodiment" means "at least one embodiment," the term "another embodiment" means "at least one additional embodiment," and the term "some embodiments" means "at least some embodiments." Other terms are defined in the following description.

[0031] It should be noted that the concepts of "first" and "second" mentioned in this disclosure are only used to distinguish different devices, modules or units, and are not used to limit the order or interdependence of the functions performed by these devices, modules or units.

[0032] It should be noted that the modifications of "one" and "multiple" mentioned in the present disclosure are illustrative rather than restrictive, and those skilled in the art should understand that unless otherwise clearly indicated in the context, they should be understood as "one or more".

[0033] The names of the messages or information exchanged between multiple devices in the embodiments of the present disclosure are only used for illustrative purposes and are not used to limit the scope of these messages or information.

[0034] Implicit geological modeling, a data-driven geological modeling method, automatically infers subsurface structures and properties by leveraging spatial relationships and physical features in geological data to generate a geological model. Compared to traditional methods, implicit geological modeling offers advantages such as high automation, rapid modeling speed, and applicability to complex geological environments. Its core concept is to transform the geological modeling problem into an optimization problem. By minimizing the residual between the geological data and the model, model parameters are optimized to best fit the geological data, thereby achieving geological modeling. However, implicit geological modeling methods are computationally complex, and traditional methods often fail to meet the accuracy, real-time, and efficiency requirements, especially when processing large-scale geological data.

[0035] One or more embodiments of the present application propose an implicit geological modeling method based on GPU acceleration to solve the problems of existing implicit geological modeling methods, such as long calculation time, low model accuracy, difficulty in handling complex constrained structural areas, and low modeling efficiency.

[0036] One or more embodiments of the present application achieve high precision by constructing an algebraic multigrid (AMG) solver based on GPU acceleration to accelerate the solution of linear equations, thereby realizing efficient modeling of large-scale geological data sets and efficient implicit geological modeling.

[0037] One or more embodiments of the present application propose an efficient implicit geological modeling algorithm. Given the significant characteristics of implicit geological modeling data, which are large in scale and highly complex, an algebraic multigrid parallel solver is constructed based on the excellent and powerful parallel processing capabilities of the GPU to address the difficulty of solving the large-scale linear equations involved in the implicit geological modeling algorithm. This reduces the time cost required to solve the linear equations and compresses the model generation time, thereby properly addressing application scenarios that urgently need to quickly obtain geological models to support emergency decision-making.

[0038] One or more embodiments of the present application construct an algebraic multigrid solver based on GPU acceleration, use algebraic multigrid to reduce the size of the linear equation system, and solve it in parallel through GPU, thereby greatly improving the accuracy and efficiency of the implicit geological modeling method.

[0039] Please refer to Figure 1 , which shows the process of an embodiment of the implicit geological modeling method based on GPU acceleration according to the present disclosure. Figure 1 The implicit geological modeling method based on GPU acceleration includes the following steps:

[0040] Step 101: Generate a grid set of the calculation area according to geological modeling parameters and geological source data.

[0041] Step 102: Construct a linear equation system Ax=b.

[0042] In this embodiment, A is the coefficient matrix of the grid set, and the matrix vector b is used to represent point constraints and direction constraints.

[0043] Step 103: construct a multi-grid based on the grid set, and determine the matrix coefficients of each grid layer.

[0044] Step 104 : Classify the grid points in the multi-grid into color class sets of different colors.

[0045] In this embodiment, there is no mutual dependence between the grid points in each color class set.

[0046] Step 105 : For each color class set, enable GPU threads with the same number of grid points as that in the color class set.

[0047] In this embodiment, each GPU thread is used to process the Gauss-Seidel iterative calculation of a coefficient matrix row corresponding to a grid point.

[0048] The color class set is determined based on a standard greedy coloring algorithm. It can be understood that using four colors can achieve that adjacent nodes belong to different color subsets.

[0049] Step 106 : When all threads in a color class set complete the specified number of iterations or meet the convergence condition, switch to the next color class set until the entire multigrid converges.

[0050] Step 107: Output the modeling result according to the converged multigrid.

[0051] In this embodiment, the algebraic multigrid method (AMG) is used to solve the linear equation system Ax = b. By utilizing the hierarchical structure of the problem and approximation techniques, the algebraic multigrid method (AMG) transforms the original problem into a series of smaller problems, thereby significantly improving solution efficiency and convergence speed.

[0052] The multi-grid algorithm utilizes the characteristics of errors and constructs multiple levels of grids to transfer and transform the errors. The low-frequency error components are transferred to the coarse grid and converted into high-frequency error components. This avoids the need for a large number of iterations to solve large-scale linear equations and performs accurate solutions at the coarsest grid layer. At this time, the grid size is relatively small and can converge quickly. Finally, the errors are corrected to the fine grid through interpolation to complete the iterative solution.

[0053] Multigrid algorithms consist of two steps: coarsening and iterative solving. Coarsening is the most important component of multigrid. Geometric multigrid uses the geometric relationship between coarse and fine grids for coarsening, while algebraic multigrid uses algebraic relationships between matrices to complete the coarsening process without grid information. Iterative solving involves performing coarse grid correction and smoothing on the multiple grids generated by coarsening to eliminate errors.

[0054] It should be noted that the majority of the time spent in the multigrid iteration process is spent on the Gauss-Seidel method iterations during the smoothing and solution process. Leveraging the powerful computing capabilities of GPUs, this iterative computation is parallelized using the GPU, building a GPU-accelerated parallel algebraic multigrid solver that significantly reduces the time spent on algebraic multigrid computations.

[0055] It should be noted that the Gauss-Seidel iteration method is sequential in nature and cannot update all interdependent nodes at the same time. In order to enable the Gauss-Seidel iteration method to perform parallel computing on the GPU, this application proposes to improve the Gauss-Seidel iteration. The improved method can be called an iterative V loop to improve its parallelism.

[0056] like Figure 4 As shown, by analyzing the three-dimensional finite difference template, it is found that when performing differential interpolation, for any grid point in the three-dimensional grid, its differential coefficient is only related to the grid points in its 18 neighborhoods, and there is no dependency on other grid points. And extended to the coefficient matrix, each row of the coefficient matrix represents a grid point in the grid, and the non-zero elements of this row only exist in the columns corresponding to the 18 neighborhood grid points. At this time, the coefficients a of the other positions in the coefficient matrix are ij Both are 0, combined with high

[0057]

[0058] The formula of the Steffen-Seidel iteration method is:

[0059] At each iteration, each row will use the latest result of this round of iteration for the j≤i solution, but at this time the coefficient a contained in the j≤i part ij , except for the coefficient a of the 18 neighborhood columns corresponding to the grid point in the row ij May not be 0, the rest a ij Are all 0, no matter which still It has no effect on the actual results, that is, when the j≤i part of the row does not contain the 18 neighborhood grid points corresponding to the corresponding grid point in the row, the sequential dependence of the Gauss-Seidel iteration method does not exist.

[0060] In some embodiments proposed in the present disclosure, the nodes in the color class set can be updated in the order of color classes (i.e., color class sets), that is, the nodes in each color class are updated in sequence, and for each color class, the coefficient matrix rows corresponding to each grid point therein can be processed iteratively in parallel, and the same number of threads as the number of grid points in the color class are opened in the GPU, and each thread is responsible for the Gauss-Seidel iterative calculation of the corresponding row in the coefficient matrix. All threads of the color class together form a package (also called a wrap), and through the wrap, the GPU can achieve merged access to the memory, which has achieved the maximum data exchange speed, such as Figure 3 shown.

[0061] When all threads in a color class reach the specified number of iterations or meet the convergence conditions, the wrap is switched and the calculation of the wrap corresponding to the next color class is started until the entire multigrid system converges.

[0062] In some embodiments, multiple GPU threads corresponding to a color class set are processed in parallel.

[0063] Optionally, each thread is responsible for the Gauss-Seidel iterative calculation of the corresponding row in the coefficient matrix. All threads of the color class together form a package (also called a wrap). Through wrap, the GPU can achieve merged access to memory and achieve the maximum data exchange speed, such as Figure 3 shown.

[0064] In some embodiments, classifying all grid points in the multi-grid into color class sets of different colors includes: classifying all grid points in the multi-grid based on vertex coloring and standard greedy coloring, wherein adjacent grid points belong to different color class sets.

[0065] Here, the vertex coloring method is used to classify all grid points in the grid, such as Figure 5 As shown, grid points that do not have mutual dependencies are classified into subsets of the same color, and all points in a color class can be updated in parallel. Using the standard greedy coloring algorithm, using 4 colors, adjacent nodes can be divided into different color subsets.

[0066] In some embodiments, the method further includes: dividing the error solution of the multiple grids into a high-frequency part and a low-frequency part according to preset judgment rules; for the high-frequency part, using Jacobi iteration or Gauss-Seidel iteration for processing; for the low-frequency part, using iterative V loop for processing.

[0067] Take the first kind of boundary conditions of the one-dimensional Poisson equation as an example:

[0068] Using central difference to simplify the above equation is:

[0069] -u i-1 +2u i -u i+1 =h 2 f

[0070] After simplifying to the matrix form Ax = b, Jacobi iteration is used to solve the linear equations and rewrite the matrix A as follows:

[0071] A=D+L+U

[0072] Where D, L, and U are the diagonal matrix, strictly lower triangular matrix, and strictly upper triangular matrix after the triangular decomposition of matrix A. Construction formula:

[0073] x m+1 =Bx m +f

[0074] Among them, the iteration matrix B = -D -1 (L+U), constant term f=D -1 b, let x be the exact solution, we can get:

[0075] x=Bx+f

[0076] The mth iteration error can be obtained as e (m) =x (m) -x, and the iterative error of the m+1th time is:

[0077] e (m+1) =Bx (m)

[0078] λ k ,v k represents the eigenvalue and eigenvector of the B matrix, k represents the eigenfrequency, and k≤N-1. After derivation, we can get:

[0079] e (m+1) =Bx (m) =a1λ1v1+a2λ2v2+...+a N λ N v N

[0080] The part where 1≤k≤N / 2 is called a low-frequency component, and the part where N / 2≤k≤N is called a high-frequency component.

[0081] It should be noted that in one or more embodiments of the present disclosure, the solution error can be decomposed into two parts. The first part is the high-frequency part relative to the initial grid. By using traditional iterative methods, such as Jacobi iteration or Gauss-Seidel iteration, this part of the error can be quickly reduced. The second part is the low-frequency part, which is very slowly reduced in the initial grid. By converting the low-frequency error into a high-frequency error on the coarse grid, it is easier to eliminate the high-frequency error. That is, by transferring and eliminating errors between grids, the system of equations can be solved efficiently and quickly.

[0082] In contrast, traditional iterative methods such as Jacobi iteration, Gauss-Seidel iteration, and conjugate gradient iteration often show a rapid decrease in error at the beginning of the iteration when solving large linear equations. However, after a few iterations, the convergence speed will become slower and slower, making it difficult for the iteration to converge within an acceptable time. This is because within a grid, high-frequency errors decay much faster than low-frequency errors, and ultimately low-frequency errors remain that are difficult to eliminate.

[0083] In some embodiments, the iterative V-loop processing includes: solving the linear equations on the finest grid using the Gauss-Seidel iterative algorithm; passing the residual of the current grid to the next layer of coarse grid through the interpolation operator, and repeating this process until the residual is restricted to the coarsest grid; accurately solving the linear equations in the coarsest grid to obtain the calculation result; returning the solution of the coarsest grid to the previous layer of fine grid through interpolation; using the new solution to perform Gauss-Seidel iteration to obtain a better solution, and repeating this process until the solution is returned to the finest grid; judging whether the termination condition is met, if the current result meets the preset termination error, the iteration ends; if the termination condition is not met, the Gauss-Seidel iterative calculation is continued with the current result as the initial value.

[0084] Optionally, the iterative V loop can be understood as iteratively solving the linear equations using a k-layer grid. The steps are as follows:

[0085] Solve for the finest grid, the formula is: u 1 =G(A 1 ,f 1 ), where G is the Gauss-Seidel iterative algorithm.

[0086] Limit the residual and pass the residual in the coarsest grid through the interpolation operator Restricted to the coarse grid H2, the formula is: Repeat this process until the residual is limited to the coarsest grid H k middle.

[0087] In the coarsest grid H k In the above example, the linear equations are directly solved to obtain the calculation result γk .

[0088] Return the solution to the fine grid with the formula: k-1 =u k-1 +P k-1 u k

[0089] , and then use the new solution to perform Gauss-Seidel iterations to obtain a better solution. This process is repeated until the solution is returned to the finest grid.

[0090] If the result satisfies the cutoff error at this time, the iteration ends. If not, the V loop iteration continues with the result as the initial value.

[0091] Therefore, by constructing multiple levels of grids, the errors are transferred and converted, and the low-frequency error components are transferred to the coarse grid and converted into high-frequency error components, thereby avoiding a large number of iterations of solving large-scale linear equations and performing accurate solutions at the coarsest grid layer. At this time, the grid scale is relatively small and can converge quickly. Finally, the errors are corrected to the fine grid through interpolation to complete the iterative solution.

[0092] In some embodiments, constructing multiple grids based on the grid set and determining the matrix coefficients of each grid layer include: determining the coarse grid point set of the coefficient matrix A according to the coefficient matrix A, and determining the corresponding interpolation matrix P; determining the matrix coefficient Ac of the coarse grid according to the coefficient matrix A and the interpolation matrix P of the fine grid.

[0093] In some embodiments, the method of determining the coarse grid point set of the coefficient matrix A and determining the corresponding interpolation matrix P based on the coefficient matrix A includes: determining the strongly connected points based on the algebraic properties of the diagonal elements and non-zero elements of each row of the fine grid coefficient matrix; dividing the grid points into independent point sets based on the strongly connected points; dividing the rows in all fine grids into respective coarse grid point sets, and determining the intermediate transfer operator based on the clustering results; processing the intermediate transfer operator based on the Jacobi smoothing operator to obtain the interpolation matrix P.

[0094] Specifically, the coarse grid points of each layer are clustered according to the algebraic properties of the diagonal elements and non-zero elements of each row of the fine grid coefficient matrix. First, the strong connection points corresponding to each row of the matrix are determined using the following formula. The non-zero elements that meet the following formula are marked as strong connection points of this row:

[0095]

[0096] Where, ε = 0.08 × 0.25 l , represents the current grid layer, a ij Represents the element in row i and column j of the coefficient matrix. The grid points are divided into independent point sets according to the following algorithm:

[0097] Initialization: Set And j=0.

[0098] For all i∈R, if Then j=j+1,

[0099]

[0100] Classify the remaining points i in R into the set that the point strongly depends on according to the following rules (This set exists): For all i∈R, if there exists K such that but R = R\i. If there are multiple sets that meet the conditions, the set with the strongest dependency is selected.

[0101] Process the remaining points.

[0102] For the point i that still remains in R, let j=j+1,

[0103] After all rows of the corresponding matrix in the fine grid are divided into their respective coarse grid point sets,

[0104]

[0105] According to clustering Define the intermediate transfer operator:

[0106] Then use a simple Jacobi smoothing operator on P l , and the final interpolation operator P is as follows:

[0107]

[0108] Where D represents A l The diagonal matrix of It is defined as follows:

[0109]

[0110] In some embodiments, determining the matrix coefficients Ac of the coarse grid based on the coefficient matrix A of the fine grid and the interpolation matrix P includes: obtaining the coefficient matrix corresponding to each grid layer by calculating the following formula:

[0111] A c =P T AP

[0112] Where A c is the coefficient matrix corresponding to the coarse grid, and A is the coefficient matrix corresponding to the fine grid of the previous layer.

[0113] In some embodiments, determining the strong connection points based on the algebraic properties of the diagonal elements and non-zero elements in each row of the fine grid coefficient matrix includes: using the following formula to determine the strong connection points corresponding to each row of the matrix, and the non-zero elements that meet the following formula are marked as the strong connection points of this row:

[0114]

[0115] Where, ε = 0.08 × 0.25 l , represents the current grid layer, a ij Represents the element in the i-th row and j-th column of the coefficient matrix.

[0116] It should be noted that by constructing multiple levels of grids, the errors are transferred and converted, and the low-frequency error components are transferred to the coarse grid and converted into high-frequency error components. This avoids a large number of iterations of solving large-scale linear equations, and accurately solves them at the coarsest grid layer. At this time, the grid scale is relatively small and can converge quickly. Finally, the error is corrected to the fine grid through interpolation to complete the iterative solution.

[0117] like Figure 7 As shown, CUDA is a parallel computing platform and programming model launched by NVIDIA that allows developers to leverage the massively parallel computing capabilities of GPUs (graphics processing units) to accelerate general computing tasks. Its core advantage lies in significantly improving the performance of compute-intensive applications through data parallelism and thread-level parallelism. It is particularly widely used in scientific computing, deep learning, image processing, and other fields. Parallel solving of multiple grids will be performed on the GPU through the CUDA framework.

[0118] In CUDA, threads are the most basic program execution unit. CUDA's parallel computing is achieved through the parallel execution of a large number of threads. In CUDA, each executed core program corresponds to a kernel core, and the kernel function in this kernel runs on a grid. Each grid consists of several thread blocks.

[0119] Within each thread block, threads can synchronize with each other and exchange data quickly through shared memory. Each thread block is heterogeneous, allowing threads to be organized in one, two, or three dimensions as needed. Each block and thread has its own ID, allowing you to find the corresponding thread and thread block.

[0120] The general process of CUDA calculation is as follows Figure 8 As shown:

[0121] First, perform initialization, obtain the parameters of the GPU device, and set the current computing device.

[0122] Allocate host and device memory. Host memory is allocated through standard library functions, and device memory is allocated through the corresponding CUDA API.

[0123] Perform data transfer to transfer the data required for parallel computing from the host memory to the device memory.

[0124] Define the kernel function, including the main task to be executed and the configuration information for the grid and thread blocks.

[0125] Execute kernel function tasks in parallel.

[0126] Perform thread synchronization until all threads' subtasks are completed.

[0127] Recycle thread data and transfer execution results from device memory back to host memory.

[0128] Reference below Figure 9 , which shows a schematic diagram of the structure of an electronic device suitable for implementing the embodiments of the present disclosure. The terminal device in the embodiments of the present disclosure may include, but is not limited to, mobile terminals such as mobile phones, laptop computers, digital broadcast receivers, PDAs (personal digital assistants), PADs (tablet computers), PMPs (portable multimedia players), in-vehicle terminals (such as in-vehicle navigation terminals), and fixed terminals such as digital TVs and desktop computers. Figure 9 The electronic device shown is only an example and should not limit the functions and scope of use of the embodiments of the present disclosure.

[0129] like Figure 9 As shown, the electronic device may include a processing device (e.g., a central processing unit, a graphics processing unit, etc.) 701, which can perform various appropriate actions and processes according to a program stored in a read-only memory (ROM) 702 or a program loaded from a storage device 708 into a random access memory (RAM) 703. Various programs and data required for the operation of the electronic device 700 are also stored in the RAM 703. The processing device 701, the ROM 702, and the RAM 703 are connected to each other via a bus 704. An input / output (I / O) interface 705 is also connected to the bus 704.

[0130] Typically, the following devices may be connected to the I / O interface 705: an input device 706 including, for example, a touch screen, a touchpad, a keyboard, a mouse, a camera, a microphone, an accelerometer, a gyroscope, etc.; an output device 707 including, for example, a liquid crystal display (LCD), a speaker, a vibrator, etc.; a storage device 708 including, for example, a magnetic tape, a hard disk, etc.; and a communication device 709. The communication device 709 may allow the electronic device to communicate with other devices wirelessly or by wire to exchange data. Although Figure 7 The electronic device is shown with various devices, but it should be understood that it is not required to implement or possess all of the devices shown. More or fewer devices may be implemented or possessed instead.

[0131] In particular, according to an embodiment of the present disclosure, the process described above with reference to the flowchart can be implemented as a computer software program. For example, an embodiment of the present disclosure includes a computer program product, which includes a computer program carried on a non-transitory computer-readable medium, and the computer program includes a program code for executing the method shown in the flowchart. In such an embodiment, the computer program can be downloaded and installed from the network through the communication device 709, or installed from the storage device 708, or installed from the ROM 702. When the computer program is executed by the processing device 701, the above-mentioned functions defined in the method of the embodiment of the present disclosure are performed.

[0132] It should be noted that the computer-readable medium mentioned above in the present disclosure may be a computer-readable signal medium or a computer-readable storage medium, or any combination of the two. A computer-readable storage medium may be, for example, but not limited to, an electrical, magnetic, optical, electromagnetic, infrared, or semiconductor system, device, or component, or any combination of the above. More specific examples of computer-readable storage media may include, but are not limited to: an electrical connection with one or more wires, a portable computer disk, a hard disk, a random access memory (RAM), a read-only memory (ROM), an erasable programmable read-only memory (EPROM or flash memory), an optical fiber, a portable compact disk read-only memory (CD-ROM), an optical storage device, a magnetic storage device, or any suitable combination of the above. In the present disclosure, a computer-readable storage medium may be any tangible medium that contains or stores a program that can be used by or in conjunction with an instruction execution system, device, or component. In the present disclosure, a computer-readable signal medium may include a data signal propagated in baseband or as part of a carrier wave, which carries computer-readable program code. Such a propagated data signal may take a variety of forms, including but not limited to electromagnetic signals, optical signals, or any suitable combination of the above. A computer-readable signal medium may also be any computer-readable medium other than a computer-readable storage medium that can transmit, propagate, or transport a program for use by or in conjunction with an instruction execution system, apparatus, or device. The program code contained on the computer-readable medium may be transmitted using any suitable medium, including but not limited to wires, optical cables, RF (radio frequency), etc., or any suitable combination thereof.

[0133] In some embodiments, the client and server can communicate using any currently known or future developed network protocol, such as HTTP (HyperText Transfer Protocol), and can be interconnected with any form or medium of digital data communication (e.g., a communication network). Examples of communication networks include a local area network ("LAN"), a wide area network ("WAN"), an internet (e.g., the Internet), and a peer-to-peer network (e.g., an ad hoc peer-to-peer network), as well as any currently known or future developed network.

[0134] The computer-readable medium may be included in the electronic device, or may exist independently without being incorporated into the electronic device.

[0135] The computer-readable medium carries one or more programs. When the one or more programs are executed by the electronic device, the electronic device: generates a grid set of the calculation area according to geological modeling parameters and geological source data; constructs a linear equation system Ax=b, where A is the coefficient matrix of the grid set and the matrix vector b is used to represent point constraints and direction constraints; constructs a multigrid based on the grid set and determines the matrix coefficients of each grid layer; classifies the grid points in the multigrid into color class sets of different colors, where there is no mutual dependence between the grid points in each color class set; for each color class set, starts GPU threads with the same number of grid points as the color class set, where each GPU thread is used to process the Gauss-Seidel iterative calculation of the coefficient matrix row corresponding to one grid point; when all threads in a color class set complete the specified number of iterations or meet the convergence conditions, switch to the next color class set until the entire multigrid converges; and outputs the modeling results based on the converged multigrid.

[0136] Computer program code for performing the operations of the present disclosure may be written in one or more programming languages, or a combination thereof, including, but not limited to, object-oriented programming languages ​​such as Java, Smalltalk, C++, and conventional procedural programming languages ​​such as "C" or similar programming languages. The program code may be executed entirely on the user's computer, partially on the user's computer, as a stand-alone software package, partially on the user's computer and partially on a remote computer, or entirely on the remote computer or server. In cases involving a remote computer, the remote computer may be connected to the user's computer through any type of network, including a local area network (LAN) or a wide area network (WAN), or may be connected to an external computer (e.g., through the Internet using an Internet service provider).

[0137] The flowcharts and block diagrams in the accompanying drawings illustrate the possible implementation architecture, functions and operations of the systems, methods and computer program products according to various embodiments of the present disclosure. In this regard, each box in the flowchart or block diagram can represent a module, program segment, or a part of code, and the module, program segment, or a part of code contains one or more executable instructions for realizing the specified logical function. It should also be noted that in some alternative implementations, the functions marked in the box can also occur in a different order than that marked in the accompanying drawings. For example, two boxes represented in succession can actually be executed substantially in parallel, and they can sometimes be executed in the opposite order, depending on the functions involved. It should also be noted that each box in the block diagram and / or flowchart, and the combination of the boxes in the block diagram and / or flowchart, can be implemented with a dedicated hardware-based system that performs the specified function or operation, or can be implemented with a combination of dedicated hardware and computer instructions.

[0138] The units involved in the embodiments described in this disclosure may be implemented in software or hardware, wherein the name of a unit does not necessarily limit the unit itself.

[0139] The functions described above herein may be performed, at least in part, by one or more hardware logic components. For example, and without limitation, exemplary types of hardware logic components that may be used include: field programmable gate arrays (FPGAs), application specific integrated circuits (ASICs), application specific standard products (ASSPs), systems on chip (SOCs), complex programmable logic devices (CPLDs), and the like.

[0140] In the context of the present disclosure, a machine-readable medium can be a tangible medium that can contain or store a program for use by or in conjunction with an instruction execution system, device or equipment. A machine-readable medium can be a machine-readable signal medium or a machine-readable storage medium. A machine-readable medium can include, but is not limited to, an electronic, magnetic, optical, electromagnetic, infrared, or semiconductor system, device or equipment, or any suitable combination of the foregoing. A more specific example of a machine-readable storage medium can include an electrical connection based on one or more lines, a portable computer disk, a hard disk, a random access memory (RAM), a read-only memory (ROM), an erasable programmable read-only memory (EPROM or flash memory), an optical fiber, a portable compact disk read-only memory (CD-ROM), an optical storage device, a magnetic storage device, or any suitable combination of the foregoing.

[0141] The above description is merely a preferred embodiment of the present disclosure and an illustration of the technical principles employed. Those skilled in the art should understand that the scope of disclosure involved in the present disclosure is not limited to the technical solutions formed by the specific combination of the above-mentioned technical features, but also includes other technical solutions formed by any combination of the above-mentioned technical features or their equivalents without departing from the above-mentioned disclosed concepts. For example, a technical solution formed by replacing the above-mentioned features with (but not limited to) technical features with similar functions disclosed in this disclosure.

[0142] In addition, although each operation is described in a specific order, this should not be understood as requiring these operations to be performed in the specific order shown or in a sequential order. Under certain circumstances, multitasking and parallel processing may be advantageous. Similarly, although some specific implementation details have been included in the above discussion, these should not be interpreted as limiting the scope of the present disclosure. Some features described in the context of a separate embodiment can also be implemented in a single embodiment in combination. On the contrary, the various features described in the context of a single embodiment can also be implemented in multiple embodiments individually or in any suitable sub-combination mode.

[0143] Although the subject matter has been described in language specific to structural features and / or methodological logical acts, it should be understood that the subject matter defined in the appended claims is not necessarily limited to the specific features or acts described above. Rather, the specific features and acts described above are merely example forms of implementing the claims.

Claims

1. An implicit geological modeling method based on GPU acceleration, characterized in that: The method comprises: Generate a grid set of the calculation area based on geological modeling parameters and geological source data; Construct a linear equation system Ax=b, where A is the coefficient matrix of the grid set and the matrix vector b is used to represent the point constraints and direction constraints; constructing a multigrid based on the grid set and determining matrix coefficients of each grid layer; Classify the grid points in the multigrid into color class sets of different colors, wherein there is no mutual dependence between the grid points in each color class set; For each color class set, start the same number of GPU threads as the number of grid points in the color class set, where each GPU thread is used to process the Gauss-Seidel iterative calculation of the coefficient matrix row corresponding to one grid point; When all threads in a color class set complete the specified number of iterations or meet the convergence conditions, they switch to the next color class set until the entire multigrid converges; Output modeling results based on the converged multigrid state; Multiple GPU threads corresponding to a color class set are processed in parallel.

2. The method according to claim 1, characterized in that The method of classifying all grid points in the multigrid into color class sets of different colors includes: Based on vertex coloring and standard greedy coloring, all grid points in the multigrid are classified, where adjacent grid points belong to different color class sets.

3. The method according to claim 1, characterized in that The method further comprises: According to the preset judgment rules, the error solution of the multi-grid is divided into high-frequency part and low-frequency part; For the high frequency part, Jacobi iteration or Gauss-Seidel iteration is used for processing; The low-frequency part is processed using an iterative V cycle.

4. The method according to claim 3, characterized in that The iterative V cycle is used for processing, including: On the finest grid, the Gauss-Seidel iterative algorithm is used to solve the linear equations; The residual of the current grid is passed to the next coarse grid through the interpolation operator, and this process is repeated until the residual is limited to the coarsest grid; In the coarsest grid, the linear equations are accurately solved to obtain the calculation results; Return the solution of the coarsest grid to the previous fine grid through interpolation; Use the new solution to perform Gauss-Seidel iterations to obtain a better solution, and repeat this process until the solution is returned to the finest grid; Determine whether the termination condition is met. If the current result meets the preset termination error, the iteration ends; If the termination condition is not met, the Gauss-Seidel iterative calculation is continued with the current result as the initial value.

5. The method according to claim 1, wherein The constructing of a multigrid based on the grid set and determining matrix coefficients of each grid layer include: According to the coefficient matrix A, a coarse grid point set of the coefficient matrix A is determined, and a corresponding interpolation matrix P is determined; The matrix coefficients Ac of the coarse grid are determined based on the coefficient matrix A of the fine grid and the interpolation matrix P.

6. The method according to claim 5, characterized in that The step of determining a coarse grid point set of the coefficient matrix A and determining a corresponding interpolation matrix P according to the coefficient matrix A includes: Determine the strong connection points based on the algebraic properties of the diagonal elements and non-zero elements of each row of the fine grid coefficient matrix; Based on the strong connection points, the grid points are divided into independent point sets; After dividing all rows in the fine grid into their respective coarse grid point sets, and determining the intermediate transfer operator based on the clustering results; The intermediate transfer operator is processed based on the Jacobi smooth operator to obtain the interpolation matrix P.

7. The method according to claim 5, characterized in that The determination of the matrix coefficients Ac of the coarse grid according to the coefficient matrix A of the fine grid and the interpolation matrix P includes: The coefficient matrix corresponding to each grid layer is obtained by calculating the following formula: , Where, is the coefficient matrix corresponding to the coarse grid, A is the coefficient matrix corresponding to the previous fine grid.

8. The method according to claim 6, characterized in that The step of determining strong connection points based on the algebraic properties of the diagonal elements and non-zero elements of each row of the fine grid coefficient matrix includes: Use the following formula to determine the strong connection points corresponding to each row of the matrix. The non-zero elements that meet the following formula are marked as strong connection points of this row: , Where, , represents the current grid layer, Represents the coefficient matrix i The element in row j and column j.

9. An electronic device, characterized in that: include: one or more processors; a storage device for storing one or more programs, When the one or more programs are executed by the one or more processors, the one or more processors implement the method according to any one of claims 1 to 8.

Images