Generation of 3D Geological Model of Geographic Sites

Abstract

A computer-implemented method for generating a three-dimensional geological model of a subsurface volume of a geographical site. The method includes processing borehole data of a plurality of boreholes of the geographical site and prior geological data associated with the geographical site using a machine learning-based model to generate a plurality of geological profiles of the subsurface volume of the geographical site. The borehole data includes data related to subsurface soil and / or rock information along at least a portion of each of the plurality of boreholes. The prior geological data includes data related to prior stratigraphic information associated with the geographical site. Each geological profile contains corresponding stratigraphic information. The method further includes forming a three-dimensional geological model for display based on the plurality of generated geological profiles.

Description

Technical Field

[0001] The present invention relates to computer-implemented three-dimensional (3D) geological modeling techniques. Background Art

[0002] In earth science, a 3D digital twin of the subsurface space includes a virtual / digital representation of the subsurface stratigraphy of a geographic site in the real world. Such a subsurface 3D digital twin can be used, for example, to facilitate the design, construction, performance monitoring, and / or maintenance of the subsurface infrastructure of a geographic site.

[0003] The creation of a subsurface 3D digital twin, especially the 3D modeling of the subsurface stratigraphy of a geographic site, is often challenging. This is because typically only limited subsurface stratigraphy data of a geographic site is available (e.g., only a small portion of the subsurface space of the geographic site has been explored or surveyed), and there is a lack of a reliable tool that can accurately and efficiently determine or infer missing subsurface stratigraphy data from the limited available data. Summary of the Invention

[0004] The present invention provides a computer-implemented tool (system, method, etc.) for 3D geological modeling. The present invention can enable or facilitate efficient and accurate 3D geological modeling. The present invention is closely related to a computer / processor environment.

[0005] In a first aspect, there is provided a computer-implemented method for generating a 3D geological model of the subsurface volume of a geographic site. The method includes: using a machine learning-based model to process borehole data of a plurality of boreholes of a geographic site and prior geological data associated with the geographic site to generate a plurality of geological profiles of the subsurface volume of the geographic site. The borehole data includes data related to subsurface soil and / or rock information (e.g., types of subsurface soil and / or rock) along at least a portion of each of the plurality of boreholes. The prior geological data includes data related to prior stratigraphic information associated with the geographic site. Each geological profile contains corresponding stratigraphic information. The method further includes forming a 3D geological model for display based on the generated geological profiles. The prior geological data represents prior geological knowledge associated with the geographic site. The 3D geological model can correspond to a subsurface digital twin of the geographic site.

[0006] Optionally, the data related to prior stratigraphic information associated with the geographic site includes data related to prior stratigraphic information of another site near the geographic site and / or another site geologically similar to the geographic site.

[0007] Optionally, the data related to the previous stratigraphic information associated with the geographical site includes data related to the stratigraphic profile associated with the geographical site.

[0008] Optionally, the stratigraphic information of each geological profile includes the corresponding stratigraphic profile (e.g., along one or more corresponding directions).

[0009] Optionally, the data related to the previous stratigraphic information associated with the geographical site can be characterized as one or more (e.g., only one, only two, etc.) 2D images of the stratigraphic profile associated with the geographical site.

[0010] Optionally, the machine learning-based model is trained using the one or more 2D images.

[0011] Optionally, the machine learning-based model includes an artificial neural network, such as a convolutional neural network. Optionally, the machine learning-based model includes a multi-layer perceptron. Optionally, the machine learning-based model includes an algorithm based on gradient boosting (i.e., an algorithm based on gradient boosting techniques or a variant thereof).

[0012] Optionally, the processing includes: applying the borehole data and the prior geological data to a 3D coordinate system, and sequentially processing multiple slices for the subsurface volume in the 3D coordinate system using the machine learning-based model to continuously generate multiple geological profiles.

[0013] In some embodiments of the first aspect, the one or more 2D images include only a single 2D image (consisting of only a single 2D image). Optionally, in these embodiments, applying the borehole data and the prior geological data to the 3D coordinate system includes: arranging the borehole data in the 3D coordinate system along the depth (e.g., vertical) direction; arranging the single 2D image to extend along the depth (e.g., vertical) direction in the 3D coordinate system.

[0014] In some embodiments of the first aspect, the one or more 2D images include (or alternatively, only include) two 2D images. Optionally, in these embodiments, applying the borehole data and the prior geological data to the 3D coordinate system includes: arranging the borehole data in the 3D coordinate system along the depth (e.g., vertical) direction; arranging the two 2D images to extend along the depth (e.g., vertical) direction in the three-dimensional coordinate system and such that the corresponding imaginary plane containing the two 2D images defines a dihedral angle in the three-dimensional coordinate system. The dihedral angle is greater than 0 degrees and less than 180 degrees, e.g., approximately 90 degrees.

[0015] Optionally, the plurality of slices includes a first plurality of slices and a second plurality of slices. The first plurality of slices are arranged to be generally parallel to each other and generally parallel to one of the two 2D images. The second plurality of slices are arranged to be generally parallel to each other and generally parallel to the other of the two 2D images. The first plurality of slices are not parallel to the second plurality of slices. The first plurality of slices and the second plurality of slices may include the same number of slices or different numbers of slices.

[0016] Optionally, the processing further includes: identifying the slice containing the most borehole data. Optionally, the sequential processing includes: first processing the slice containing the most borehole data using a machine learning-based model; then, (i) if the first processed slice is one of the first plurality of slices, alternately processing one slice of the second plurality of slices and one slice of the first plurality of slices using the machine learning-based model until all of the plurality of slices are processed; and (ii) if the first processed slice is one of the second plurality of slices, alternately processing one slice of the first plurality of slices and one slice of the second plurality of slices using the machine learning-based model until all of the plurality of slices are processed.

[0017] Optionally, the geological profile generated from the processed slice is used as additional borehole data for the processing of the next slice.

[0018] Optionally, the sequential processing includes: first processing the slice containing the most borehole data; then processing the remaining slices in decreasing order of the amount of borehole data and additional borehole data contained in the slices until all of the plurality of slices are processed. Optionally, the sequential processing further includes: if two or more of the remaining slices contain the same amount of borehole data and additional borehole data, randomly selecting one of the two or more slices for processing.

[0019] Optionally, the processing of each slice of the plurality of slices includes, for each slice: if there is no borehole data along the boundary, performing an extrapolation operation using a machine learning-based model to determine the formation information of the boundary of the slice; and performing an interpolation operation using the machine learning-based model to determine the formation information of other parts of the slice where there is no borehole data.

[0020] Optionally, performing the extrapolation operation includes: determining a template based on the borehole data included in the slice, the template including a plurality of spaced unit columns parallel to the depth direction, the unit columns including a first column arranged along a set of borehole data, a second column arranged along another set of borehole data, and a third column arranged along a boundary without borehole data, which is used for extrapolation (i.e., extrapolating the borehole data to the third column); applying the template to a 2D image parallel to the slice to extract formation relationships and / or statistical data; training a machine learning-based model using the extracted formation relationships and / or statistical data; and applying the template to the slice along one or more paths to iteratively process the slice using the trained machine learning-based model to determine the formation information of the boundary.

[0021] Optionally, performing the interpolation operation includes: determining a plurality of templates based on the borehole data included in the slice, each of the plurality of templates including a corresponding plurality of spaced unit columns parallel to the depth direction, the unit columns including a corresponding first column arranged along a corresponding set of borehole data, a corresponding second column arranged along a corresponding another set of borehole data, and a third column in the space between the first column and the second column and without borehole data, which is used for interpolation (i.e., interpolating the borehole data to the third column); applying the template to a 2D image parallel to the slice to extract formation relationships and / or statistical data; training a machine learning-based model using the extracted formation relationships and / or statistical data; and applying the template to the slice along one or more paths to iteratively process the slice using the trained machine learning-based model to determine the formation information of other parts of the slice without borehole data.

[0022] Optionally, forming the 3D geological model includes combining the generated geological profiles.

[0023] Optionally, the 3D geological model has a 3 spatial resolution of 1 m or less.

[0024] Optionally, the 3D geological model includes a 3D point cloud, and each point in the 3D point cloud contains subsurface soil and / or rock information (e.g., the type of subsurface soil and / or rock) modeled as a categorical variable. Each point in the 3D point cloud may correspond to a voxel.

[0025] Optionally, the computer-implemented method further includes quantifying the formation uncertainty associated with the 3D model.

[0026] Optionally, quantifying the formation uncertainty associated with the 3D model includes: for each point in the 3D point cloud, determining a value representing the confidence level associated with the subsurface soil and / or rock information contained in the point.

[0027] Optionally, quantifying formation uncertainty associated with a 3D model includes performing a statistical analysis on a plurality of equally probable 3D geological realizations or random samples generated based on Monte Carlo simulations.

[0028] Optionally, the computer-implemented method further includes displaying a map of the 3D geological model and / or the quantified formation uncertainty.

[0029] Optionally, the computer-implemented method further includes: obtaining new borehole data of a plurality of boreholes at a geographical site and / or new prior geological data associated with the geographical site; after obtaining the new borehole data and / or the new prior geological data, automatically performing the following operations: using a machine learning-based model to process at least one of the new borehole data and / or the new prior geological data, and the borehole data and the prior geological data, to generate a plurality of updated geological profiles of the subsurface volume, each of which contains corresponding formation information (e.g., patterns) associated with the geographical site. The computer-implemented method further includes forming an updated 3D geological model based on the generated updated geological profiles.

[0030] Optionally, the computer-implemented method further includes: detecting the availability of new borehole data and / or new prior geological data; and automatically obtaining the new borehole data and / or the new prior geological data when the availability of the new borehole data and / or the new prior geological data is detected. Optionally, the detection is performed according to a predetermined schedule (e.g., periodically). Optionally, the detection is continuous.

[0031] In a second aspect, there is provided a computer-readable medium, such as a non-transitory computer-readable storage medium, storing one or more programs (computer programs) configured to be executed by one or more processors, the one or more programs including instructions for performing or facilitating the performance of the method of the first aspect.

[0032] In a third aspect, there is provided a system for generating a 3D geological model of a subsurface volume of a geographical site. The system includes one or more processors and a memory. The memory stores one or more programs (computer programs) configured to be executed by one or more processors. The one or more programs include instructions for performing or facilitating the performance of the method of the first aspect. The system may further include one or more displays operably coupled to the one or more processors for displaying a map of the 3D geological model and / or the quantified formation uncertainty.

[0033] In a fourth aspect, there is provided a computer-implemented method for operating a computer-implemented geographic information platform. The method includes: receiving a request to generate a 3D geological model of a geographical site; in response to the request, performing the method of the first aspect; and outputting a 3D geological model of the subsurface volume of the geographical site for display.

[0034] In a fifth aspect, there is provided a computer-readable medium, such as a non-transitory computer-readable storage medium, storing one or more programs (computer programs) configured to be executed by one or more processors, the one or more programs including instructions for performing or facilitating the performance of the method of the fourth aspect.

[0035] In a sixth aspect, there is provided a computer-implemented geographic information platform, comprising: one or more processors and a memory. The memory stores one or more programs (computer programs) configured to be executed by one or more processors. The one or more programs include instructions for performing or facilitating the performance of the method of the fourth aspect. The computer-implemented geographic information platform can be implemented using hardware and software.

[0036] In a seventh aspect, there is provided a system, particularly a computer system, having suitable means for performing or facilitating the performance of the method of the first aspect and / or the fourth aspect.

[0037] By considering the detailed description and the drawings, other features and aspects of the present invention will become apparent. Where appropriate and applicable, any one or more features described herein with respect to one aspect or embodiment can be combined with any one or more other features described herein with respect to any one or more other aspects or embodiments.

[0038] Degree terms or terms of correlation related to quantity or condition (such as "generally", "about", "approximately", "substantially", etc.) are used to account for (depending on the context) at least one of the following: manufacturing tolerances, degradation, assembly, use, trends, tendencies, imperfect realities, etc. BRIEF DESCRIPTION OF THE DRAWINGS

[0039] Embodiments of the present invention will now be described by way of example with reference to the drawings, in which:

[0040] Figure 1 is a schematic diagram showing a two-dimensional (2D) simulation-based framework for generating a 3D geological model of a subsurface volume of a geographic site according to an embodiment of the present invention;

[0041] Figure 2A is a schematic diagram showing an arrangement of two 2D training images according to an embodiment of the present invention, the arrangement being used in determining the position and order of slices for 2D simulation;

[0042] Figure 2B is a schematic diagram showing an arrangement of site-specific measurements (drilling data) according to an embodiment of the present invention, the arrangement being used in determining the position and order of slices for 2D simulation;

[0043] Figure 3AIt is a schematic diagram showing the extraction of large-scale formation features and corresponding spatial predictions using a modified CNN model in an XGBoost-based algorithm according to an embodiment of the present invention;

[0044] Figure 3B It is a schematic diagram showing the comprehensive prediction operation in an XGBoost-based algorithm according to an embodiment of the present invention;

[0045] Figure 4 It is a schematic diagram showing the delimitation of a single simulation slice using an XGBoost-based algorithm according to an embodiment of the present invention;

[0046] Figure 5A It is a diagram of a simulated 3D training geological domain (volume) according to an embodiment of the present invention;

[0047] Figure 5B It is a diagram of a simulated 3D test geological domain (volume) according to an embodiment of the present invention;

[0048] Figure 5C It is Figure 5B a perspective view of the simulated 3D test geological domain (volume);

[0049] Figure 5D It is a diagram showing two substantially vertical training images extracted from Figure 5A the simulated 3D training geological domain;

[0050] Figure 5E It is a diagram showing Figure 5B 25 site-specific boreholes extracted from the simulated 3D test geological domain;

[0051] Figure 6 It is a graph showing the relationship between the percentage change in the most likely prediction and the number of 3D realizations or random samples (up to 100) according to an embodiment of the present invention;

[0052] Figure 7A It is a diagram showing the most likely 3D geological model derived from 100 3D realizations according to an embodiment of the present invention;

[0053] Figure 7B It is a diagram showing the dispersion of 100 3D realizations regarding the most likely prediction according to an embodiment of the present invention;

[0054] Figure 8A It is a graph showing the relationship between the prediction accuracy (compared with the real geological domain) of the most likely prediction for different soil types and the number of 3D realizations (up to 100) according to an embodiment of the present invention;

[0055] Figure 8Bis a diagram showing the accuracy of the determined 3D geological model (compared to the true geological domain);

[0056] Figure 9A is a graph showing the most likely prediction of a selected slice (X = 0) in the Figure 7A most likely 3D geological model;

[0057] Figure 9B is a graph showing the dispersion diagram of a selected slice (X = 0) in the Figure 7A most likely 3D geological model;

[0058] Figure 9C is a graph showing the accuracy diagram of a selected slice (X = 0) in the Figure 7A most likely 3D geological model;

[0059] Figure 9D is a graph showing the most likely prediction of a selected slice (X = 10) in the Figure 7A most likely 3D geological model;

[0060] Figure 9E is a graph showing the dispersion diagram of a selected slice (X = 10) in the Figure 7A most likely 3D geological model;

[0061] Figure 9F is a graph showing the accuracy diagram of a selected slice (X = 10) in the Figure 7A most likely 3D geological model;

[0062] Figure 9G is a graph showing the most likely prediction of a selected slice (X = 40) in the Figure 7A most likely 3D geological model;

[0063] Figure 9H is a graph showing the dispersion diagram of a selected slice (X = 40) in the Figure 7A most likely 3D geological model;

[0064] Figure 9I is a graph showing the accuracy diagram of a selected slice (X = 40) in the Figure 7A most likely 3D geological model;

[0065] Figure 10A is a graph showing the most likely prediction of a selected slice (Y = 60) in the Figure 7A most likely 3D geological model;

[0066] Figure 10B is a graph showing the dispersion diagram of a selected slice (Y = 60) in the Figure 7A most likely 3D geological model;

[0067] Figure 10C is a diagram showing the accuracy map of a selected slice (Y = 60) in the most likely 3D geological model of Figure 7A ;

[0068] Figure 10D is a diagram showing the most likely prediction of a selected slice (Y = 90) in the most likely 3D geological model of Figure 7A ;

[0069] Figure 10E is a diagram showing the dispersion map of a selected slice (Y = 90) in the most likely 3D geological model of Figure 7A ;

[0070] Figure 10F is a diagram showing the accuracy map of a selected slice (Y = 90) in the most likely 3D geological model of Figure 7A ;

[0071] Figure 10G is a diagram showing the most likely prediction of a selected slice (Y = 100) in the most likely 3D geological model of Figure 7A ;

[0072] Figure 10H is a diagram showing the dispersion map of a selected slice (Y = 100) in the most likely 3D geological model of Figure 7A ;

[0073] Figure 10I is a diagram showing the accuracy map of a selected slice (Y = 100) in the most likely 3D geological model of Figure 7A ;

[0074] Figure 11A is a diagram showing 25 unevenly arranged boreholes in an embodiment of the present invention;

[0075] Figure 11B is a diagram showing a single training image extracted from the simulated 3D training geological domain of Figure 5A in an embodiment of the present invention;

[0076] Figure 12A is a diagram showing the most likely prediction of three selected slices (X = 20, Y = 20, and Z = 5) when the boreholes have equal spacing in an embodiment of the present invention;

[0077] Figure 12B is a diagram showing the dispersion of three selected slices (X = 20, Y = 20, and Z = 5) when the boreholes have equal spacing in an embodiment of the present invention;

[0078] Figure 12CA graph showing the accuracy of three selected slices (X = 20, Y = 20, and Z = 5) when the drill holes have equal spacing in one embodiment of the present invention:

[0079] Figure 12D A graph showing the most likely prediction of three selected slices (X = 20, Y = 20, and Z = 5) when the drill holes have unequal spacing in one embodiment of the present invention;

[0080] Figure 12E A graph showing the dispersion of three selected slices (X = 20, Y = 20, and Z = 5) when the drill holes have unequal spacing in one embodiment of the present invention;

[0081] Figure 12F A graph showing the accuracy of three selected slices (X = 20, Y = 20, and Z = 5) when the drill holes have unequal spacing in one embodiment of the present invention:

[0082] Figure 13A A graph showing the most likely prediction of three selected slices (X = 0, Y = 0, and Z = 0) when using two training images (X = 50, Y = 50) in one embodiment of the present invention;

[0083] Figure 13B A graph showing the dispersion of three selected slices (X = 0, Y = 0, and Z = 0) when using two training images (X = 50, Y = 50) in one embodiment of the present invention;

[0084] Figure 13C A graph showing the accuracy of three selected slices (X = 0, Y = 0, and Z = 0) when using two training images (X = 50, Y = 50) in one embodiment of the present invention;

[0085] Figure 13D A graph showing the most likely prediction of three selected slices (X = 0, Y = 0, and Z = 0) when using one training image (X = 50) in one embodiment of the present invention;

[0086] Figure 13E A graph showing the dispersion of three selected slices (X = 0, Y = 0, and Z = 0) when using one training image (X = 50) in one embodiment of the present invention;

[0087] Figure 13F A graph showing the accuracy of three selected slices (X = 0, Y = 0, and Z = 0) when using one training image (X = 50) in one embodiment of the present invention;

[0088] Figure 14AA chart showing the most likely prediction of a selected slice (Y = 60) in the most likely 3D geological model obtained using a single training image (X = 50) in one embodiment of the present invention;

[0089] Figure 14B A chart showing the dispersion map of a selected slice (Y = 60) in the most likely 3D geological model obtained using a single training image (X = 50) in one embodiment of the present invention;

[0090] Figure 14C A chart showing the accuracy map of a selected slice (Y = 60) of the most likely 3D geological model obtained using a single training image (X = 50) in one embodiment of the present invention;

[0091] Figure 14D A chart showing the most likely prediction of a selected slice (Y = 90) in the most likely 3D geological model obtained using a single training image (X = 50) in one embodiment of the present invention;

[0092] Figure 14E A chart showing the dispersion map of a selected slice (Y = 90) in the most likely 3D geological model obtained using a single training image (X = 50) in one embodiment of the present invention;

[0093] Figure 14F A chart showing the accuracy map of a selected slice (Y = 90) of the most likely 3D geological model obtained using a single training image (X = 50) in one embodiment of the present invention;

[0094] Figure 14G A chart showing the most likely prediction of a selected slice (Y = 100) in the most likely 3D geological model obtained using a single training image (X = 50) in one embodiment of the present invention;

[0095] Figure 14H A chart showing the dispersion map of a selected slice (Y = 100) in the most likely 3D geological model obtained using a single training image (X = 50) in one embodiment of the present invention;

[0096] Figure 14I A chart showing the accuracy map of a selected slice (Y = 100) of the most likely 3D geological model obtained using a single training image (X = 50) in one embodiment of the present invention;

[0097] Figure 15 A functional block diagram of an example information processing system that can be used to execute one or more method embodiments of the present invention;

[0098] Figure 16is a flowchart showing a method for generating a 3D geological model of the subsurface volume of a geographical site according to an embodiment of the present invention;

[0099] Figure 17 is a flowchart showing a method for operating a computer-implemented geographical information platform according to an embodiment of the present invention; and

[0100] Figure 18 is an example screenshot of a user interface of a computer-implemented geographical information platform according to an embodiment of the present invention. Detailed Embodiments

[0101] Figure 16 shows a method 1600 for generating a 3D geological model of the subsurface volume of a geographical site according to some embodiments of the present invention. The method 1600 is a computer-implemented method that operates using computer hardware and / or software.

[0102] The method 1600 begins at step 1602, where borehole data and prior geological data are processed using a machine learning-based model to generate a geological profile of the subsurface volume of the geographical site. The borehole data is obtained from multiple boreholes at the geographical site. The borehole data includes data related to subsurface soil and / or rock information (e.g., the type of subsurface soil and / or rock) along at least a portion of each borehole. The prior geological data relates to data reflecting prior (known) geological knowledge associated with the geographical site, e.g., learned or determined by a geotechnical engineer. The prior geological data includes data related to prior stratigraphic information (e.g., a stratigraphic section) associated with the geographical site. In one example, the prior stratigraphic information associated with the geographical site is the prior stratigraphic information of the geographical site. In another example, the prior stratigraphic information associated with the geographical site is the prior stratigraphic information of another geographical site near the geographical site and / or another geographical site that is geologically similar to the geographical site (e.g., when no local reliable data is available for the geographical site to be modeled). In some embodiments, the data related to the stratigraphic information associated with the geographical site can be characterized as a two-dimensional (2D) image of a stratigraphic section associated with the geographical site. Each generated geological profile contains corresponding stratigraphic information, e.g., a stratigraphic section.

[0103] In one example, the processing in step 1602 includes applying the borehole data and the prior geological data to a 3D coordinate system. This helps to orient the data to the same reference frame for further processing. After applying the borehole data and the prior geological data to the 3D coordinate system, the processing includes sequentially processing multiple slices of the subsurface volume in the 3D coordinate system using the machine learning-based model to continuously (i.e., one after another) generate geological profiles.

[0104] In one example, data related to formation information associated with a geographical site can be represented by a single 2D image. In this example, applying borehole data and prior geological data to a 3D coordinate system includes arranging the borehole data in the 3D coordinate system along a depth (e.g., vertical) direction and arranging the single 2D image in the 3D coordinate system along a depth (e.g., vertical) direction. The single 2D image also extends along another direction in the 3D coordinate system.

[0105] In another example, data related to formation information associated with a geographical site can be represented by two 2D images. In this example, applying borehole data and prior geological data to a 3D coordinate system includes arranging the borehole data in the 3D coordinate system along a depth (e.g., vertical) direction and extending the two 2D images respectively in the 3D coordinate system along a depth (e.g., vertical) direction. The two 2D images also each extend along another direction in the 3D coordinate system such that the two 2D images define two imaginary planes. These imaginary planes can define a dihedral angle in the 3D coordinate system. The dihedral angle is greater than 0 degrees and less than 180 degrees, e.g., between about 45 degrees and about 135 degrees, at about 90 degrees, etc. In the sequential processing of the plurality of slices described above, the slices can include a plurality of slices arranged generally parallel to each other and parallel to one of the 2D images, and a plurality of slices arranged generally parallel to each other and parallel to the other 2D image. The slices arranged parallel to one of the 2D images and the slices arranged parallel to the other 2D image are not parallel to each other and can have different numbers or the same number of slices.

[0106] In some embodiments, the processing further includes identifying, from the slices, the slice that contains the largest amount of borehole data. The sequential processing can include first processing, using a machine learning-based model, the slice that contains the largest amount of borehole data, and then, (i) if the first processed slice is a slice parallel to one of the 2D images, alternately processing, using the machine learning-based model, the slices parallel to the other 2D image and the slices parallel to the one 2D image until all slices are processed, and (ii) if the first processed slice is a slice parallel to the other 2D image, alternately processing, using the machine learning-based model, the slices parallel to one 2D image and the slices parallel to the other 2D image until all slices are processed.

[0107] The geological profile generated from the processed (already processed) slice can be used as additional borehole data for the processing of the next slice (slice to be processed). In some embodiments, sequential processing may include first processing the slice containing the largest amount of borehole data, and then processing the remaining slices in decreasing order of the amount of borehole data and additional borehole data contained in the slice until all slices are processed. During this process, if the remaining slices contain the same amount of borehole data and additional borehole data, one of them is selected (e.g., randomly) for processing.

[0108] For each slice, the processing may include: if there is no borehole data along the slice boundary, performing an extrapolation operation using a machine learning-based model to determine the formation information of the slice boundary, and performing an interpolation operation using a machine learning-based model to determine the formation information of other parts of the slice without borehole data. As an example, the extrapolation operation includes determining a template (processing template) based on the borehole data contained in the slice. The template includes a plurality of spaced unit columns parallel to the depth direction, where the first column is arranged along a set of borehole data, the second column is arranged along another set of borehole data, and the third column is arranged along the boundary without borehole data for extrapolation. In this example, the extrapolation operation further includes applying the template to a 2D image parallel to the slice to extract formation relationships and / or statistical data, training a machine learning-based model using the extracted formation relationships and / or statistical data, and applying the template to the slice along one or more paths to iteratively process the slice using the trained machine learning-based model to determine the formation information of the boundary. As an example, the interpolation operation includes determining a plurality of templates based on the borehole data contained in the slice. Each template includes corresponding spaced unit columns parallel to the depth direction, which include a corresponding first column along a corresponding set of borehole data, a corresponding second column along another corresponding set of borehole data, and a third column in the space between the first column and the second column and without borehole data for interpolation. In this example, the interpolation operation further includes applying the templates to a 2D image parallel to the slice to extract formation relationships and / or statistical data, training a machine learning-based model using the extracted formation relationships and / or statistical data, and applying the templates to the slice along one or more paths to iteratively process the slice using the trained machine learning-based model to determine the formation information of other parts of the slice without borehole data.

[0109] In step 1602, a machine learning-based model can be trained using a 2D image. The machine learning-based model can be, for example, an artificial neural network (e.g., a convolutional neural network), a multi-layer perceptron, an algorithm based on gradient boosting (i.e., an algorithm based on gradient boosting techniques or its variants), etc.

[0110] After step 1602, method 1600 proceeds to step 1604 and step 1606. In step 1604, a 3D geological model of the subsurface volume of the geographical site is formed based on the generated geological profile. In step 1606, the 3D geological model is displayed. Step 1604 and step 1606 can be executed separately (formed then displayed) or simultaneously (displayed while forming). The 3D geological model can be formed by combining the generated geological profiles. In some embodiments, the 3D geological model has a spatial resolution of 1m 3 or less (relative to the geographical site). The formed 3D geological model can be a 3D point cloud with multiple points (e.g., voxels), where each point contains subsurface soil and / or rock information (e.g., the type of subsurface soil and / or rock) modeled as a categorical variable. The 3D geological model of the subsurface volume of the geographical site can correspond to the subsurface digital twin of the geographical site.

[0111] In some embodiments, method 1600 further includes quantifying the formation uncertainty associated with the 3D geological model. For example, the quantification can include determining a value for each point in the 3D point cloud that indicates the confidence level associated with the subsurface soil and / or rock information contained in that point. The quantification can include performing a statistical analysis on a plurality of equally probable 3D geological realizations or random samples generated based on Monte Carlo simulations. A map of the quantified formation uncertainty can be displayed, e.g., separately from or side by side with the 3D geological model.

[0112] The following disclosure (refer to Figures 1 to 15 ) provides some example embodiments of method 1600.

[0113] The inventors of the present invention have developed a data-driven machine learning algorithm called iterative convolutional XGBoost, which is used to delineate 2D subsurface geological profiles, also known as "2D IC-XGBoost", and is disclosed in the literature of Shi et al., "Development of subsurface geological cross-section from limited site-specific boreholes and prior geological knowledge using iterative convolution, Journal of Geotechnical and Geoenvironmental Engineering, 147(9), 04021082", the entire content of which is incorporated herein by reference. The 2D IC-XGBoost algorithm combines a training image with a convolutional neural network (CNN) framework to adaptively extract formation connectivity at different scales from a single training image.

[0114] Through further research, experiments, and tests, the inventors of the present invention have now designed that 2D IC-XGBoost is efficient and can be modified and improved for generating 3D subsurface geological models and subsurface digital twins.

[0115] The following disclosure relates to a data-driven formation modeling and uncertainty quantification technique for spatially delineating a 3D subsurface formation from limited site-specific boreholes and prior geological knowledge (or prior geological data). In some embodiments, the technique relies on sequential simulation of 2D slices, performed using an enhanced IC-XGBoost algorithm called "IC-XGBoost3D" (modified from the aforementioned 2D IC-XGBoost). In some embodiments, the technique employs a sequential simulation strategy that avoids the need for a complete 3D training image for stochastic modeling and addresses the challenge of the lack of 3D training images in engineering practice. In some embodiments, the reliable formation statistics of a geographical site can be represented by a 2D training image, which can reveal potential anisotropy. In some embodiments, the technique utilizes sequential prediction at unsampled locations to generate a quasi-3D representation of the subsurface geology.

[0116] 3D Subsurface Geological Modeling Method

[0117] The inventors of the present invention have previously recognized that machine learning-based algorithms such as IC-XGBoost can be used to delineate 2D subsurface geological profiles. The basic principle is to combine site-specific measurements with 2D formation relationships learned from a single training image for conditional simulation. In this type of algorithm, the training image can be regarded as a collection of prior geological knowledge of the surveyed site.

[0118] The inventors of the present invention have now thought that, based on the same principle, theoretically, a 3D geological domain can be modeled from limited site-specific measurements and based on a complete 3D training image. However, the inventors have found that the problem is that in practical applications, it is difficult, if not impossible, to obtain a qualified 3D training image, and site-specific measurements are usually too sparse (insufficient) to directly construct a detailed 3D geological model with the required accuracy. Although ideally a 3D training image should be used to represent and capture 3D geological complexity, in geotechnical practice, 3D training images are usually non-existent or difficult to obtain because they require explicit information (usually non-existent) about the potential 3D spatial variation of soil types. On the other hand, a common practice in geotechnical engineering is to design based on key or representative 2D profiles. These 2D profiles can provide a rich database for the selection of training images. The inventors of the present invention have therefore designed that, as an alternative to using 3D training images, 2D training images can be used to construct a 3D subsurface geological domain from sparse site-specific data for geotechnical practice.

[0119] The following disclosure provides a sequential simulation method of one embodiment, which can be used to solve the difficulties associated with the lack of 3D training data and sparse site-specific measurements (i.e., limited real data) in 3D geological modeling.

[0120] Figure 1 A framework 100 for shaping or generating a 3D subsurface geological model (i.e., delineating a 3D subsurface geological domain) of a geographical site based on a 2D training image according to one embodiment of the present invention is shown. The 3D geological domain of interest can be visualized as a 3D point cloud.

[0121] In this embodiment, in 102, the training image and site-specific measurements are shown. The training image is two 2D training images that are generally arranged perpendicular to each other. The site-specific measurements are measurements / data specific to the geographical site, which in this example include borehole data obtained from boreholes BH1 and BH2 at the site. In 104, the simulation order is determined. Specifically, the entire 3D domain can be split into a series of 2D vertical slices along two perpendicular directions. A random simulation path can be determined to access each simulation slice based on the general principle of first accessing or simulating the slice with the most site-specific data. In 106, the 2D profiles of the slices are sequentially predicted or determined. Specifically, the 2D simulation slices (i.e., geological profiles) are accessed in sequence, and the relevant stratigraphic changes are delineated using a data-driven algorithm (described in further detail below). Then, each simulated geological profile is added back to the geological domain and used as a measurement for subsequent conditional simulations. In 108, the 3D realization is reconstructed. Specifically, after shaping the geological profiles along all 2D slices, the complete 3D realization is reconstructed by combining all the simulated 2D slices.

[0122] Since a single 2D training image only contains spatial statistics in one direction, in order to reveal the potential formation anisotropy in different directions, a combination of multiple 2D training images in different directions is required. In this example, the potential formation connectivity of each 2D slice is learned from the training images parallel to the relevant 2D slice. A pair of 2D perpendicular training images can be used to represent the anisotropic formation connectivity in different directions in the 3D geological domain. In one example, if the geology of the site of interest is isotropic, a single training image can also be used to delineate the formation connectivity in different directions and generate an isotropic 3D geological domain from the single training image. In practice, the determination of isotropic and anisotropic geological features depends on the subjective judgment of the engineer.

[0123] In this example, to address the issues related to sparse or limited site-specific measurements, a sequential simulation program is applied to continuously model a series of 2D vertical slices. The 2D slice from the previous simulation step is used as additional site-specific measurement for the conditional simulation of the subsequent 2D slice. In other words, at each simulation step, the training image parallel to the simulated slice is used and combined with the site-specific measurements (including the original measurements and the additional data obtained from the previous simulation) for stochastic simulation. Using this method, the modeling of the 3D geological domain can be transformed into a sequential simulation of a series of 2D vertical slices, which can be achieved through machine learning of 2D training images and sparse site-specific measurements. As the spatial prediction progresses, the formation connectivity extracted from the 2D vertical training images is gradually added to the 3D point cloud, thus replicating the 3D spatial anisotropy pattern. Therefore, the difficulties in the digital representation of the 3D geological domain due to the lack of 3D training images and sparse site-specific measurements can be solved by the method of this embodiment (utilizing a pair of 2D training images and the sequential simulation process).

[0124] In this example, the existing 2D IC-XGBoost algorithm is modified and improved to construct a 3D geological domain by learning the formation profiles from a pair of training images following the sequential simulation order. Additionally, in this example, the existing 2D IC-XGBoost is improved by combining extrapolation and pre-training strategies, which aim to improve the relevant computational efficiency and thus facilitate the generation of the 3D geological model. The main parts of the framework for 3D geological modeling in this embodiment are further described below.

[0125] One of the main parts of the framework is the training image.

[0126] The method in this embodiment constructs a 3D subsurface geological domain conditioned on site-specific measurements and a pair of training images, which are conceptual statistical models of the geology that are considered (credible) to exist in the region of interest (the site). The training images provide a quantitative framework that can be used to perform geological modeling using valuable prior geological knowledge.

[0127] In Figure 1 , 102 shows an example of two vertical training images and two site-specific boreholes (BH1 and BH2). In this example, the 3D subsurface geological model includes four soil types. The two training images TI XZ and TI YZ reflect the trust or credible spatial variability of the subsurface stratigraphic profiles in two vertical directions, i.e., the X direction and the Y direction. Each training image depicts the geometric features of each soil layer expected to exist in the formation and the complex stratigraphic connectivity between different soil layers in the same direction. The combination of the two training images depicts horizontal anisotropy and facilitates a quasi-3D representation of the subsurface strata. In this embodiment, the two training images are arranged vertically. However, it should be noted that in some other embodiments, the two training images do not have to be perpendicular to each other and can intersect at an angle other than 90 degrees (greater than 0 degrees and less than 180 degrees). Without loss of generality, in the following disclosure, two vertical training images are used unless otherwise stated.

[0128] As described above, in some embodiments, when the geological domain is isotropic and a single training image adopted can describe the formation characteristics in other directions, a single training image can be used for spatial modeling. In this case, the term isotropic means that the stratigraphic connectivity reflected in different directions (e.g., the X direction and the Y direction) can be correctly reflected by or correspond to a single set of training images. In earth sciences, potential training images can be generated directly from outcrop structures or indirectly using object-based and process-based techniques. For geological and geotechnical applications, the shaping of 2D profiles is a prerequisite for actual projects. Therefore, there is a rich database of training images. In addition, training images can be easily borrowed from nearby sites or previous projects with similar geological environments. Different techniques can be used to evaluate the compatibility of potential training images with site-specific data (such as boreholes), and these techniques include, for example, running distributions, multi-point density functions (MPDF), or stratigraphic edge direction histograms. The application of training images in spatial modeling is based on the assumption that stratigraphic profiles will repeat in regions with similar geological origins.

[0129] Another major part of the framework is related to determining the simulation order and profile positions.

[0130] In this embodiment, a 3D geological domain is constructed by sequentially generating 2D geological profiles and recombining them into a 3D point cloud. The order and position of the 2D simulation slices affect the accuracy of the 3D geological domain modeling. Essentially, the simulation order is established based on the principle that each 2D simulation slice contains as much conditional data as possible. It should be noted that the conditional data includes the original site-specific boreholes and the results obtained from previous simulation steps.

[0131] Figure 2A and Figure 2B Generally illustrates the process of determining the order and position of 2D simulation slices in one embodiment, where Figure 2A shows the layout of the training image, while Figure 2B shows the layout of the site-specific measurements. In this embodiment, a 3D point cloud of dimension N X ×N Y ×N Z is projected onto a 2D horizontal layout map (Z = 0), and the vertical 2D profiles are represented by lines in Figure 2A and Figure 2B . In Figure 2A , two vertical training images (TI YZ and TI XZ ) are represented by thick lines. In Figure 2B , four site-specific boreholes (i.e., BH1, BH2, BH3, and BH4) are represented by solid squares.

[0132] The simulation order starts from the simulation slice (SS) with the most site-specific boreholes. For example, there are two candidate simulation slices, namely SS XZ,Y＝a and SS YZ,X＝b , where a and b represent the coordinates along the Y and X directions respectively. SS XZ,Y＝a contains three boreholes. SS YZ,X＝b contains two boreholes. Since SS XZ,Y＝a has more site-specific measurements, it ranks first in the simulation order. Mathematically, the above process for determining the simulation slice order can be described by the following equation.

[0133]

[0134] where N BH [·] is the total number of boreholes within a given simulation slice. N BH [SS YZ,X＝b represents the total number of boreholes within the profile SS YZ,X＝b (including the original boreholes and the results obtained from previous simulation steps), which is a profile parallel to the YZ plane and has an X value of b. Similarly, N BH [SS XY,Y＝a represents the profile SSXY,Y＝a The total number of drill holes within XY,Y＝a . If a set (two or more) of 2D simulation slices has the same number of measurements, one 2D simulation slice is randomly selected from the set for spatial prediction. The remaining 2D simulation slices in the set are simulated in alternating directions until all 2D simulation slices have been simulated. In other words, if a randomly selected simulation slice parallel to the X direction is predicted first, the next simulation will move to a randomly selected simulation slice parallel to the Y direction, and then back to another randomly selected simulation slice parallel to the X direction. Similarly, if a randomly selected simulation slice parallel to the Y direction is predicted first, the next simulation will move to a randomly selected simulation slice parallel to the X direction, and then back to another randomly selected simulation slice parallel to the Y direction. After identifying the first simulation slice in the simulation order, spatial prediction of the unsampled profiles within the simulation slice is performed using the algorithm described in more detail below.

[0135] In this embodiment, only the training images parallel to the simulation slices are used for spatial prediction of the unsampled regions in the simulation slices. After shaping the geological profile simulation slices, the entire simulation slice is considered a measurement (i.e., a set of additional drill holes) and is added back Figure 2B to the layout map in Figure 2B to determine the next simulation slice. For example, if the size of the simulation slice is 50 (depth) × 100 (horizontal) and the resolution is 2m, after adding the shaped slice (i.e., SS) back Figure 2B to the layout map in Figure 2B , the total number of additional drill holes is 50 (= 100 / 2). Thus, the adjusted total number of drill holes for the next simulation slice includes the site-specific drill holes and the additional drill holes calculated from the previous step. The simulation slices with the maximum N BH along the X direction and along the Y direction are determined and compared. The simulation slice with the most drill holes (and thus the most drill hole data) is designated as the next simulation slice. If the number of the most drill holes (and thus the most drill hole data) is the same in both directions, a random draw is performed to select the next simulation slice. The above steps are repeated until all 2D slices within the 3D geological domain have been simulated.

[0136] Another major part of the framework involves spatial prediction of the 2D simulation slices using an XGBoost-based algorithm.

[0137] In this embodiment, after determining the order and position of the 2D simulation slices, an enhanced IC-XGBoost machine learning algorithm is used for spatial prediction with the drilled and parallel training images as conditions. The enhanced IC-XGBoost machine learning algorithm is an improvement of the 2D IC-XGBoost in the literature publicly disclosed by Shi et al. and incorporated herein by reference. 2D IC-XGBoost can be regarded as a regression algorithm, which is used to determine the unknown soil types between known site-specific measurements. The basic mathematics of the 2D IC-XGBoost algorithm can be found in the literature publicly disclosed by Shi et al. and incorporated herein by reference. For the sake of brevity, only some main parts of the 2D IC-XGBoost algorithm are described below.

[0138] Figure 3A and Figure 3B shows some main features of the 2D IC-XGBoost algorithm, where Figure 3A shows the extraction of large-scale formation features and corresponding spatial prediction using a modified CNN model Figure 3B shows the iterative prediction operation.

[0139] 2D IC-XGBoost is rooted in the framework of a convolutional neural network (CNN). It uses multi-scale templates to extract potential formation relationships from a single training image through iterative convolution. The extracted formation profiles can be used as input for training an XGBoost model, which is used for subsequent spatial prediction.

[0140] Figure 3A shows the extraction of large-scale formation features and corresponding spatial prediction. First, large-scale simulation slices adapted to the original measurements are determined, i.e., three columns with known and unknown formations. Subsequently, the simulation slices are transferred to the training section to scan the training profiles. Convolution in the CNN is used to extract the general formation connectivity in the training profiles. The extracted formation profiles can be used as input data to train an XGBoost model, which is used for forward spatial prediction of large-scale grids. The spatial prediction of large-scale patterns follows a random simulation path (dashed line in Figure 3A ), until all unknown points at the current grid scale have been predicted. In this example, the above process can only predict the soil type at a single grid scale.

[0141] Figure 3BShows an iterative prediction process using gradually decreasing simulation slices. The results obtained from the previous grid scale are treated as additional measurements and combined with the original measurements for the next round of spatial prediction. The above process is repeated until predictions have been made for all unknown points in the simulation profile. It should be noted that the 2D IC-XGBoost model cannot be used for extrapolation. Additionally, if directly used for 3D geological modeling, the iterative convolution and prediction processes in the 2D IC-XGBoost model will be time-consuming.

[0142] To address these issues, some embodiments of the present invention provide an enhanced IC-XGBoost algorithm, also known as IC-XGBoost3D, which will be described in more detail below. Briefly, however, in these embodiments of the enhanced IC-XGBoost algorithm, new multi-scale templates are designed to enable extrapolation, and a pre-training strategy is applied to improve computational efficiency, which is particularly useful for shaping 3D geological models. By using the enhanced IC-XGBoost algorithm in these embodiments, unsampled points within the simulation slice can be determined more accurately and efficiently. Additionally, in some embodiments of the enhanced IC-XGBoost algorithm, for a given simulation slice, extrapolation is first performed, followed by interpolation. Multi-scale templates suitable for the simulation slice are determined. These templates can be used to scan parallel training images to extract relevant stratigraphic profiles, which are then used to train the XGBoost model and subsequent conditional simulations.

[0143] Another major part of the framework involves the quantification of stratigraphic uncertainty.

[0144] A complete 3D simulation domain is called a realization. Multiple random realizations are generated conditioned on the same set of site-specific measurements and training images. Each realization represents a possible representation of the 3D subsurface geological model of the site under investigation. The most likely prediction of Z mp (x) (where x is a spatial coordinate and belongs to 3D space, i.e., ) can be obtained by assigning to each point x in 3D space the soil class with the highest occurrence frequency. This most likely prediction is used as the final prediction result of the method in this embodiment. Additionally, the stratigraphic uncertainty associated with the most likely prediction can be explicitly quantified through statistical analysis of multiple realizations. In this embodiment, the uncertainty of the soil type at x is measured by the dispersion Dp(x), which is defined as the proportion of soil types that do not match (compared to Z r )(x)) among N mp realizations. Mathematically, the dispersion Dp(x) is calculated by the following equation.

[0145]

[0146] where I is an indicator function (with a value of 1 if the condition within the parentheses is true and 0 otherwise); Z r (z) is the r-th realization result. The dispersion represents the confidence in the spatial prediction for each point x in 3D space. The dispersion ranges between 0 and 1, and the larger the value, the lower the confidence in the prediction. A 3D dispersion map can be constructed by combining the dispersions estimated at each point x in 3D space.

[0147] The principle of random simulation in some method embodiments of the present invention is the same as those of geostatistical simulations (such as SIS, MPS, and SGS). The uncertainty in the methods in some embodiments of the present invention may come from two sources: random simulation and neural networks. In this embodiment, the last activation function of the neural network employed is the softmax function, which normalizes the output of the network to a probability distribution of the predicted output classes. The uncertainty from these two sources is caused by limited local data. Equation (2) only calculates the dispersion of the model output, rather than the "uncertainty" being higher or lower in terms of the probability that the model corresponds to the "true" geology.

[0148] In some embodiments, when the most likely prediction derived does not change significantly as the number of realizations increases, the random simulation process terminates. In this embodiment, the percentage change in the most likely prediction Pc for every k additional realizations is used to regulate the simulation process.

[0149]

[0150] where N X 、N Y and N Z are the total number of points in the X, Y, and Z directions respectively; and represent the soil types at x r -k and N r derived from the N i realizations respectively. For example, in the illustrative example of the present disclosure, k is taken as 10 and a Pc value of 0.5% is specified as the termination criterion. For the following illustrative example, the underlying ground truth geology domain is available. Therefore, the generated 3D model can be compared with the ground truth model to verify the method embodiments of the present invention. The deviation of multiple realizations from the ground truth geology model can be measured by accuracy (Acc):

[0151]

[0152] where, when the most likely prediction Z mp (x i ) and the underlying ground truth geology domain Z T (xi )At the same position x i When the soil types are the same, the value of the indicator function I is 1. The above formula can be used to calculate the accuracy of each soil type by simply replacing the 3D dimensions (N X ×N Y ×N Z ) with the dimensions of the corresponding soil types. It should be understood that in actual on-site surveys, the basic reference truth 3D geological domain is not available, and here the method of this embodiment is only verified by accuracy.

[0153] IC-XGBoost3D algorithm

[0154] As described above, some embodiments of the present invention provide an IC-XGBoost3D algorithm for 3D geological modeling. The IC-XGBoost3D algorithm can be regarded as an improved version of the above 2D IC-XGBoost algorithm. Specifically, two improvements are made to the 2D IC-XGBoost algorithm to improve its capabilities and computational efficiency. One of the improvements is to handle extrapolation (when there are no site-specific measurements available on the model boundary). In some embodiments, as long as there are no site-specific measurements along the model boundary, extrapolation of the model boundary is first performed. After extrapolation, the 2D IC-XGBoost algorithm is used to delineate the unsampled positions within the 2D geological profile.

[0155] Figure 4 Shows the process of shaping a 2D geological profile when there are no site-specific measurements along the model boundary. Specifically, Figure 4 Shows the delineation of a single simulation slice using the IC-XGBoost3D algorithm in this embodiment.

[0156] See Figure 4 , in 402, a multi-scale template adapted to the simulation slice is determined, and in 404, all potential multi-scale templates are collected. In 406, the template is applied to exhaust all multi-scale stratigraphic profiles from the training image. In 408, all the collected stratigraphic profiles are used to train or pre-train the XGBoost model, and in 410, spatial prediction of the 2D simulation slice is performed.

[0157] Referring to 402, the 2D simulation slice SS Xz,Y＝a contains three boreholes BH1, BH2, and BH3. In 406, a multi-scale template suitable for SS Xz,Y＝a is used to scan the parallel training image TI XZ, and extract the relevant formation connectivity. As shown in 402 and 404, in this example, each template includes nine cells and is divided into three columns with different horizontal spacings. For example, templates B and C are adapted to simulate slices, and there is an unknown central column on the center line between two columns of known soil formations. Then, the templates determined in 404 are transferred to parallel training images to extract formation profiles, and the extracted formation profiles are used to train the XGBoost model and subsequent conditional simulation in the 2D IC-XGBoost algorithm. In this embodiment, template A is designed for extrapolation. The unknown column of this template is located at the right boundary of SS XZ,Y＝a and the other two spaced columns are aligned with the two nearest boreholes BH1 and BH2. The new template A is specifically designed for the extrapolation operation and is based on the assumption that the soil type along the unknown profile is related to the formations revealed by the nearest neighbors. In this example, among all potential multi-scale templates (A, B, C, and D) from the simulated slices, compared with other templates, template A is preferred for feature extraction from the training image and subsequent conditional simulation. This is because templates B, C, and D are designed for interpolation, and their performance depends on the calculation results of extrapolation using template A. For example, the right column of template D is aligned with the right model boundary, and the interpolation of the central unknown column depends on the result of extrapolation of the right column using template A. It should be noted that the three columns of the template do not necessarily need to have equal spacing. In some embodiments, the spacing can be adapted to the simulation results updated at each simulation step. For example, as new simulation results are added, the spacing between adjacent boreholes gradually decreases. A multi-scale template with a reduced spacing between adjacent columns is used to extract a formation profile with a reduced scale.

[0158] Another improvement is to determine the multi-scale simulation templates for pre-training, or to determine all multi-scale simulation templates of different scales suitable for the simulated slices at the beginning of the simulation. For example, in Figure 4 , the direct templates adapted to the original site-specific measurements include templates A, B, and C. Template D is classified as an intermediary because its performance depends on the results calculated in previous simulation steps. All direct and intermediary templates are predetermined at the beginning of the simulation. The collected multi-scale templates are used to exhaust all formation relationships embedded in the training image. All extracted spatial statistics are used to pre-train the XGBoost model, which is continuously used for the sequential delineation of 2D simulated slices. As shown in 404 and 406, the three templates A, B, and C and other templates with reduced spacing are collected together and applied to a single training profile TI XZ, to extract all potential stratigraphic relationships. All the extracted spatial statistics are used as inputs for training an XGBoost model, which is then used for extrapolation and interpolation. In this embodiment, the pre-training strategy can significantly reduce the computational cost associated with the method of this embodiment. This is because at the beginning of 3D geological modeling, the features extracted from the training images are only performed twice: once for extrapolation and once for interpolation. In one example, due to the implementation of the pre-training strategy, the simulation time for delineating a single 2D simulation slice can be reduced from more than 5 minutes to less than 15 seconds (using a computer equipped with Core TM i7-4790 CPU@3.60GHz and 8.00GB RAM), greatly improving the computational efficiency. In this example, a Python script is developed to run the numerical simulation.

[0159] Illustrative example

[0160] Two 3D geological domains with dimensions of 100m (X-distance) × 100m (Y-distance) × 50m (Z-depth) are simulated to demonstrate the operability and performance of the method in some embodiments of the present invention. In this example, each complete 3D geological domain is discretized into a point cloud of 100×100×50 with a resolution of 1m, for a total of 500,000 voxels. The resolution and model size adopted are generally consistent with typical geotechnical site investigations for excavation and slopes. In this example, the stratigraphic boundaries separating different soil layers are modeled as 3D surfaces and can be described by the following equation:

[0161] AX 2 +BX+CY 2 +DY+EZ 2 +FZ+D = 0 (5)

[0162] where X, Y, Z represent the directions of the standard Cartesian coordinate system, and A, B, C, D, E, F, G are the coefficients of this 3D surface. Voxels located above the 3D surface and voxels located below the 3D surface are assigned different soil types. In this example, there are a total of six 3D surfaces. Table 1 summarizes the coefficient statistics of all 3D stratigraphic boundaries.

[0163] Table 1 – Summary of the coefficients of the 3D stratigraphic function

[0164] AX 2 +BX+CY 2 +DY+EZ 2 +FZ+G = 0

[0165]

[0166]

[0167] Figure 5A and Figure 5B shows the generated 3D training and test geological domains with four soil types. Figure 5C Shows a perspective view of the 3D model. It can be seen that Soil 4 mainly consists of two parts. The shallower part is mainly located at a depth of about 20 m while the deeper part is distributed between depths of 30 m and 50 m. As an example, as Figure 5D shown, two vertical slices are taken from the 3D training geological domain to be used as training images. The two vertical slices are perpendicular to each other and are located on the centerline of the 3D training geological domain. In this example, the planar area of the site is 100 m × 100 m. If the borehole spacing is 20 m, the total number of boreholes is 25 (= 5 × 5). Figure 5E Shows 25 site-specific boreholes extracted from the 3D test geological domain. These 25 boreholes account for 0.25% of all voxels in the 3D point cloud. Although the boreholes are evenly distributed in this embodiment, it should be noted that some embodiments of the present invention are also applicable to cases of site-specific measurements with irregular intervals.

[0168] Further simulation results are obtained using the method of this embodiment.

[0169] After the above simulation process, multiple 3D realizations are generated. The total simulation time for generating the 3D subsurface geological model is approximately 14 minutes (using a personal computer equipped with an Intel(R) Core(TM) i7-10750H CPU @ 2.60 GHz and 16.0 GB RAM), which is generally acceptable for practical engineering applications. Using a more advanced computer or efficient programming techniques can reduce the simulation time.

[0170] In this example, the stochastic simulation process terminates when the percentage change Pc in the most likely prediction drops to 0.5%.

[0171] Figure 6 Shows the percentage change in the most likely prediction for every 10 realizations. It can be seen that the percentage change shows a monotonically decreasing trend, and when the number of realizations increases to more than 20, the percentage becomes less than 0.5%. In this example, 100 3D realizations are generated conditioned on a pair of vertical training images and 25 site-specific boreholes.

[0172] Figure 7A and Figure 7B Shows the results of the method of this embodiment conditioned on 25 boreholes, where Figure 7A shows the most likely 3D prediction while Figure 7B shows the dispersion plot.

[0173] Specifically, Figure 7AShows the most likely 3D geological model derived from 100 realizations. The overall accuracy calculated using Equation (4) is 95.4%. The relevant formation connectivity revealed from the three outer surfaces at X = 100, Y = 100, and Z = 0 correctly replicates Figure 5B the true spatial pattern in

[0174] . The dispersion of 100 3D realizations relative to the most likely prediction was calculated using Equation (2). Figure 7 shows the relevant 3D dispersion map. It should be noted that the large dispersion value bands are mainly concentrated around the predicted soil layer boundaries.

[0175] Figure 8A and Figure 8B show the comparison between the most likely prediction result and the true geological domain, where Figure 8A shows the variation of the prediction accuracy of the most likely prediction result with the number of realizations, while Figure 8B shows the accuracy map.

[0176] Specifically, Figure 8A shows the variation of the prediction accuracy of different soil types with the number of realizations. As Figure 8A shown, when the number of realizations is less than 40, the prediction accuracy of each soil type gradually increases. Although all soil types have different thicknesses and spatial connectivities, the prediction accuracy of Soil 2 remains the lowest among all four soil types regardless of the number of realizations. This is because the training formation profile of Soil 2 is relatively limited in the training image compared to other soil types. When the number of realizations increases to more than 40, all prediction accuracies approximately stabilize.

[0177] Figure 8B shows the accuracy of the 3D model calculated using Equation (4). It can be seen that the areas with low accuracy or large prediction errors are mainly concentrated near the true formation boundaries, which is similar to the spatial pattern revealed by the Figure 7B dispersion map in

[0178] . Based on the good correlation between the dispersion map and the accuracy map, the dispersion map can be used alone to quantify the prediction uncertainty.

[0179] Figures 9A to 9I shows the prediction results for the selected slices at X = 0, 10, and 40. It should be noted that the slice (X = 0, 40) does not contain any initial site-specific boreholes. In contrast, the slice at X = 10 contains five boreholes, providing more information than the other two slices. Figure 9AShows the most likely prediction derived from 100 realizations at X = 0. Although the results were obtained by both extrapolation and interpolation, an accuracy of approximately 94.3% can be achieved based on Equation (4). It should be noted that the true soil layer boundaries extracted at the same location in Figure 5B have also been superimposed as dashed lines for comparison. Clearly, the most likely prediction reasonably captures the overall stratigraphic relationship. Figure 9B The dispersion plot in shows that the largest stratigraphic uncertainties are mainly concentrated around the true soil layer boundaries. For this example, the base ground truth geological profile is available. Figure 9C The accuracy plot in shows that most of the prediction errors fluctuate around the true soil layer boundaries, which is consistent with the Figure 9B dispersion results in. The consistent results between the dispersion and accuracy mean that the dispersion can be used alone to quantify the associated prediction uncertainties. For the slice at X = 10, its most likely prediction ( Figure 9D ) achieved an accuracy of 97.3%, which is the highest among the three slices because it contains more site-specific boreholes than the other two slices. The large dispersion area ( Figure 9E ) and low accuracy (or large prediction errors, Figure 9F ) are mainly concentrated around the true soil layer boundaries. For the most likely prediction at X = 40 ( Figure 9G ), an accuracy of 90.4% can be achieved. It should be noted that the central horizontal band of Soil 4 does not fully replicate the true stratigraphic connectivity. This prediction error is reflected as high dispersion values in the Figure 9H dispersion plot. Essentially, the dispersion plot represents the confidence in the prediction using the method of this embodiment. Figure 9I The accuracy plot in is similar to the stratigraphic profile in the dispersion plot.

[0180] Figures 10A to 10I Shows the prediction results for the selected slices at Y = 60, 90, and 100. It should be noted that all three most likely predictions have an accuracy of over 92%. The true soil layer boundaries have also been superimposed as dashed lines. Clearly, the method of this embodiment correctly replicates the stratigraphic relationships of the three slices (see Figure 10A , Figure 10D and Figure 10G ). The large dispersion areas are mainly concentrated around the true soil layer boundaries (see Figure 10B , Figure 10E and Figure 10H ). Figure 10C , Figure 10F and Figure 10I The accuracy plots in share the same spatial pattern as the corresponding dispersion plots, which further justifies the application of the dispersion plots to uncertainty quantification.

[0181] Perform a sensitivity analysis further.

[0182] The 3D geological modeling method in this embodiment relies on the sequential shaping of 2D geological profiles using IC-XGBoost. The inventors of the present invention have designed that, in some embodiments, for a 2D geological profile discretized into a total of 100 to 150 grids in the horizontal direction, four boreholes (borehole data) can be considered sufficient to obtain reasonable modeling results. Therefore, in some embodiments, if any simulation slice contains less than 4 boreholes, the performance of the 3D modeling method may be reduced.

[0183] The effects of irregular borehole spacing and the number of training images on the 3D subsurface geological model obtained using the method of this embodiment have also been studied.

[0184] Figure 11A and Figure 11B shows the results of a sensitivity study using different combinations of training images and site-specific boreholes, where Figure 11A shows the case of irregular borehole spacing while Figure 11B shows the case of using a single training image in this method. The two vertical training images in the first case (see Figure 11A ) are in the same positions as those in Figure 5D , but the spacing of the site-specific boreholes varies between 5m and 40m, and this variation is arranged in a geometric sequence along the X and Y directions, i.e., 5, 10, 20, 40, 80. In contrast, Figure 11B shown in the second case has the same number and positions of site-specific boreholes as those in Figure 5E , but only one training image extracted from the training geological domain (see Figure 5A ) is used to reflect the prior geological knowledge along the X and Y directions.

[0185] To study the effect of irregular borehole spacing on the method of the embodiment, 100 realizations were generated until the percentage change in the most likely prediction decreased to less than 0.5%.

[0186] Figures 12A to 12F shows the effect of irregular borehole spacing on the prediction results, especially the comparison of the prediction results conditioned on 25 boreholes with different spacings. For ease of comparison, three slices (X = 20, Y = 20, and Z = 5) were selected. Figure 12A shows the most likely prediction when the 25 boreholes are equally spaced. The most likely prediction results along the two vertical slices (X = 20, Y = 20) reasonably capture the spatial connectivity. Figure 12B The dispersion diagram in Figure 12C shows that the formation changes are mainly concentrated around the predicted soil layer boundaries. Similarly, Figure 12B The accuracy diagram in Figure 12DShows the most likely prediction for the condition of 25 boreholes with different spacings. For this figure, the prediction accuracy calculated using formula (4) is 91.6%, which is less than Figure 12A 95.4% in Figure 11A . It should be noted that since the borehole positions are arranged in a geometric order with different densities, and only one borehole is located in the area where X>40m and Y>40m (see Figure 12D ), the accuracy of the predicted 3D geological model is reduced, especially in the areas with larger borehole spacings (i.e., lower borehole densities). Due to the high borehole density in the areas where X<40m and Y<40m, the formation connectivity within the two vertical slices (X = 20, Y = 20) is better constrained, and a finer formation profile (e.g., for soil 4) can be reasonably predicted. The largest difference occurs in the horizontal plane at Z = 5, especially near the lateral boundary intersecting with X = 100. This can be explained by the fact that the boundary line is derived using extrapolation based on the boreholes at a distance (e.g., with intervals of 20m and 60m), which results in a high formation uncertainty in the prediction of soil types. In addition, the large formation uncertainty propagates to the areas near the outer boundary because the internal areas near the boundary line are interpolated based on the results extrapolated at the boundary line. In Figure 12D , in the low borehole density area along Z = 5, due to the high formation uncertainty, the predicted soil layer boundaries (e.g., soil 4) are not as smooth as those in Figure 12A . A similar phenomenon is also observed in terms of the dispersion map (see Figure 12E ). For the two vertical slices (X = 20, Y = 20), since the random simulation is better constrained by the high-density boreholes, the bands of large dispersion values are less concentrated, especially in the areas where X<40m and Y<40m. On the contrary, the dispersion distribution along the horizontal plane (Z = 5) varies greatly. A similar trend is also noticed in the accuracy map in Figure 12F . In Figure 12F , the areas with low accuracy or large prediction errors are much wider, especially in the areas near the lateral model boundary. Therefore, in some embodiments, the presence of some boreholes during the in-situ investigation near the site boundary (and thus the availability of the relevant borehole data) is beneficial for mitigating the potential propagation of formation errors during the shaping of the 3D subsurface geological model.

[0187] Figures 13A to 13FShows the impact of using a single training image on the prediction results in this method. In this example, a complete 3D geological domain is constructed conditioned on site-specific measurements and prior geological knowledge reflected in a single training image. It should be noted that the ground conditions in this example can be classified as "isotropic" because all six formation connectivity patterns (i.e., soil 1 / soil 2, soil 1 / soil 4, soil 2 / soil 3, soil 2 / soil 4, soil 3 / soil 4, and soil 4 / soil 3) shown in the existing site-specific measurements can find corresponding ones in Figure 11B the single training image. Thus, the single training image is considered to contain typical formation relationships along the X and Y directions. Using the 3D modeling method of this embodiment, 100 3D realizations are generated. The relevant most probable predictions, dispersions, and accuracies are obtained using Equation (2) and Equation (3) respectively. Figures 13A to 13C Shows the prediction results of three selected slices (X = 0, Y = 0, and Z = 0) conditioned on a single training image and 25 site-specific boreholes.

[0188] Figure 13A The most probable prediction in Figure 5B reasonably replicates the Figure 13B formation connectivity in Figure 13C . The dispersion in Figures 13D to 13F indicates that formation uncertainties are mainly concentrated around the predicted formation boundaries. Regions with large prediction uncertainties overlap with regions of low accuracy or large prediction errors in Figure 13D . When only a single training image is combined with site-specific measurements for conditional simulation, the calculated results are shown in Figure 13A . The most probable prediction in Figure 13E and Figure 13F reaches an accuracy of 93.8%, slightly less than the 95.4% in Figure 13B and Figure 13C . In contrast, the dispersion and accuracy maps in

[0189] show that regions with large dispersions and prediction errors are more scattered than those in

[0190] Figures 14A to 14I Shows the prediction results along the selected slices (Y = 60, 90, 100) using a single training image (X = 50). Figure 14AShows the most likely prediction along Y = 60. For better comparison, the true formation boundaries are also superimposed as dashed lines. Although an accuracy of 93.5% is obtained, during the simulation, the soil 4 formation boundary at a depth of 10 m to 20 m is not well constrained. This can be explained by the fact that all formation profiles are learned from a single training image rather than parallel training images. Figure 14B The dispersion in Figure 14B indicates that the formation variations vary greatly around the true formation boundaries. Figure 14C The regions of low accuracy or large prediction errors in Figure 14C are similar in pattern to the regions in Figure 14B Figure 14B . Figure 14D and Figure 14G The accuracies of the most likely prediction results along Y = 90 and 100 in Figure 14G are approximately 95.7% and 89.7% respectively. Both of these accuracies are slightly less than the accuracies (96.8% and 92.4%) calculated using a pair of vertical training images in Figure 10D and Figure 10G Figure 10G . Figure 14E and Figure 14H The dispersion in Figure 14E and Figure 14H is more spread out around the true soil boundaries. There is a good correlation between the accuracy maps ( Figure 14F and Figure 14I Figure 14I ) and the dispersion maps ( Figure 14E and Figure 14H Figure 14H ). It is worth mentioning that even when using a single training image for conditional simulation, the 3D spatial distribution of the formation profiles can be reasonably replicated. This demonstrates the feasibility of using a single training image for 3D geological modeling of isotropic ground conditions in some embodiments.

[0191] Among other things, the above disclosure has provided a 3D subsurface formation modeling and uncertainty quantification method (IC-XGBoost3D) that can effectively generate a 3D subsurface geological model from sparse site-specific measurements and prior geological knowledge. In some embodiments, the method relies on sequential simulation of 2D slices, which are subsequently recombined into a complete 3D geological domain. In some implementations, an enhanced machine learning algorithm called iterative convolutional XGBoost (IC-XGBoost, e.g., 2D IC-XGBoost) can be used to implement the simulation of each 2D slice. In some implementations, the potential formation connectivity is learned from a pair of vertical training images, which are a collection of recognized (credible) spatial statistics in the region of interest. In some embodiments, the 3D subsurface geological modeling method is purely data-driven. The above results show that the method in the above embodiments of the present invention can accurately reconstruct a 3D geological domain conditioned on a pair of training images (e.g., vertical) and limited site-specific boreholes. In some embodiments, the method not only infers the most likely 3D subsurface geological model but also provides a quantitative assessment of the associated 3D formation uncertainty.

[0192] Figure 15An example information processing system 1500 of an embodiment of the present invention is shown, which can be used as a server or other types of information processing systems. The information processing system 1500 can be used to execute one or more method embodiments (partial or all) of the present invention. The information processing system 1500 generally includes appropriate components required to receive, store, and execute appropriate computer instructions, commands, and / or codes. The main components of the information processing system 1500 are a processor 1502 and a memory (memory) 1504. The processor 1502 can include one or more of: one or more CPUs, one or more MCUs, one or more logic circuits, one or more Raspberry Pi chips, one or more digital signal processors (DSPs), one or more application specific integrated circuits (ASICs), one or more field programmable gate arrays (FPGAs), and / or any one or more digital or analog circuits / circuits configured to interpret and / or execute program instructions and / or process signals and / or information and / or data. The memory 1504 can include one or more volatile memories (such as RAM, DRAM, SRAM), one or more non-volatile memories (such as ROM, PROM, EPROM, EEPROM, FRAM, MRAM, FLASH, SSD, NAND, NVDIMM), or any combination thereof. Appropriate computer instructions, commands, codes, information, and / or data can be stored in the memory 1504. The computer instructions for executing or facilitating the execution of the method embodiments of the present invention can be stored in the memory 1504. The processor 1502 and the memory (memory) 1504 can be integrated or separated (and operably connected). Optionally, the information processing system 1500 further includes one or more input devices 1506. Examples of the input device 1506 include: keyboard, mouse, stylus, image scanner, microphone, tactile / touch input device (such as a touch-sensitive screen), image / video input device (such as a camera), etc. Optionally, the information processing system 1500 further includes one or more output devices 1508. Examples of the output device 1508 include: display (such as a monitor, screen, projector, etc.), speaker, headphones, earphones, printer, additive manufacturing machine (such as a 3D printer), etc. The display can include an LCD display, an LED / OLED display, or other suitable displays, which can be or can not be touch-sensitive. The information processing system 1500 can also include one or more disk drives 1512, which can include one or more of the following: solid state drive, hard disk drive, optical drive, flash drive, tape drive, etc. A suitable operating system can be installed in the information processing system 1500, for example, in the disk drive 1512 or the memory 1504. The memory 1504 and the disk drive 1512 can be operated by the processor 1502.Optionally, the information processing system 1500 further includes a communication device 1510 for communicating with devices such as servers, personal computers, terminals, tablets, telephones, watches, Internet of Things (IoT) devices, or other wireless computing devices. The communication device 1510 may include one or more of the following: a modem, a network interface card (NIC), an integrated network interface, an NFC transceiver, a ZigBee transceiver, a Wi-Fi transceiver, a Bluetooth transceiver, a radio frequency transceiver, a cellular (2G, 3G, 4G, 5G, 5G+, etc.) transceiver, an optical port, an infrared port, a USB connection, or other wired or wireless communication interfaces. The transceiver may be implemented by one or more devices (integrated transmitter and receiver, separate transmitter and receiver, etc.). The communication link may be wired or wireless and is used to transmit commands, instructions, information, and / or data. In one example, the processor 1502, the memory 1504 (optionally, the input device 1506, the output device 1508, the communication device 1510, and the disk drive 1512, if present) are directly or indirectly interconnected with each other through a bus, a peripheral component interconnect (PCI) such as PCI Express, a universal serial bus (USB), an optical bus, or other similar bus structures. In one embodiment, at least some of these components may be wirelessly connected, for example, through a network such as the Internet or a cloud computing network. Those skilled in the art will understand that... Figure 15 The information processing system 1500 shown in [reference] is only an example, and in other embodiments, the information processing system 1500 may have different configurations (e.g., including additional components, having fewer components, etc.).

[0193] In some embodiments, the method of the present invention may be integrated with or implemented on a computer-implemented geographic information platform such as a Geographic Information System (GIS) platform. The computer-implemented geographic information platform may run as an application (desktop, mobile APP, web APP, etc.).

[0194] Figure 17 A method 1700 for operating a geographic information platform according to an embodiment of the present invention is shown. The method 1700 may be implemented using suitable hardware and software. The geographic information platform may be implemented using suitable hardware and / or software. The method 1700 may be considered to implement the method 1600 on or in conjunction with the geographic information platform. Therefore, steps 1706, 1708, and 1710 in the method 1700 generally correspond to the corresponding steps 1602, 1604, and 1606 in the method 1600.

[0195] In step 1702 of method 1700, it is determined whether new borehole data and / or new prior geological data are available. This determination can be performed by an information processing system operating a geographic information platform. This determination can be performed automatically. In one example, new borehole data and / or new prior geological data are uploaded, transmitted, input, etc. to the information processing system or another system operably coupled to the information processing system, and the determination can be made by actively monitoring such inputs or based on notifications received or retrieved related to the availability of new data. If it is determined that no new borehole data and new prior geological data are available, method 1700 stops or remains at step 1702 to monitor the availability of new data. If it is determined that new data is available, method 1700 proceeds to step 1704 to obtain the new borehole data and / or new prior geological data (those that are available). The new data obtained can be stored in a memory, internal memory, database, server, etc. for subsequent use.

[0196] In step 1705, the information processing system or the geographic information platform receives a request to generate a 3D geological model of a geographic site via a user interface. If the request is received, method 1700 proceeds to step 1706 to generate a geological profile of the subsurface volume of the geographic site using the borehole data and prior geological data (including the new borehole data and / or new prior geological data obtained in step 1704, if any). The details of step 1706 are substantially the same as those of step 1602 in method 1600 and will not be repeated here. Similarly, the details of steps 1708 and 1710 are substantially the same as those of steps 1604 and 1606 in method 1600 and will not be repeated here. The 3D geological model generated in step 1710 is presented to the user via the user interface of the geographic information platform. Figure 18 An example user interface 1800 of the geographic information platform is shown.

[0197] In one embodiment, if new borehole data and / or new prior geological data are available and are obtained while the 3D geological model is being formed or displayed, the 3D geological model is updated based on the new data. In one embodiment, the 3D geological model can thus be dynamically and automatically updated.

[0198] The above-described embodiments of the present invention have one or more unique features. For example, the methods of some embodiments bypass the need for 3D training images (which are required in some existing 3D geological modeling methods), thus addressing the challenge of the absence of 3D training images in geological or geotechnical applications. For example, the methods of some embodiments construct a 3D geological domain by sequentially generating 2D geological profiles conditioned on 2D training images and site-specific borehole data and recombining them into a 3D point cloud. The simulation sequence is established based on the principle that each current 2D vertical slice contains as much conditional data as possible. The conditional data includes the original site-specific boreholes and the results obtained from the previous simulation steps. For example, in some embodiments, credible formation statistics (regarded as "truth") can be briefly and implicitly stored in a pair of vertical 2D training images that reveal potential horizontal anisotropy and can be obtained from previous projects with similar geological settings. A 3D representation of the subsurface geology is generated through sequential prediction at unsampled locations.

[0199] For example, in the methods of some embodiments, the sequential simulation of each 2D vertical slice is performed using an enhanced IC-XGBoost algorithm ("IC-XGBoost3D") with improved capabilities and computational efficiency. The 2D IC-XGBoost algorithm previously developed by the inventors relied on a simulation template (i.e., a 3×3 cell grid with the unknown cell in the central column) that was used to extract formation features from 2D training images and perform spatial interpolation. Due to such a template, the 2D IC-XGBoost algorithm could not be used for extrapolation. In IC-XGBoost3D, one enhancement is the incorporation of extrapolation (which is often encountered in engineering practice when there are no site-specific boreholes available along the model boundary). A new simulation template is designed for pattern extraction and extrapolation, with the unknown cell in the side column. As long as there are no site-specific boreholes along the model boundary, the extrapolation prediction of the model boundary is first performed. After extrapolation, the 2D IC-XGBoost algorithm is used to demarcate the unsampled locations within the 2D geological profile or spatial interpolation. Another enhancement is the determination of all multi-scale simulation templates of different sizes suitable for the simulation slices at the beginning of the simulation. The collected multi-scale templates are applied to exhaust all formation profiles embedded in the 2D training images. All extracted spatial statistics are used to pre-train an XGBoost model, which is continuously used for the sequential demarcation of 2D vertical slices. In some embodiments, the implementation of the pre-training strategy can significantly improve the simulation efficiency of the demarcation of individual 2D vertical slices. This may be important for shaping a 3D subsurface geological model using IC-XGBoost3D, as a large number of 2D simulations need to be performed to combine the 3D model.

[0200] Some embodiments of the present invention may provide an automated method for 3D subsurface formation modeling and uncertainty quantification, which facilitates the establishment and update of a 3D geological model of a subsurface digital twin conditioned on a limited number of site-specific boreholes and 2D training images reflecting prior geological knowledge. Modeling of the 3D geological domain may include sequential simulation of many 2D vertical slices following a random simulation path and recombining them into a complete 3D geological domain. Formation uncertainty may be quantified by statistical analysis of multiple equally-probable 3D geological realizations generated under a Monte Carlo simulation framework. The 2D training images may utilize existing 2D geological profiles from nearby sites or previous projects with similar geological environments. The random 2D simulation sequence may be determined based on the principle that each current vertical slice should contain the most site-specific data, and a random draw is made if there are two or more simulation slices containing the same amount of site-specific data. Prediction of the 2D simulation slices may be performed using an enhanced 2D IC-XGBoost algorithm called IC-XGBoost3D. The conditional data for simulating the 2D vertical slices may include existing site-specific boreholes and simulation results from previous simulations. The IC-XGBoost3D algorithm may achieve extrapolation prediction and may implement a pre-training strategy to achieve the computational efficiency required for 3D modeling. If there are unsampled points along the boundary of the 3D geological model, extrapolation prediction may always be performed before interpolation prediction. The pre-training strategy may be implemented by determining all multi-scale templates (i.e., 3×3 cell grids) based on site-specific borehole data and then applying the templates to exhaustively embed all formation relationships in the training image to train a classification algorithm (i.e., XGBoost or multi-layer perceptron).

[0201] Example functions of some embodiments of the present invention include reconstructing or constructing a 3D subsurface geological model from sparse site-specific boreholes and 2D training images reflecting prior geological knowledge. Some embodiments of the present invention can be used to build and dynamically update a 3D subsurface geological domain in a digital twin of subsurface critical infrastructure for smart city applications. In the past decade, different cities around the world have developed specific strategies to promote the development of smart cities. One of the tasks in developing a smart city is to establish a 3D digital twin, which is a virtual representation of a 3D physical entity in the real world. Digital twins can help practitioners streamline the design, construction, performance monitoring, and maintenance of subsurface infrastructure. Therefore, there is a need for an effective method to accurately determine 3D subsurface formations and quantitatively evaluate the associated formation uncertainties involved in the subsurface digital twin. The methods and 3D geological models in some embodiments of the present invention can be used in combination with Internet of Things (IoT) technologies for real-time monitoring of civil construction processes for data-driven engineering analysis and decision-making. The methods of the present invention can be packaged and coded, for example, in a Python package, and / or incorporated into or otherwise applied to a Geographic Information System (GIS) platform, such as ArcGIS.

[0202] Some embodiments of the present invention may include one or more of the following advantages. First, some embodiments utilize valuable prior geological knowledge that is typically overlooked in traditional statistical methods and represent it using 2D training images. Second, some embodiments rely on an image-based stochastic simulation framework to infer the best estimate of the 3D subsurface geological domain conditioned on 2D training images and site-specific data and quantify the associated formation uncertainties in a data-driven manner. For traditional statistical methods, such as Markov-based models (e.g., hidden Markov chains and Markov random fields), the parametric functional form that controls the parameter distribution must be specified in advance. For these traditional statistical methods, the specification and determination of the relevant parameters are challenging when only limited site-specific measurements are available. Current engineering practice largely relies on engineering judgment to manually / humanly digitize the 3D subsurface geological domain for engineering analysis and design. Some embodiments can help contractors and designers build a digital twin of the subsurface critical civil infrastructure, which can be used for intelligent geotechnical site investigations, geotechnical engineering design, construction site analysis, etc.

[0203] Those skilled in the art will understand that the present invention as shown in specific embodiments can be subject to different changes and / or modifications to provide other embodiments of the present invention. Therefore, the described embodiments of the present invention should be considered illustrative rather than restrictive in all respects. Examples of optional features of some aspects of the present invention are set forth in the Summary of the Invention section above. Some embodiments of the present invention may include one or more of these optional features (some of which are not specifically shown in the drawings). Some embodiments of the present invention may lack one or more of these optional features (some of which are not specifically shown in the drawings). One or more features in one embodiment and one or more features in another embodiment can be combined to provide further embodiments of the present invention.

Claims

1. A computer-implemented method for generating a three-dimensional geological model of a subsurface volume of a geographical site, and comprising: processing, using a machine learning-based model: (a) borehole data of a plurality of boreholes of the geographical site, the borehole data including data related to subsurface soil and / or rock information along at least a portion of each of the plurality of boreholes, and (b) prior geological data related to the geographical site, wherein the prior geological data includes data related to prior stratigraphic information associated with the geographical site; wherein the data related to the prior stratigraphic information associated with the geographical site can be characterized as one or more two-dimensional images of a stratigraphic section associated with the geographical site; to generate a plurality of geological sections of the subsurface volume of the geographical site; each of the plurality of geological sections includes corresponding stratigraphic information; applying the borehole data and the prior geological data to a three-dimensional coordinate system; sequentially processing, using the machine learning-based model, a plurality of slices for the subsurface volume in the three-dimensional coordinate system to continuously generate the plurality of geological sections; and forming the three-dimensional geological model for display based on the generated plurality of geological sections; wherein the processing of each of the plurality of slices includes, for each slice: if there is no borehole data along the boundary of the slice, performing an extrapolation operation using the machine learning-based model to determine the stratigraphic information of the boundary of the slice; and performing an interpolation operation using the machine learning-based model to determine the stratigraphic information of other parts of the slice where there is no borehole data; wherein performing the extrapolation operation includes: determining a template based on the borehole data included in the slice, the template including a plurality of spaced unit columns parallel to the depth direction, the plurality of spaced unit columns including a first column arranged along a set of borehole data, a second column arranged along another set of borehole data, and a third column arranged along the boundary where there is no borehole data for extrapolation; applying the template to a two-dimensional image parallel to the slice to extract stratigraphic relationships and / or statistical data; training the machine learning-based model using the extracted stratigraphic relationships and / or the statistical data; and applying the template to the slice along one or more paths to iteratively process the slice using the trained machine learning-based model to determine the stratigraphic information of the boundary.

2. The computer-implemented method according to claim 1, wherein, The machine learning-based model includes a gradient boosting-based algorithm.

3. The computer-implemented method according to claim 1, wherein, The data related to the prior stratigraphic information associated with the geographical site includes data related to the stratigraphic information of another site near the geographical site and / or another site geologically similar to the geographical site.

4. The computer-implemented method according to claim 1, wherein, The machine learning-based model is trained based on the one or more two-dimensional images.

5. The computer-implemented method according to claim 1, Among them, the one or more two-dimensional images include only a single two-dimensional image; and wherein applying the borehole data and the prior geological data to the three-dimensional coordinate system includes: Arrange the borehole data in the three-dimensional coordinate system along the depth direction; and Arrange the single two-dimensional image to extend along the depth direction in the three-dimensional coordinate system.

6. The computer-implemented method according to claim 1,[[]] Among them, The one or more two-dimensional images include two two-dimensional images; And wherein applying the borehole data and the prior geological data to the three-dimensional coordinate system includes: Arrange the borehole data in the three-dimensional coordinate system along the depth direction; and Arrange the two two-dimensional images to extend along the depth direction in the three-dimensional coordinate system and such that the corresponding imaginary planes containing the two two-dimensional images define a dihedral angle in the three-dimensional coordinate system.

7. The computer-implemented method according to claim 6, wherein, The dihedral angle is approximately 90 degrees.

8. The computer-implemented method according to claim 6, wherein, The plurality of slices includes: A first plurality of slices, which are arranged to be generally parallel to each other and generally parallel to one of the two two-dimensional images; and A second plurality of slices, which are arranged to be generally parallel to each other and generally parallel to the other of the two two-dimensional images; and The first plurality of slices is not parallel to the second plurality of slices.

9. The computer-implemented method according to claim 8, wherein, The processing further includes: Identify the slice containing the most borehole data; First process the slice containing the most borehole data using the machine learning-based model; and Then, (i) If the slice processed first is one of the first plurality of slices, use the machine learning-based model to alternately process one of the second plurality of slices and one of the first plurality of slices until all of the plurality of slices are processed; and (ii) If the slice processed first is one of the second plurality of slices, use the machine learning-based model to alternately process one of the first plurality of slices and one of the second plurality of slices until all of the plurality of slices are processed.

10. The computer-implemented method according to claim 9, wherein, The geological profile generated by the processed slice is used as additional borehole data in the processing of the next slice.

11. The computer-implemented method according to claim 10, wherein, The sequential processing includes: First process the slice containing the most borehole data; and Then process the remaining slices in decreasing order of the amount of borehole data and additional borehole data contained in the slices until all of the plurality of slices are processed.

12. The computer-implemented method according to claim 11, wherein, The sequential processing includes: If two or more of the remaining slices contain the same amount of borehole data and additional borehole data, randomly select one of the two or more slices for processing.

13. The computer-implemented method according to claim 1, further comprising: Display the three-dimensional geological model.

14. A computer-implemented method for generating a three-dimensional geological model of the subsurface volume of a geographical site, and comprising: Using a machine learning-based model, process: (a) The borehole data of a plurality of boreholes of the geographical site, the borehole data including data related to subsurface soil and / or rock information along at least a portion of each of the plurality of boreholes, and (b) Prior geological data related to the geographical site, wherein the prior geological data includes data related to previous stratigraphic information associated with the geographical site; wherein the data related to the previous stratigraphic information associated with the geographical site can be characterized as one or more two-dimensional images of the stratigraphic profile associated with the geographical site. to generate a plurality of geological profiles of the subsurface volume of the geographical site; each of the plurality of geological profiles includes corresponding stratigraphic information; applying the borehole data and the prior geological data to a three-dimensional coordinate system; sequentially processing, in the three-dimensional coordinate system, a plurality of slices for the subsurface volume using the machine learning-based model to continuously generate the plurality of geological profiles; and forming the three-dimensional geological model for display based on the generated plurality of geological profiles; wherein processing each of the plurality of slices includes, for each slice: if there is no borehole data along the boundary of the slice, performing an extrapolation operation using the machine learning-based model to determine the stratigraphic information of the boundary of the slice; and performing an interpolation operation using the machine learning-based model to determine the stratigraphic information of other parts of the slice where there is no borehole data; wherein performing the interpolation operation includes: determining a plurality of templates based on the borehole data included in the slice, each of the plurality of templates including a corresponding plurality of spaced unit columns parallel to the depth direction, the unit columns including a corresponding first column arranged along a corresponding set of borehole data, a corresponding second column arranged along a corresponding other set of borehole data, and a third column in the space between the first column and the second column and without borehole data for interpolation; applying the template to a two-dimensional image parallel to the slice to extract stratigraphic relationships and / or statistical data; training the machine learning-based model using the extracted stratigraphic relationships and / or statistical data; and applying the template to the slice along one or more paths to iteratively process the slice using the trained machine learning-based model to determine the stratigraphic information of other parts of the slice where there is no borehole data.

15. A system for generating a three-dimensional geological model of a subsurface volume of a geographical site, comprising: one or more processors; and a memory storing one or more programs configured to be executed by the one or more processors, the one or more programs including instructions for performing the steps of the method according to claim 1 or 14.

16. A non-transitory computer storage medium storing one or more programs configured to be executed by one or more processors, the one or more programs including instructions for performing the steps of the method according to claim 1 or 14.

Images