Determination method for spatial variability of shearing strength parameter of deep and thick covering layer

Abstract

The invention discloses a method for determining the spatial variability of shear strength parameters of a deep and thick covering layer, particularly relates to the field of geotechnical engineering and engineering geology, and is used for determining the spatial variability of the shear strength parameters of the deep and thick covering layer, and the method comprises the following steps: S1, constructing a three-dimensional geologic model of the covering layer, and initializing a shear strength parameter space field; s2, operating a rock-soil numerical model based on the parameter space field to obtain a simulated displacement field; s3, actually measured displacement data under the corresponding working condition are obtained; s4, comparing the simulated displacement with the actually measured displacement, and if the difference exceeds a threshold value, updating the parameter space field through an optimization algorithm and carrying out iterative calculation until the difference meets a convergence condition; and S5, outputting the converged shear strength parameter space field as a final result. According to the method, an initial shear strength parameter space field is constructed, a rock-soil numerical model is driven to generate a simulated displacement field, the simulated displacement field is compared with actually measured displacement data under the corresponding working condition, and the parameter field is iteratively updated by using an optimization algorithm until the difference between the simulated displacement and the actually measured displacement meets a convergence condition; and finally outputting the shear strength parameter space distribution consistent with the actual deformation behavior.

Description

Technical Field

[0001] This invention relates to the fields of geotechnical engineering and engineering geology, and more specifically, to a method for determining the spatial variability of shear strength parameters in deep overburden layers. Background Technology

[0002] In geotechnical engineering design involving thick overburden layers (such as alluvial and diluvial strata, often exceeding 50 meters in thickness), the spatial distribution of shear strength parameters (including cohesion and internal friction angle) directly affects slope stability, foundation bearing capacity, and the reliability of deformation prediction. Traditional methods primarily rely on limited borehole sampling for indoor or in-situ tests, extrapolating point test results to the entire region through empirical interpolation, or simplifying them to homogeneous or layered homogeneous assumptions for numerical analysis.

[0003] However, deep overburden layers typically exhibit significant heterogeneity and spatial variability. Their complex depositional structures and diverse origins lead to non-stationary variations in shear strength parameters in both the lateral and vertical directions. When using such simplified parameter fields for numerical simulations, the predicted displacement response often deviates systematically from actual monitoring data, thus affecting the accuracy of engineering safety assessments. Therefore, there is an urgent need for a method that can integrate field-measured displacement responses with physical mechanism models to obtain a spatial distribution of shear strength parameters that conforms to actual deformation behavior. To address this issue, a method for determining the spatial variability of shear strength parameters in deep overburden layers is proposed. Summary of the Invention

[0004] To overcome the aforementioned deficiencies of the prior art, this invention constructs an initial shear strength parameter spatial field and drives a geotechnical numerical model to generate a simulated displacement field. This simulated displacement field is then compared with the measured displacement data under the corresponding working conditions. An optimization algorithm is used to iteratively update the parameter field until the difference between the simulated and measured displacements meets the convergence condition. Finally, a spatial distribution of shear strength parameters consistent with the actual deformation behavior is output.

[0005] To achieve the above objectives, the present invention provides the following technical solution: a method for determining the spatial variability of shear strength parameters of deep overburden layers, characterized in that the determination method includes the following steps: S1. Construct a three-dimensional geological model of the overburden and initialize the spatial field of shear strength parameters; S2. Run the geotechnical numerical model based on the parameter spatial field to obtain the simulated displacement field; S3. Obtain the measured displacement data under the corresponding working conditions; S4. Compare the simulated displacement with the measured displacement. If the difference exceeds the threshold, update the parameter space field through the optimization algorithm and iterate the calculation until the difference meets the convergence condition. S5. Output the converged shear strength parameter spatial field as the final result.

[0006] In a preferred embodiment: the construction of the three-dimensional geological model of the overburden in step S1 includes establishing stratigraphic interfaces based on borehole columnar sections, engineering geological profiles, and geophysical exploration data, and dividing the stratigraphic units into multiple geological units with different genetic types or lithological characteristics. The initialization of the spatial field of shear strength parameters includes assigning initial mean values ​​and coefficients of variation of cohesion and internal friction angle to each geological unit, and generating a spatially correlated parameter distribution using a random field method. The spatial correlation distance of the random field is determined according to the depositional environment, grain composition, or consolidation history of each unit, and non-stationary covariance structure modeling is achieved through variogram analysis to reflect the variability differences of the thick overburden in the lateral and vertical directions.

[0007] In a preferred embodiment: the geotechnical numerical model in step S2 is a finite element model or a finite difference model. The boundary conditions of the model are set according to the actual terrain, groundwater level and external load. The material constitutive relation adopts the Mohr-Coulomb model, HardeningSoil model or modified Cambridge model. The stress path evolution caused by construction sequence or water level change is considered in the calculation process. The mesh division of the geotechnical numerical model adopts an adaptive refinement strategy, setting finer element size in potential deformation concentration areas, and using special element types for the foundation interface of the structure, weak interlayer or soil-rock contact surface to accurately simulate slip, separation or shear behavior.

[0008] In a preferred embodiment: the measured displacement data in step S3 includes surface deformation data obtained by synthetic aperture radar interferometry, displacement time histories recorded by global navigation satellite system monitoring points, settlement data obtained by leveling, deep horizontal displacement profiles measured by inclinometers, or compression of each soil layer recorded by layered settlement markers. The time nodes of the measured displacement data correspond to the loading stages of the simulated working conditions in the numerical model. Before use, the original monitoring data is spatiotemporally registered, outlier removed, and multi-source data consistency verified to ensure that the displacements obtained by different monitoring methods are comparable under the same coordinate system and time reference.

[0009] In a preferred embodiment: Step S4, comparing the simulated displacement with the measured displacement, includes interpolating the simulated displacement field to the spatial location of the measured displacement data, calculating the displacement difference between the two at the corresponding monitoring points, wherein the difference exceeds a threshold, meaning that the weighted root mean square of the displacement difference is greater than a preset allowable error, and the optimization algorithm includes a genetic algorithm, a particle swarm optimization algorithm, a Bayesian optimization method, or a gradient-based inversion algorithm, used to adjust the parameter values ​​of each grid node or region in the shear strength parameter spatial field, wherein the weights are dynamically allocated according to the measurement accuracy, data stability, and spatial representativeness of the monitoring points, and parameters of some high-confidence borehole locations are retained as anchor points during the optimization process to prevent overall distortion of the parameter field.

[0010] In a preferred embodiment: during the iterative calculation in step S4, physical rationality constraints are applied to the updated shear strength parameter spatial field, including that the cohesion value is not less than zero, the internal friction angle ranges from zero to forty-five degrees, and the parameter variation between adjacent spatial units does not exceed a preset geological smoothing threshold. The geological smoothing threshold is preset according to the stratigraphic continuity, sedimentary rhythm, or soil structural characteristics. For areas with clear geological interfaces, parameter abrupt changes are allowed, while within homogeneous soil layers, parameter gradient constraints are forcibly implemented to maintain spatial continuity.

[0011] In a preferred embodiment: the spatial field of shear strength parameters output in step S5 is stored in the form of a three-dimensional grid. Each grid cell contains the corresponding cohesion value and internal friction angle value, and can be exported as a data file that can be recognized by general geological modeling format or mainstream geotechnical numerical analysis software. The output results also contain parameter uncertainty information, which is added to each grid cell in the form of standard deviation or confidence interval, and can be directly called by subsequent slope stability analysis modules, foundation bearing capacity verification programs or engineering digital twin platforms.

[0012] The technical effects and advantages of this invention are as follows: This invention constructs an initial shear strength parameter spatial field and drives a geotechnical numerical model to generate a simulated displacement field. The simulated displacement field is then compared with the measured displacement data under the corresponding working conditions. An optimization algorithm is used to iteratively update the parameter field until the difference between the simulated and measured displacements meets the convergence condition. Finally, the spatial distribution of shear strength parameters consistent with the actual deformation behavior is output.

[0013] This process combines the physical constraints of the numerical model with the objective feedback from field monitoring data. Based on the initialization of the random field, it introduces geological unit partitioning and non-stationary covariance structure to give the parameter field reasonable spatial correlation. At the same time, through weight allocation, anchor point retention and physical rationality constraints, it takes into account both data fitting accuracy and geological interpretability during the optimization process. The final output parameter field not only reflects the spatial variation trend, but also adds uncertainty information for subsequent analysis.

[0014] Therefore, the determination of the parameter field no longer depends on the subjective extrapolation of sparse point experiments, but is driven by the inversion of the global displacement response, which systematically improves its spatial continuity, physical consistency and engineering applicability. Attached Figure Description

[0015] Figure 1 This is a flowchart of the steps of the present invention.

[0016] Figure 2 This is a flowchart of the inversion and iteration process for the shear strength parameters of the present invention. Detailed Implementation

[0017] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0018] The exemplary implementation method for determining the spatial variability of shear strength parameters in deep overburden will now be described more fully with reference to the accompanying drawings. The method includes the following steps: Step S1 involves constructing a three-dimensional geological model of the overburden, which includes establishing stratigraphic interfaces based on borehole columnar sections, engineering geological profiles, and geophysical exploration data, and dividing the stratigraphic units into multiple geological units with different genetic types or lithological characteristics. The initialization of the spatial field of shear strength parameters includes assigning initial mean values ​​and coefficients of variation of cohesion and internal friction angle to each geological unit, and generating a parameter distribution with spatial correlation using a random field method. The spatial correlation distance of the random field is determined based on the depositional environment, grain composition, or consolidation history of each unit, and non-stationary covariance structure modeling is achieved through variogram analysis to reflect the variability differences of the deep overburden in the lateral and vertical directions. First, borehole columnar sections, engineering geological profiles, and geophysical exploration results (such as high-density electrical resistivity tomography and seismic refraction / reflection data) are collected within the region. A three-dimensional stratigraphic interface model of the overburden is then established using stratigraphic correlation and interpolation methods (such as Delaunay triangulation or Kriging interpolation). Subsequently, the model is divided into several geological units, each with relatively consistent genetic types (such as alluvial, diluvial, and lacustrine deposits) or lithological characteristics (such as silty clay, medium sand, and gravel). For each geological unit, a cohesion setting is then established. and internal friction angle Initial statistical parameters, including the mean With coefficient of variation Based on this, a spatially continuous parameter field is generated using random field theory.

[0019] Specifically, for any spatial location Its parameter value can be expressed as: ; in, represent or The mean, Standard deviation, For a standardized Gaussian random field, its covariance function Determined by fitting the variation function; Considering the different scales of variation in the horizontal and vertical directions of thick overburden, an anisotropic index model is adopted:

[0020] in These are the spatially related distances in the horizontal and vertical directions, values ​​taken empirically based on the sedimentary environment (e.g., alluvial deposits). For non-stationary situations, different partition settings can be applied. Values ​​are used to model non-stationary covariance.

[0021] In step S2, the geotechnical numerical model is a finite element model or a finite difference model. The boundary conditions of the model are set according to the actual terrain, groundwater level and external load. The material constitutive relation adopts the Mohr-Coulomb model, HardeningSoil model or modified Cambridge model. In the calculation process, the stress path evolution caused by construction sequence or water level change is considered. The mesh division of the geotechnical numerical model adopts an adaptive refinement strategy, setting finer element size in the potential deformation concentration area, and using special element types for the foundation interface of the structure, weak interlayer or soil-rock contact surface to accurately simulate slip, separation or shear behavior. After obtaining the initial shear strength parameter spatial field, it is mapped to the grid element properties of the geotechnical numerical model. This invention can be calculated using finite element software (such as Plaxis, ABAQUS) or finite difference programs (such as FLAC3D). Model boundary conditions are set with displacement constraints based on the actual terrain; groundwater level is applied with pore water pressure according to survey data; external loads include self-weight, building loads, reservoir water level, or seismic inertial forces.

[0022] The material constitutive model can be selected based on the soil type: the Mohr-Coulomb model (suitable for coarse-grained soils), the HardeningSoil model (suitable for normally consolidated clays), or the modified Cambridge model (suitable for overconsolidated soils). During the simulation, loading is performed in stages according to the actual construction or water impoundment sequence to accurately reflect the influence of stress path evolution on deformation.

[0023] To improve computational accuracy, the model mesh employs a local refinement strategy (e.g., reducing element size to 1 / 3 to 1 / 5 of the original size) around potential slip zones, weak interlayers, or the foundation of structures. For the soil-structure interface, Goodman contact elements or interface elements are used, and the friction angle and bond strength are set to simulate possible slippage or delamination behavior. After the calculation is completed, the nodal displacements at the corresponding locations of all monitoring points are extracted to form a simulated displacement field: .

[0024] The measured displacement data in step S3 includes surface deformation data obtained by synthetic aperture radar interferometry, displacement time histories recorded by global navigation satellite system monitoring points, settlement data obtained by leveling, deep horizontal displacement profiles measured by inclinometers, or compression of each soil layer recorded by layered settlement markers. The time nodes of the measured displacement data correspond to the loading stages of the simulated working conditions in the numerical model. Before use, the original monitoring data is spatiotemporally registered, outlier is removed, and multi-source data consistency is verified to ensure that the displacements obtained by different monitoring methods are comparable under the same coordinate system and time reference. The measured displacement data comes from various monitoring methods, including but not limited to: Using Sentinel-1 or TerraSAR-X satellite data, the annual average surface deformation rate was inverted using PS-InSAR technology; three-dimensional displacement time series recorded by GNSS monitoring stations; elevation changes of settlement markers obtained from first- or second-order leveling surveys; and horizontal displacement profiles recorded by inclinometers at different depths. The compression of each soil layer is monitored by the stratified settlement magnetic ring beacon. All measured data need to be preprocessed: first, they are unified to the same coordinate system (such as CGCS2000) and elevation datum; second, the time series is interpolated or averaged to correspond to a specific working condition in the numerical model (such as "30 days after water is impounded to the normal high water level"); finally, obvious outliers (such as data points with abrupt changes exceeding 3 times the mean square error) are removed, and the consistency of multi-source data is checked (such as the difference between InSAR and GNSS at common points does not exceed the monitoring accuracy limit).

[0025] The step S4, comparing the simulated displacement with the measured displacement, includes interpolating the simulated displacement field to the spatial position of the measured displacement data, calculating the displacement difference between the two at the corresponding monitoring points, and the difference exceeding the threshold means that the weighted root mean square of the displacement difference is greater than the preset allowable error. The optimization algorithm includes a genetic algorithm, a particle swarm optimization algorithm, a Bayesian optimization method, or a gradient-based inversion algorithm, which is used to adjust the parameter values ​​of each grid node or region in the shear strength parameter spatial field. The weights are dynamically allocated according to the measurement accuracy, data stability, and spatial representativeness of the monitoring points, and some parameters of high-confidence borehole locations are retained as anchor points during the optimization process to prevent overall distortion of the parameter field. The simulated displacement field is interpolated to the spatial coordinates of the measured monitoring points using trilinear interpolation or radial basis function interpolation to obtain the simulated displacement values ​​at the corresponding locations; the objective function (i.e., the difference index) is defined as the weighted root mean square error:

[0026] in To determine the effective number of monitoring points, For the first The weight of a point is set according to the accuracy of its monitoring method (e.g., GNSS weight is 1.0, InSAR weight is 0.7, and inclinometer weight is 0.9); if (With a preset tolerance, such as 5mm), the optimization algorithm is then started to update the parameter field; the optimization process can use a gradient-free method (such as genetic algorithm GA): the parameter field is discretized into several control points, and each individual is encoded as a group. The values ​​are used to generate a new population through crossover and mutation, and a forward model is called to evaluate fitness (i.e., ...). (value); to accelerate convergence, Bayesian optimization can also be used to construct a surrogate model based on historical sampling to predict the optimal parameter combination.

[0027] Force physical constraints to be applied after each update: ensure all , And mutual The parameter difference between adjacent elements satisfies:

[0028] in , The geological smoothing threshold can be set according to the soil type (e.g., a smaller value for cohesive soil, and a larger gradient allowed for gravel layers); in addition, the parameter values ​​are fixed at existing high-confidence borehole locations as anchor points to prevent the inversion results from deviating from the geological reality; iteration continues until... Or it may reach the maximum number of iterations.

[0029] The shear strength parameter spatial field output in step S5 is stored in the form of a three-dimensional grid. Each grid cell contains the corresponding cohesion value and internal friction angle value, and supports exporting to a general geological modeling format or a data file that can be recognized by mainstream geotechnical numerical analysis software. The output results also contain parameter uncertainty information, which is added to each grid cell in the form of standard deviation or confidence interval, and can be directly called by subsequent slope stability analysis modules, foundation bearing capacity verification programs or engineering digital twin platforms. After the inversion converges, the final shear strength parameter spatial field is stored in the form of a structured three-dimensional mesh, with each mesh element containing the corresponding cohesion value and internal friction angle value. The result can be exported to a general format (such as VTK, DXF, GOCAD or Plaxis material table) for direct use in subsequent geotechnical engineering analysis. To further support reliability design, the system can also output parameter uncertainty information. For example, when using Bayesian optimization or multi-startpoint inversion, the statistical distribution of each convergence result is recorded, the standard deviation or 95% confidence interval of each grid cell parameter is calculated, and stored along with the principal values. This uncertainty field can be used for Monte Carlo slope stability analysis, foundation bearing capacity probabilistic assessment, or risk warning modules in digital twin platforms.

[0030] Finally, in a specific embodiment, this method was applied to the deep overburden layer (approximately 90 meters thick) of the dam foundation of a hydropower station in Southwest China. A three-dimensional model was constructed using 32 boreholes and seismic refraction profiles, dividing the area into four geological units. The initial parameter field was generated using anisotropic random fields. The measured data included 18 Sentinel-1 InSAR deformation maps and data from 5 inclinometer boreholes. The inversion was performed using an improved particle swarm optimization algorithm, and the objective function converged to 3.2 mm (with an allowable error of 5 mm) after 11 iterations. The final inversion results showed a significant lateral decreasing trend in the internal friction angle of the sand layer, which was in good agreement with the results of subsequent verification borehole tests.

[0031] The above description is merely a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the protection scope of the present invention.

Claims

1. A method for determining the spatial variability of shear strength parameters of deep overburden layers, characterized in that: The determination method includes the following steps: S1. Construct a three-dimensional geological model of the overburden and initialize the spatial field of shear strength parameters; S2. Run the geotechnical numerical model based on the parameter spatial field to obtain the simulated displacement field; S3. Obtain the measured displacement data under the corresponding working conditions; S4. Compare the simulated displacement with the measured displacement. If the difference exceeds the threshold, update the parameter space field through the optimization algorithm and iterate the calculation until the difference meets the convergence condition. S5. Output the converged shear strength parameter spatial field as the final result.

2. The method for determining the spatial variability of shear strength parameters of deep overburden layers according to claim 1, characterized in that: The construction of the three-dimensional geological model of the overburden in step S1 includes establishing stratigraphic interfaces based on borehole columnar sections, engineering geological profiles, and geophysical exploration data, and dividing the stratigraphic units into multiple geological units with different genetic types or lithological characteristics. The initialization of the spatial field of shear strength parameters includes assigning initial mean values ​​and coefficients of variation of cohesion and internal friction angle to each geological unit, and generating a spatially correlated parameter distribution using a random field method. The spatial correlation distance of the random field is determined according to the depositional environment, grain composition, or consolidation history of each unit, and non-stationary covariance structure modeling is achieved through variogram analysis to reflect the variability differences of the thick overburden in the lateral and vertical directions.

3. The method for determining the spatial variability of shear strength parameters of deep overburden layers according to claim 1, characterized in that: In step S2, the geotechnical numerical model is a finite element model or a finite difference model. The boundary conditions of the model are set according to the actual terrain, groundwater level and external load. The material constitutive relation adopts the Mohr-Coulomb model, HardeningSoil model or modified Cambridge model. In the calculation process, the stress path evolution caused by construction sequence or water level changes is considered. The mesh division of the geotechnical numerical model adopts an adaptive refinement strategy, setting finer element sizes in areas with concentrated potential deformation. Special element types are used for the foundation interface of the structure, weak interlayers or soil-rock contact surface to accurately simulate slip, separation or shear behavior.

4. The method for determining the spatial variability of the shear strength parameter of a deep overburden layer according to claim 1, characterized in that: The measured displacement data in step S3 includes surface deformation data obtained through synthetic aperture radar interferometry, displacement time histories recorded by monitoring points of the Global Navigation Satellite System, settlement data obtained by leveling, deep horizontal displacement profiles measured by inclinometers, or compression of each soil layer recorded by layered settlement markers. The time nodes of the measured displacement data correspond to the loading stages of the simulated working conditions in the numerical model. Before use, the original monitoring data undergoes spatiotemporal registration, outlier removal, and multi-source data consistency verification to ensure that the displacements obtained by different monitoring methods are comparable under the same coordinate system and time reference.

5. The method for determining the spatial variability of shear strength parameters of deep overburden layers according to claim 1, characterized in that: The step S4, comparing the simulated displacement with the measured displacement, includes interpolating the simulated displacement field to the spatial position of the measured displacement data and calculating the displacement difference between the two at the corresponding monitoring points. The difference exceeding the threshold means that the weighted root mean square of the displacement difference is greater than the preset allowable error. The optimization algorithm includes a genetic algorithm, a particle swarm optimization algorithm, a Bayesian optimization method, or a gradient-based inversion algorithm, which is used to adjust the parameter values ​​of each grid node or region in the shear strength parameter spatial field. The weights are dynamically allocated according to the measurement accuracy, data stability, and spatial representativeness of the monitoring points. During the optimization process, some parameters of high-confidence borehole locations are retained as anchor points to prevent overall distortion of the parameter field.

6. The method for determining the spatial variability of the shear strength parameter of a deep overburden layer according to claim 1, characterized in that: In step S4, during the iterative calculation process, physical constraints are applied to the updated shear strength parameter spatial field, including that the cohesion value is not less than zero, the internal friction angle ranges from 0 to 45 degrees, and the parameter variation between adjacent spatial units does not exceed a preset geological smoothing threshold. The geological smoothing threshold is preset based on the stratigraphic continuity, sedimentary rhythm, or soil structural characteristics. For areas with clear geological interfaces, parameter abrupt changes are allowed, while within homogeneous soil layers, parameter gradient constraints are forcibly implemented to maintain spatial continuity.

7. The method for determining the spatial variability of shear strength parameters of deep overburden layers according to claim 1, characterized in that: The shear strength parameter spatial field output in step S5 is stored in the form of a three-dimensional grid. Each grid cell contains the corresponding cohesion value and internal friction angle value, and can be exported as a data file that can be recognized by general geological modeling format or mainstream geotechnical numerical analysis software. The output results also contain parameter uncertainty information, which is added to each grid cell in the form of standard deviation or confidence interval, and can be directly called by subsequent slope stability analysis modules, foundation bearing capacity verification programs or engineering digital twin platforms.

Images