A method, system and storage medium for optimizing a controlled drilling trajectory in thick alluvium

By introducing key control points and a three-dimensional geomechanical model into drilling trajectory optimization, and combining formation lithology and obstacle probability density functions, the drilling trajectory is optimized, solving the problems of rigid trajectory and lack of risk integration in existing methods. This achieves globally optimal drilling trajectory planning, reduces complexity, and improves safety and efficiency.

CN122452241APending Publication Date: 2026-07-24河南省地质研究院
View PDF 0 Cites 0 Cited by

Patent Information

Authority / Receiving Office
CN · China
Patent Type
Applications(China)
Current Assignee / Owner
河南省地质研究院
Filing Date
2026-05-07
Publication Date
2026-07-24

AI Technical Summary

Technical Problem

Existing drilling trajectory optimization methods cannot flexibly adapt to subtle changes in formation conditions, fail to effectively incorporate risk factors, leading to risks such as stuck drill pipe and wellbore collapse. Furthermore, they fail to establish real-time coupling between the trajectory path and the three-dimensional non-uniform geomechanical model, thus failing to achieve a globally optimal solution.

Method used

By setting key control points, a smooth trajectory is generated using a three-dimensional geomechanical model and spline interpolation function. Hard constraints are established by combining strata lithological parameters and the probability density function of underground obstacles. The total trajectory length and torque drag comprehensive evaluation value are optimized, and the optimal solution is found by using a particle swarm optimization algorithm.

Benefits of technology

It achieves global optimal drilling trajectory optimization under complex geological conditions, reduces the solution complexity of the optimization problem, ensures the drillability and safety of the optimized trajectory, and maximizes the overall drilling benefits.

✦ Generated by Eureka AI based on patent content.

Smart Images

  • Figure CN122452241A_ABST
    Figure CN122452241A_ABST
Patent Text Reader

Abstract

The application provides a thick alluvium controlled drilling trajectory optimization method, system and storage medium, the method comprises the following steps: generating a smooth candidate drilling trajectory by a first derivative continuous spline interpolation function, establishing a three-dimensional geomechanics model coupled with stratum parameters, and setting the wellbore stress meeting the rock strength criterion as a first hard constraint; performing probability modeling on geological risks to delimit boundaries, setting a safety distance as a second hard constraint, constructing a weighted spatial anisotropic drilling comprehensive cost objective function, fusing a trajectory length, a torque drag and a wellbore instability risk index, adopting a global optimization algorithm to minimize the comprehensive cost under the double constraints, solving optimal control point parameters, and interpolating to generate an optimized drilling trajectory.
Need to check novelty before this filing date? Find Prior Art

Description

Technical Field

[0001] This application belongs to the field of optimization, and in particular relates to a method, system and storage medium for controlled drilling trajectory optimization in thick alluvial layers. Background Technology

[0002] In engineering fields such as oil and gas exploration and development, geothermal utilization, and underground space development, the drilling trajectory is the only channel connecting the surface and underground target points. This is especially true when drilling operations are conducted in areas with complex geological conditions, such as thick alluvial layers, dramatic lithological variations, uneven stress distribution, unstable strata, and unknown obstacles. These conditions place extremely high demands on the optimization of the drilling trajectory. Existing drilling trajectory planning methods produce rigid trajectories that cannot flexibly adapt to subtle changes in formation conditions. Furthermore, the optimization objective function is often too simplistic, such as pursuing only the shortest path or minimum curvature, without considering risk factors like wellbore stability and torque drag. This can lead to risks such as stuck pipe, wellbore collapse, and excessive drill string fatigue during actual construction.

[0003] Some studies employ functions such as B-splines or Bézier curves to generate smoother, more continuous trajectories, enhancing their geometric flexibility. Expansion has also been made in terms of optimization objectives and constraints, such as incorporating torque-driven prediction models and collision avoidance scanning results. However, in terms of geomechanical coupling, most methods have not yet established real-time coupling between the trajectory path and a three-dimensional non-uniform geomechanical model. Wellbore stability analysis is often performed as an independent verification step after initial trajectory planning, rather than being integrated into the optimization iteration process as a hard constraint. This prevents the trajectory from actively avoiding areas of stress concentration or low rock strength, significantly reducing the global applicability of the optimization. Furthermore, the inherent uncertainty of geological exploration data, coupled with the inability of existing technologies to probabilistically assess and control risks, may lead to overly conservative or risky avoidance strategies. The methods fail to demonstrate that the relative importance of influencing factors varies across different formations and well sections, lacking the ability to identify the spatial anisotropy of drilling costs, thus failing to find a truly globally optimal solution. Summary of the Invention

[0004] To address the problem that existing methods fail to comprehensively consider the complex mechanical issues and risk factors arising from the interaction with the formation during drilling, and do not establish real-time coupling between the trajectory path and the three-dimensional non-uniform geomechanical model.

[0005] In the first aspect, the present invention proposes a method for optimizing controlled drilling trajectories in thick alluvial layers, comprising the following steps: At least three key control points are set, and the three-dimensional spatial coordinates of each key control point and the tangential vector of the well trajectory at the control point are used as optimization variables. A smooth candidate drilling trajectory is generated based on the optimization variables through a preset spline interpolation function that guarantees at least the continuity of the first derivative. A three-dimensional geomechanical model coupled with formation lithology and physical parameters is established. The first hard constraint condition is that the induced stress state around the well at any point on the candidate trajectory does not exceed the rock strength envelope defined based on the preset failure criterion. The spatial uncertainty of underground obstacles and high-risk geological bodies is modeled as a three-dimensional probability density function. The boundary of the probability exclusion body is delineated according to the preset risk tolerance. The second hard constraint condition is that the candidate drilling trajectory is kept at a distance greater than the preset safety distance from the boundary. The weighted sum of the total trajectory length, the comprehensive evaluation value of torque drag, and the comprehensive evaluation value of wellbore instability risk is used as the drilling comprehensive cost objective function. The weight coefficients of each term are determined according to the formation unit attributes and local drilling direction along the trajectory. With the goal of minimizing the drilling comprehensive cost objective function, the optimal set of optimization variables is searched and determined in the solution space that simultaneously satisfies the first and second hard constraints. Based on the optimal optimization variables, the optimized drilling trajectory is generated through the spline interpolation function.

[0006] Optionally, generating a smooth candidate drilling trajectory based on the optimization variables using a preset spline interpolation function that guarantees at least the continuity of the first derivative includes: The spline interpolation function is a third-order Bessel-Hermitian spline interpolation function. It uses dimensionless parameters as parameterized variables, takes the three-dimensional spatial coordinates of each key control point as the interpolation point, and takes the tangential vector of the interpolation point as the first derivative of the dimensionless parameters at the interpolation point. It interpolates the well section between two adjacent key control points, and splices all the interpolated well sections to generate a globally continuous smooth drilling trajectory.

[0007] Optionally, the establishment of a three-dimensional geomechanical model coupled with formation lithological parameters, using the condition that the well-circumferential induced stress state at any point on the candidate trajectory does not exceed the rock strength envelope defined based on a preset failure criterion, as the first hard constraint, includes: The first hard constraint condition includes three-dimensional geological meshing of the well area, assigning formation rock mechanics parameters obtained from well logging interpretation or seismic inversion to each mesh cell, including elastic modulus, Poisson's ratio, cohesion, and internal friction angle, calculating the stress-strain distribution of the three-dimensional geological body under the action of the in-situ stress field using the finite element method, calculating the induced stress field of the wellbore rock mass caused by drilling activity for any point on the candidate trajectory, calculating the critical collapse pressure using the Mohr-Coulomb failure criterion, and calculating the critical rupture pressure using the maximum tensile stress criterion, ensuring that the wellbore fluid column pressure corresponding to the stress field is completely within the stability window defined by the critical collapse pressure and rupture pressure.

[0008] Optionally, the weighted sum of the total trajectory length, the comprehensive evaluation value of torque drag, and the comprehensive evaluation value of wellbore instability risk is used as the objective function for overall drilling cost. The weight coefficients of each term are determined based on the formation unit attributes along the trajectory and the local drilling direction, including: The drilling comprehensive cost objective function is defined as a weighted sum of global macro-level evaluation terms and local micro-level risk evaluation terms, and its calculation formula is as follows: ,in: The objective function is the overall drilling cost. The total length of the trajectory after dimensionless transformation; The dimensionless global torque drag comprehensive evaluation value is calculated based on the vector iteration from the drill bit to the wellhead; N is the total number of micro-segments in the discrete division of the drilling trajectory; The dimensionless wellbore instability risk assessment value calculated based on the local stress field for the i-th micro-element segment; The global weight coefficient is the total length of the trajectory. This is the global weighting coefficient for the overall evaluation value of global torque drag. The local dynamic weighting coefficient for the wellbore instability risk of the i-th micro-element segment; The local dynamic weight coefficient The specific values ​​are determined by looking up a table using a pre-set engineering experience matrix, based on the lithology, permeability, and clay content of the formation unit traversed by the i-th trajectory micro-element, as well as the angle between the well trajectory at that micro-element and the formation bedding plane normal and the direction of the maximum horizontal principal stress; and for the entire well section, the weights of each item satisfy the normalization condition: .

[0009] Optionally, the step of modeling the spatial location uncertainty of underground obstacles and high-risk geological bodies as a three-dimensional probability density function, delineating the boundary of the probability exclusion body according to a preset risk tolerance, and maintaining a candidate drilling trajectory at a distance greater than a preset safety distance from the boundary as a second hard constraint condition includes: The second hard constraint condition includes modeling the spatial location uncertainty of each underground obstacle or high-risk geological body as a three-dimensional Gaussian distribution function using seismic interpretation data, setting a risk tolerance, defining the boundary of the integral region where the integral probability of the Gaussian distribution function in the spatial region reaches the tolerance value as the boundary of the probability exclusion body, and requiring that the minimum Euclidean distance from any point on the planned drilling trajectory to the boundary is greater than a preset safety distance value.

[0010] Optionally, the step of searching and determining an optimal set of optimization variables within the solution space that simultaneously satisfies the first and second hard constraints, with the objective function of minimizing the overall drilling cost, includes: The particle swarm optimization algorithm is adopted. The three-dimensional coordinates and the Cartesian component of the tangential vector of each key control point are used to form a high-dimensional optimization vector as a particle in the algorithm. In each iteration, for the candidate drilling trajectory represented by each particle, it is checked whether the candidate drilling trajectory satisfies the first and second hard constraints. If either condition is not satisfied, the cost objective function value of the particle is assigned a preset maximum value by the penalty function method. The velocity and position of the particle are updated iteratively until the maximum number of iterations is reached or the objective function value converges. The optimization variable corresponding to the globally optimal particle is output.

[0011] Optionally, the calculation of the torque drag comprehensive evaluation value includes discretizing the drilling trajectory into N micro-segments, using a soft rope model, and performing mechanical analysis segment by segment upwards from the drill bit. In each micro-segment, the weight, buoyancy, normal pressure generated by contact with the well wall, and friction force calculated based on the formation friction coefficient of the micro-segment drill string are comprehensively considered. The total tension and total torque at the wellhead are obtained through vector iteration calculation. The calculated total tension and total torque are then divided by preset reference tension and reference torque to perform dimensionless calculation, and the torque drag comprehensive evaluation value is generated by combining them according to a preset weighting method.

[0012] In another aspect, the present invention also proposes a controlled drilling trajectory optimization system for thick alluvial deposits, comprising the following modules: The generation module is used to set at least three key control points, and use the three-dimensional spatial coordinates of each key control point and the tangential vector of the well trajectory at the control point as optimization variables. A smooth candidate drilling trajectory is generated based on the optimization variables through a preset spline interpolation function that guarantees at least the continuity of the first derivative. A module is established to create a three-dimensional geomechanical model coupled with formation lithological and physical parameters. The first hard constraint condition is that the induced stress state around the well at any point on the candidate trajectory does not exceed the rock strength envelope defined based on a preset failure criterion. The spatial uncertainty of underground obstacles and high-risk geological bodies is modeled as a three-dimensional probability density function. The boundary of the probability exclusion body is delineated according to a preset risk tolerance. The second hard constraint condition is that the candidate drilling trajectory is kept at a distance greater than a preset safety distance from the boundary. The construction module is used to take the weighted sum of the total trajectory length, the comprehensive evaluation value of torque drag, and the comprehensive evaluation value of wellbore instability risk as the drilling comprehensive cost objective function. The weight coefficients of each item are determined according to the formation unit attributes and local drilling direction along the trajectory. With the goal of minimizing the drilling comprehensive cost objective function, the module searches and determines the optimal set of optimization variables in the solution space that simultaneously satisfies the first and second hard constraints. Based on the optimal optimization variables, the module generates an optimized drilling trajectory through the spline interpolation function.

[0013] Preferably, the step of generating a smooth candidate drilling trajectory based on the optimization variables using a preset spline interpolation function that guarantees at least the continuity of the first derivative includes: The spline interpolation function is a third-order Bessel-Hermitian spline interpolation function. It uses dimensionless parameters as parameterized variables, takes the three-dimensional spatial coordinates of each key control point as the interpolation point, and takes the tangential vector of the interpolation point as the first derivative of the dimensionless parameters at the interpolation point. It interpolates the well section between two adjacent key control points, and splices all the interpolated well sections to generate a globally continuous smooth drilling trajectory.

[0014] Preferably, the establishment of a three-dimensional geomechanical model coupled with formation lithological parameters, wherein the well-circumferential induced stress state at any point on the candidate trajectory does not exceed the rock strength envelope defined based on a preset failure criterion, is used as the first hard constraint condition, including: The first hard constraint condition includes three-dimensional geological meshing of the well area, assigning formation rock mechanics parameters obtained from well logging interpretation or seismic inversion to each mesh cell, including elastic modulus, Poisson's ratio, cohesion, and internal friction angle, calculating the stress-strain distribution of the three-dimensional geological body under the action of the in-situ stress field using the finite element method, calculating the induced stress field of the wellbore rock mass caused by drilling activity for any point on the candidate trajectory, calculating the critical collapse pressure using the Mohr-Coulomb failure criterion, and calculating the critical rupture pressure using the maximum tensile stress criterion, ensuring that the wellbore fluid column pressure corresponding to the stress field is completely within the stability window defined by the critical collapse pressure and rupture pressure.

[0015] Preferably, the weighted sum of the total trajectory length, the comprehensive evaluation value of torque drag, and the comprehensive evaluation value of wellbore instability risk is used as the objective function for overall drilling cost. The weight coefficients of each term are determined based on the formation unit attributes along the trajectory and the local drilling direction, including: The drilling comprehensive cost objective function is defined as a weighted sum of global macro-level evaluation terms and local micro-level risk evaluation terms, and its calculation formula is as follows: ,in: The objective function is the overall drilling cost. The total length of the trajectory after dimensionless transformation; The dimensionless global torque drag comprehensive evaluation value is calculated based on the vector iteration from the drill bit to the wellhead; N is the total number of micro-segments in the discrete division of the drilling trajectory; The dimensionless wellbore instability risk assessment value calculated based on the local stress field for the i-th micro-element segment; The global weight coefficient is the total length of the trajectory. This is the global weighting coefficient for the overall evaluation value of global torque drag. The local dynamic weighting coefficient for the wellbore instability risk of the i-th micro-element segment; The local dynamic weight coefficient The specific values ​​are determined by looking up a table using a pre-set engineering experience matrix, based on the lithology, permeability, and clay content of the formation unit traversed by the i-th trajectory micro-element, as well as the angle between the well trajectory at that micro-element and the formation bedding plane normal and the direction of the maximum horizontal principal stress; and for the entire well section, the weights of each item satisfy the normalization condition: .

[0016] Preferably, the step of modeling the spatial location uncertainty of underground obstacles and high-risk geological bodies as a three-dimensional probability density function, delineating the boundary of the probability exclusion body according to a preset risk tolerance, and maintaining a candidate drilling trajectory at a distance greater than a preset safety distance from the boundary as a second hard constraint condition includes: The second hard constraint condition includes modeling the spatial location uncertainty of each underground obstacle or high-risk geological body as a three-dimensional Gaussian distribution function using seismic interpretation data, setting a risk tolerance, defining the boundary of the integral region where the integral probability of the Gaussian distribution function in the spatial region reaches the tolerance value as the boundary of the probability exclusion body, and requiring that the minimum Euclidean distance from any point on the planned drilling trajectory to the boundary is greater than a preset safety distance value.

[0017] Preferably, the step of searching and determining an optimal set of optimization variables within the solution space that simultaneously satisfies the first and second hard constraints, with the objective function of minimizing the overall drilling cost, includes: The particle swarm optimization algorithm is adopted. The three-dimensional coordinates and the Cartesian component of the tangential vector of each key control point are used to form a high-dimensional optimization vector as a particle in the algorithm. In each iteration, for the candidate drilling trajectory represented by each particle, it is checked whether the candidate drilling trajectory satisfies the first and second hard constraints. If either condition is not satisfied, the cost objective function value of the particle is assigned a preset maximum value by the penalty function method. The velocity and position of the particle are updated iteratively until the maximum number of iterations is reached or the objective function value converges. The optimization variable corresponding to the globally optimal particle is output.

[0018] Preferably, the calculation of the torque drag comprehensive evaluation value includes discretizing the drilling trajectory into N micro-segments, using a soft rope model, and performing mechanical analysis segment by segment upwards from the drill bit. In each micro-segment, the weight, buoyancy, normal pressure generated by contact with the well wall, and friction force calculated based on the formation friction coefficient of the micro-segment drill string are comprehensively considered. The total tension and total torque at the wellhead are obtained through vector iteration calculation. The calculated total tension and total torque are then divided by preset reference tension and reference torque to perform dimensionless calculation, and the torque drag comprehensive evaluation value is generated by combining them according to a preset weighting method.

[0019] This invention establishes a refined geomechanical model coupled with formation parameters and integrates wellbore stability and obstacle avoidance issues with location uncertainties. The obstacle avoidance problem is used as a rigid constraint in the optimization process, ensuring the drillability and safety of the optimized trajectory. A comprehensive drilling cost objective function considering formation spatial anisotropy is constructed, representing and comprehensively evaluating trajectory length, engineering risks, and operational difficulty. This ensures that the optimization result is no longer limited to the geometric shortest but maximizes the overall drilling benefits. The use of a key control point parameterization method reduces the solution complexity of the optimization problem. Combined with a global optimization algorithm, it can obtain the globally optimal drilling trajectory under complex constraints. Attached Figure Description

[0020] Figure 1 A flowchart of a method for optimizing controlled drilling trajectories in thick alluvial deposits provided by this invention; Figure 2 This is a schematic diagram of the cumulative variation curve of the pull force and torque along the well section of the soft rope model provided by the present invention. Detailed Implementation

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

[0022] In a first aspect, the present invention proposes a method for optimizing controlled drilling trajectories in thick alluvial deposits, such as... Figure 1 As shown, it includes the following steps: S1. Set at least three key control points, and use the three-dimensional spatial coordinates of each key control point and the tangential vector of the well trajectory at the control point as optimization variables. Then, generate a smooth candidate drilling trajectory based on the optimization variables using a preset spline interpolation function that guarantees at least the continuity of the first derivative. Define the optimization variable as a vector containing the three-dimensional coordinates of all key control points. and tangential vector For N key control points, the dimension of the optimization variables is 6N. In one embodiment, a piecewise cubic Hermitian interpolation polynomial is selected as the spline interpolation function for any two adjacent key control points. and Using position coordinates and tangential vector and This generates the three-dimensional spatial curve parametric equations for that segment. Connecting all the segmented curves sequentially yields a smooth candidate drilling trajectory with a continuous global first derivative.

[0023] In some embodiments, generating a smooth candidate drilling trajectory based on the optimization variables using a preset spline interpolation function that guarantees at least continuity of the first derivative includes: The spline interpolation function is a third-order Bessel-Hermitian spline interpolation function. It uses dimensionless parameters as parameterized variables, takes the three-dimensional spatial coordinates of each key control point as the interpolation point, and takes the tangential vector of the interpolation point as the first derivative of the dimensionless parameters at the interpolation point. It interpolates the well section between two adjacent key control points, and splices all the interpolated well sections to generate a globally continuous smooth drilling trajectory.

[0024] For two adjacent critical control points and For any defined well section, input a local dimensionless parameter u, with a value range of [0,1]. The geometry of this well section is determined by the coordinates of its endpoints. , and endpoint tangential vector , The location vector R(u) of any point on the well section is uniquely determined. It is calculated using the following third-order Bessel-Hermitian spline interpolation formula: .in, Let i be the three-dimensional coordinates of the i-th key control point. This is the tangential vector at that point. The magnitude of this vector affects the curvature of the curve segment, and the preferred value is usually related to the well section length, for example, it can be set to 0.5 to 1.5 times the distance between the two control points. For example, for a well from... =(0,0,1000) to For the well section with coordinates (200, 100, 1200), the tangential vector can be set as... =(150,75,150), =(250,125,250). Using this formula, when u=0, , When u=1, , All N-1 well segments generated by N key control points are sequentially spliced ​​together. Since the positions and first derivatives of two adjacent well segments are equal at each key control point, the global positional and tangential continuity of the entire drilling trajectory is ensured, that is, the trajectory is smooth and without inflections.

[0025] S2. Establish a three-dimensional geomechanical model coupled with the lithological and physical parameters of the formation. The first hard constraint condition is that the induced stress state around the well at any point on the candidate trajectory does not exceed the rock strength envelope defined based on the preset failure criterion. The spatial uncertainty of underground obstacles and high-risk geological bodies is modeled as a three-dimensional probability density function. The boundary of the probability exclusion body is delineated according to the preset risk tolerance. The second hard constraint condition is that the candidate drilling trajectory is kept at a distance greater than the preset safety distance from the boundary. Using finite element analysis software such as Abaqus or ANSYS, a three-dimensional geomechanical mesh model incorporating in-situ stress, pore pressure, and rock mechanical parameter fields is established based on pre-drilling geological exploration data. The candidate well trajectory is discretized into a series of evaluation points. For each evaluation point, the far-field stress tensor is extracted from the geomechanical model. For homogeneous formation elements, the Kirsch equation is used to calculate the circumferential, radial, and axial induced stresses around the wellbore at that point. The Kirsch equation is based on the elastic mechanics assumptions of homogeneous isotropic elastic rock mass, circular wellbore, and plane strain. First, the far-field three-dimensional in-situ stress, vertical principal stress, and maximum / minimum horizontal principal stress are converted into stress components in the local coordinate system of the wellbore, combined with the well inclination angle and azimuth angle. Then, parameters such as wellbore radius, pore pressure, and wellbore fluid column pressure are substituted, and the radial, circumferential, and axial induced stresses at the wellbore location are directly calculated using analytical formulas. For non-homogeneous formation elements, a formation non-homogeneity correction factor needs to be introduced based on the Kirsch equation. The value ranges from 0.7 to 1.2, and is determined by the formation elastic modulus gradient. The greater the spatial difference in formation elastic modulus, the stronger the heterogeneity. The further it deviates from the neutral value of 1. In one embodiment, ,in Let be the magnitude of the elastic modulus gradient vector along the normal or tangential direction of a micro-element of a specific drilling trajectory. This is the average elastic modulus within the grid containing the infinitesimal element. As an empirical calibration constant, it needs to be obtained by multinomial regression or inversion fitting between the location of lithological abrupt change interfaces in regional historical logging data and the frequency of actual wellbore collapse accidents. Its value is typically between 0.1 and 0.5. By substituting the modified Kirch equation, the calculation results of radial, circumferential, and axial induced stresses around the wellbore can be corrected, so that the stress solution results are consistent with the real mechanical response of non-uniform formations, ensuring the accuracy of wellbore stability constraint judgment.

[0026] The Mohr-Coulomb criterion or the maximum tensile stress criterion is used as the rock failure criterion to determine whether the calculated wellbore induced stress is within the rock strength envelope. If the stress at any point exceeds the envelope, the trajectory violates the constraint. The center position and uncertainty of each obstacle or high-risk geological body are represented by a probability density function of a three-dimensional Gaussian distribution, defined by the expected position vector and the covariance matrix. Based on a preset risk tolerance, such as a 99.9% collision avoidance probability, an equivalent confidence ellipsoid is determined using a chi-square distribution as the boundary of the probabilistic repulsion body. The GJK algorithm is used to calculate the shortest distance between the candidate drilling trajectory curve and the surface of the ellipsoid. The continuous one-dimensional drilling trajectory generated by spline interpolation is discretized into a series of extremely short straight line segments at fixed steps (e.g., every 5 or 10 meters) less than a preset safety distance. Using each straight line segment as a central axis, and considering the drill string's outer diameter and geometric safety margin, it is enclosed in a series of convex spatial bounding boxes, such as the polygonal approximate convex hull of a cylinder, or a directed bounding box (OBB). The probabilistic exclusionary body boundary of the high-risk underground geological body is also meshed into a set of convex polyhedra. The trajectory convex bounding boxes and the geological convex polyhedra are substituted into the GJK algorithm for iterative collision detection to calculate the minimum safe Euclidean distance. It is then determined whether this shortest distance is greater than the preset safety distance value; if it is less than or equal to this value, the trajectory violates the constraints.

[0027] In some embodiments, establishing a three-dimensional geomechanical model coupled with formation lithological parameters, and using the condition that the well-circumferential induced stress state at any point on the candidate trajectory does not exceed the rock strength envelope defined based on a preset failure criterion as the first hard constraint, includes: The first hard constraint condition includes three-dimensional geological meshing of the well area, assigning formation rock mechanics parameters obtained from well logging interpretation or seismic inversion to each mesh cell, including elastic modulus, Poisson's ratio, cohesion, and internal friction angle, calculating the stress-strain distribution of the three-dimensional geological body under the action of the in-situ stress field using the finite element method, calculating the induced stress field of the wellbore rock mass caused by drilling activity for any point on the candidate trajectory, calculating the critical collapse pressure using the Mohr-Coulomb failure criterion, and calculating the critical rupture pressure using the maximum tensile stress criterion, ensuring that the wellbore fluid column pressure corresponding to the stress field is completely within the stability window defined by the critical collapse pressure and rupture pressure.

[0028] The target work area is divided into a fine three-dimensional geological grid, and specific rock mechanics parameters are assigned to each grid cell. For example, for a sandstone stratum cell, the elastic modulus E can be set to 15 GPa and the Poisson's ratio to... The internal friction angle is 0.25, the cohesion C is 8 MPa, and the internal friction angle is... The angle is 35°; for mudstone stratigraphic units, E can be set to 5 GPa. The value is 0.35, and C is 3 MPa. The angle is 22°. Boundary conditions are set, namely the in-situ stress field, including the vertical principal stresses. Maximum horizontal principal stress and minimum horizontal principal stress The pore pressure field and the initial stress tensor of each grid cell were calculated using finite element analysis software on the entire three-dimensional geological model. For a candidate well trajectory, the trajectory was discretized into a series of calculation points. For each calculation point, based on the geomechanical parameters and stress state of the grid cell containing the point, as well as the well inclination and azimuth angles, the radial stress at the wellbore was first obtained by modifying the Kirch equation. Circumferential stress Axial stress The effective stress components are obtained by subtracting pore pressure; the effective normal stress and shear stress are then substituted into the Mohr-Coulomb criterion. The minimum wellbore fluid column pressure required for shear collapse is the critical collapse pressure. Then, the circumferential effective stress is used to achieve the tensile strength of the rock. The determination condition is based on the maximum tensile stress criterion. The maximum wellbore fluid column pressure at which tensile fracturing occurs is determined by reverse calculation; this is the critical fracturing pressure. The actual fluid column pressure in the wellbore needs to be between and Only by maintaining a certain range can wellbore stability be guaranteed. For any point on the trajectory, the safe pressure window... It must be greater than a preset threshold, such as 1.5 MPa, to ensure that drilling operations have sufficient safety margin.

[0029] In some embodiments, modeling the spatial location uncertainty of underground obstacles and high-risk geological bodies as a three-dimensional probability density function, delineating the boundary of the probability exclusion body according to a preset risk tolerance, and maintaining a candidate drilling trajectory at a distance greater than a preset safety distance from the boundary as a second hard constraint condition includes: The second hard constraint condition includes modeling the spatial location uncertainty of each underground obstacle or high-risk geological body as a three-dimensional Gaussian distribution function using seismic interpretation data, setting a risk tolerance, defining the boundary of the integral region where the integral probability of the Gaussian distribution function in the spatial region reaches the tolerance value as the boundary of the probability exclusion body, and requiring that the minimum Euclidean distance from any point on the planned drilling trajectory to the boundary is greater than a preset safety distance value.

[0030] For a known underground obstacle, such as the trajectory of an adjacent well, the inclinometer data exhibits cumulative error. The location uncertainty is modeled as a three-dimensional Gaussian distribution. The center of the probability density function... The nominal position calculated from the inclinometer data, and the covariance matrix. The statistical characteristics of inclinometer errors are used to determine, for example, assigning values ​​to horizontal and vertical errors respectively, to ensure that the uncertainty modeling closely reflects reality. A risk tolerance is then set. The preferred range is to The boundary of the probabilistic repulsive body is an isoprobability density surface. Within the volume enclosed by this surface, the cumulative probability of the obstacle's presence is... For a three-dimensional Gaussian distribution, the boundary is an ellipsoid, and the equation is: Where v is the coordinate of a spatial point, and the k value is determined by a chi-square distribution with 3 degrees of freedom at a confidence level of 1. Time to determine. Set a physical safety distance. This value is typically taken as several times the maximum outer diameter of the drill string assembly. During optimization, the shortest Euclidean distance from the entire candidate drilling trajectory curve to the boundary of the aforementioned probabilistic repulsion ellipsoid is calculated based on the GJK algorithm. Only when the minimum distance from all points on the entire trajectory to the repulsion ellipsoid boundary of all obstacles is strictly greater than... Only when this condition is met is the trajectory considered to satisfy the second hard constraint.

[0031] S3, the weighted sum of the total trajectory length, torque drag comprehensive evaluation value and wellbore instability risk comprehensive evaluation value is used as the drilling comprehensive cost objective function. The weight coefficients of each item are determined according to the formation unit attributes and local drilling direction along the trajectory. With the goal of minimizing the drilling comprehensive cost objective function, the optimal set of optimization variables, namely the three-dimensional spatial coordinates and tangential vectors of each key control point, are searched and determined in the solution space that simultaneously satisfies the first and second hard constraints. Based on the optimal optimization variables, the optimized drilling trajectory is generated through the spline interpolation function.

[0032] The total arc length was calculated using numerical integration via the Gauss-Gande quadrature method for the candidate trajectories. A flexible rope model was applied, discretizing the trajectory into elements. Taking into account wellbore curvature, inclination angle, formation friction coefficient, and drill string assembly parameters, the total torque and drag force were iteratively calculated. The flexible rope model is a simplified mechanical model for drilling torque drag analysis. It assumes the drill string is a flexible rope with no bending stiffness, transmitting only axial tension, torque, and friction, neglecting the drill string's own bending resistance. The contact normal force and friction between the drill string and the wellbore were iteratively calculated segment by segment along the wellbore trajectory to solve for the total tension and torque at the wellhead. The total torque and drag force were then normalized to a global torque drag comprehensive evaluation value. A wellbore stability safety factor, which is the ratio of rock strength to the maximum stress around the well, was calculated along the trajectory line. A weight database corresponding to the 3D geological model mesh is established. Based on the lithology of the formation where the micro-segment on the trajectory is located, such as abrasiveness and stability, as well as the well inclination and azimuth angles of that micro-segment, the unique local instability penalty weight coefficient for that micro-segment is queried and interpolated. This local instability penalty weight coefficient is multiplied by the reciprocal of the corresponding safety factor, and then integrated along the global trajectory to obtain the comprehensive evaluation value of wellbore instability risk. Simultaneously, the weight coefficients for the total trajectory length and torque drag are preset according to engineering design requirements. A genetic algorithm or particle swarm optimization algorithm is used as the global optimization algorithm. In each iteration of the algorithm, for each generated candidate solution (i.e., a set of optimization variables), a hard constraint check is first performed. If the constraint is violated, a maximum objective function value matching the magnitude of the dimensionless objective function is given using the penalty function method. If the constraint is satisfied, the normalized total arc length and the comprehensive evaluation value of global torque drag are multiplied by their corresponding weight coefficients, and then added to the comprehensive evaluation value of wellbore instability risk obtained from the aforementioned integration to calculate the scalarized comprehensive drilling cost objective function value. The algorithm continuously searches for the solution that minimizes the objective function value through operations such as selection, crossover, mutation, and particle velocity and position updates until the maximum number of iterations or the convergence criterion is reached. The optimal solution output at termination is the optimal set of key control point positions and tangential vectors. Substituting this set of variables into the aforementioned piecewise cubic Hermitian interpolation polynomial generates the optimized drilling trajectory.

[0033] In some embodiments, the weighted sum of the total trajectory length, the comprehensive evaluation value of torque drag, and the comprehensive evaluation value of wellbore instability risk is used as the objective function for overall drilling cost. The weight coefficients of each term are determined based on the formation unit attributes along the trajectory and the local drilling direction, including: The drilling comprehensive cost objective function is defined as a weighted sum of global macro-level evaluation terms and local micro-level risk evaluation terms, and its calculation formula is as follows: ,in: The objective function is the overall drilling cost. The total length of the trajectory after dimensionless transformation; The dimensionless global torque drag comprehensive evaluation value is calculated based on the vector iteration from the drill bit to the wellhead; N is the total number of micro-segments in the discrete division of the drilling trajectory; The dimensionless wellbore instability risk assessment value calculated based on the local stress field for the i-th micro-element segment; The global weight coefficient is the total length of the trajectory. This is the global weighting coefficient for the overall evaluation value of global torque drag. The local dynamic weighting coefficient for the wellbore instability risk of the i-th micro-element segment; The local dynamic weight coefficient The specific values ​​are determined by looking up a table using a pre-set engineering experience matrix, based on the lithology, permeability, and clay content of the formation unit traversed by the i-th trajectory micro-element, as well as the angle between the well trajectory at that micro-element and the formation bedding plane normal and the direction of the maximum horizontal principal stress; and for the entire well section, the weights of each item satisfy the normalization condition: .

[0034] Dimensionless processing is performed on each sub-item: Total length of the global trajectory Normalization is achieved by dividing the total length of the generated trajectory by the straight-line distance from the wellhead to the target point; global torque drag comprehensive evaluation value. The total torque and total tension calculated iteratively up to the wellhead are divided by the maximum rated torque and maximum rated load of the drilling rig, respectively, and then normalized by a weighted combination according to a preset ratio; the local wellbore instability risk assessment value of the i-th micro-element segment. The stability of this micro-element segment is characterized by the reciprocal of its wellbore stability safety factor, and is dimensionless. The smaller the safety factor, the better. The higher the value, the higher the risk of instability in that local stratum.

[0035] Global weight coefficient and These are the baseline values ​​set based on overall macro-constraints of the drilling project, such as schedule requirements and the upper limit of drilling rig performance; while the local dynamic weighting coefficients... This relies on a pre-defined five-dimensional engineering experience matrix, whose dimensions include: lithology, permeability grade, angle between the drilling trajectory and the formation normal, angle between the drilling trajectory and the maximum horizontal principal stress, and clay content. Each element of the matrix outputs a baseline risk weight characterizing local geological sensitivity. For example, when a micro-segment of the trajectory passes through mudstone formations with high clay content and the drilling direction has a large angle with the formation normal, the risk of local wellbore instability is extremely high, and the baseline risk weight obtained from the table is extremely large; conversely, if the trajectory drills at a near-vertical angle in stable limestone, the baseline risk weight obtained from the table is extremely small. The global baseline value and the risk weights obtained from the table for all micro-segments are linearly scaled together to ensure that the final weights of all items strictly satisfy the condition that the sum is 1. The value integrates the optimization of global macro costs (trajectory length, total friction) and the dynamic punishment and avoidance of local wellbore instability risks in high-risk formations.

[0036] In some embodiments, the step of searching and determining an optimal set of optimization variables within the solution space that simultaneously satisfies the first and second hard constraints, with the objective function of minimizing the overall drilling cost, includes: The particle swarm optimization algorithm is adopted. The three-dimensional coordinates and the Cartesian component of the tangential vector of each key control point are used to form a high-dimensional optimization vector as a particle in the algorithm. In each iteration, for the candidate drilling trajectory represented by each particle, it is checked whether the candidate drilling trajectory satisfies the first and second hard constraints. If either condition is not satisfied, the cost objective function value of the particle is assigned a preset maximum value by the penalty function method. The velocity and position of the particle are updated iteratively until the maximum number of iterations is reached or the objective function value converges. The optimization variable corresponding to the globally optimal particle is output.

[0037] Assuming four key control points are defined, the positions of the two middle control points and the tangential vectors of all four points are fixed, while the positions of the wellhead and the target point are fixed. If the tangential vectors are represented using Cartesian coordinate components, the total number of optimization variables is 18. Each particle is an 18-dimensional vector representing a complete set of trajectory definition parameters. The particle swarm optimization (PSO) algorithm's operating parameters are set as follows: particle population size is set to 100, maximum number of iterations is set to 500, inertia weight w decreases linearly from 0.9 to 0.4, and the learning factor... and All values ​​are set to 2.0. In each iteration, for each particle in the population, a unique candidate drilling trajectory is generated based on the value of the particle's 18-dimensional vector. Constraint checks are performed on this trajectory. The cost objective function value of the candidate drilling trajectory is calculated. If the constraint check passes, the comprehensive drilling cost is calculated normally. If any constraint is not satisfied, a penalty function is applied, and... The value is assigned to a maximum value, for example According to the standard update rules of the particle swarm optimization algorithm, update the velocity and position of each particle: ; .in It is the particle's own historical optimal position. This is the globally optimal position for the entire population. This iterative process continues until the number of iterations reaches 500, or the globally optimal solution shows no substantial improvement in 50 consecutive iterations, at which point the algorithm is considered converged. The globally optimal particle... The 18 optimization variables included were identified as the optimal key control point parameters, and the optimized drilling trajectory was generated using these parameters.

[0038] In some embodiments, the calculation of the torque drag comprehensive evaluation value includes discretizing the drilling trajectory into N micro-segments, using a soft rope model, and performing mechanical analysis segment by segment upwards from the drill bit. In each micro-segment, the weight, buoyancy, normal pressure generated by contact with the well wall, and friction force calculated based on the formation friction coefficient of the micro-segment drill string are comprehensively considered. The total tension and total torque at the wellhead are obtained through vector iteration calculation. The calculated total tension and total torque are then divided by preset reference tension and reference torque to perform dimensionless calculation, and the torque drag comprehensive evaluation value is generated by combining them according to a preset weighting method.

[0039] Discretize the candidate drilling trajectory with a total length of L into N paths with lengths of L and L. The calculation starts from the smallest infinitesimal segment N where the drill bit is located and iterates upwards segment by segment to the wellhead. In the i-th segment, the input is the tension at the bottom of that segment. and torque For drilling operations, the tension at the top of this section... The calculation formula is: ,in This represents the unit weight of the drilling tool in the drilling fluid. This represents the average wellbore inclination angle for this section. The coefficient of friction, This is positive pressure. Positive pressure Determined by the components of wellbore curvature and drill string weight, it can be approximately calculated as follows: ,in The curvature angle of this segment. The torque at the top of this segment. The calculation formula is: , where r is the drill string joint radius. The initial condition for iteration is: at the drill bit... =0, =0. The total pulling force at the wellhead is obtained through iteration from i=N to i=1. and total torque ,like Figure 2 As shown, dimensionless processing is performed, and a reference tension is set. and reference torque The dimensionless torque drag comprehensive evaluation value TD is: .

[0040] On the other hand, the present invention also provides a controlled drilling trajectory optimization system for thick alluvial layers, comprising the following modules: The generation module is used to set at least three key control points, and use the three-dimensional spatial coordinates of each key control point and the tangential vector of the well trajectory at the control point as optimization variables. A smooth candidate drilling trajectory is generated based on the optimization variables through a preset spline interpolation function that guarantees at least the continuity of the first derivative. A module is established to create a three-dimensional geomechanical model coupled with formation lithological and physical parameters. The first hard constraint condition is that the induced stress state around the well at any point on the candidate trajectory does not exceed the rock strength envelope defined based on a preset failure criterion. The spatial uncertainty of underground obstacles and high-risk geological bodies is modeled as a three-dimensional probability density function. The boundary of the probability exclusion body is delineated according to a preset risk tolerance. The second hard constraint condition is that the candidate drilling trajectory is kept at a distance greater than a preset safety distance from the boundary. The construction module is used to take the weighted sum of the total trajectory length, the comprehensive evaluation value of torque drag, and the comprehensive evaluation value of wellbore instability risk as the drilling comprehensive cost objective function. The weight coefficients of each item are determined according to the formation unit attributes and local drilling direction along the trajectory. With the goal of minimizing the drilling comprehensive cost objective function, the module searches and determines the optimal set of optimization variables in the solution space that simultaneously satisfies the first and second hard constraints. Based on the optimal optimization variables, the module generates an optimized drilling trajectory through the spline interpolation function.

[0041] The various embodiments in this specification are described in a progressive manner. Each embodiment focuses on the differences from other embodiments. The various embodiments can be combined as needed, and the same or similar parts can be referred to each other.

[0042] The above description of the disclosed embodiments enables those skilled in the art to make or use this application. Various modifications to these embodiments will be readily apparent to those skilled in the art, and the general principles defined herein may be implemented in other embodiments without departing from the spirit or scope of this application. Therefore, this application is not to be limited to the embodiments shown herein, but is to be accorded the widest scope consistent with the principles and novel features disclosed herein.

Claims

1. A method for optimizing controlled drilling trajectories in thick alluvial deposits, characterized in that, Includes the following steps: At least three key control points are set, and the three-dimensional spatial coordinates of each key control point and the tangential vector of the well trajectory at the control point are used as optimization variables. A smooth candidate drilling trajectory is generated based on the optimization variables through a preset spline interpolation function that guarantees at least the continuity of the first derivative. A three-dimensional geomechanical model coupled with formation lithology and physical parameters is established. The first hard constraint condition is that the induced stress state around the well at any point on the candidate trajectory does not exceed the rock strength envelope defined based on the preset failure criterion. The spatial uncertainty of underground obstacles and high-risk geological bodies is modeled as a three-dimensional probability density function. The boundary of the probability exclusion body is delineated according to the preset risk tolerance. The second hard constraint condition is that the candidate drilling trajectory is kept at a distance greater than the preset safety distance from the boundary. The weighted sum of the total trajectory length, the comprehensive evaluation value of torque drag, and the comprehensive evaluation value of wellbore instability risk is used as the drilling comprehensive cost objective function. The weight coefficients of each term are determined according to the formation unit attributes and local drilling direction along the trajectory. With the goal of minimizing the drilling comprehensive cost objective function, the optimal set of optimization variables is searched and determined in the solution space that simultaneously satisfies the first and second hard constraints. Based on the optimal optimization variables, the optimized drilling trajectory is generated through the spline interpolation function.

2. The method according to claim 1, characterized in that, The step of generating a smooth candidate drilling trajectory based on the optimized variables using a preset spline interpolation function that guarantees at least the continuity of the first derivative includes: The spline interpolation function is a third-order Bessel-Hermitian spline interpolation function. It uses dimensionless parameters as parameterized variables, takes the three-dimensional spatial coordinates of each key control point as the interpolation point, and takes the tangential vector of the interpolation point as the first derivative of the dimensionless parameters at the interpolation point. It interpolates the well section between two adjacent key control points, and splices all the interpolated well sections to generate a globally continuous smooth drilling trajectory.

3. The method according to claim 1, characterized in that, The establishment of a three-dimensional geomechanical model coupled with formation lithological parameters, using the condition that the well-circumferential induced stress state at any point on the candidate trajectory does not exceed the rock strength envelope defined based on a preset failure criterion, as the first hard constraint, includes: The first hard constraint condition includes three-dimensional geological meshing of the well area, assigning formation rock mechanics parameters obtained from well logging interpretation or seismic inversion to each mesh cell, including elastic modulus, Poisson's ratio, cohesion, and internal friction angle, calculating the stress-strain distribution of the three-dimensional geological body under the action of the in-situ stress field using the finite element method, calculating the induced stress field of the wellbore rock mass caused by drilling activity for any point on the candidate trajectory, calculating the critical collapse pressure using the Mohr-Coulomb failure criterion, and calculating the critical rupture pressure using the maximum tensile stress criterion, ensuring that the wellbore fluid column pressure corresponding to the stress field is completely within the stability window defined by the critical collapse pressure and rupture pressure.

4. The method according to claim 1, characterized in that, The weighted sum of the total trajectory length, torque drag comprehensive evaluation value, and wellbore instability risk comprehensive evaluation value is used as the drilling comprehensive cost objective function. The weight coefficients of each term are determined based on the formation unit attributes along the trajectory and the local drilling direction, including: The drilling comprehensive cost objective function is defined as a weighted sum of global macro-level evaluation terms and local micro-level risk evaluation terms, and its calculation formula is as follows: ,in: The objective function is the overall drilling cost. The total length of the trajectory after dimensionless transformation; The dimensionless global torque drag comprehensive evaluation value is calculated based on the vector iteration from the drill bit to the wellhead; N is the total number of micro-segments in the discrete division of the drilling trajectory; The dimensionless wellbore instability risk assessment value calculated based on the local stress field for the i-th micro-element segment; The global weight coefficient is the total length of the trajectory. This is the global weighting coefficient for the overall evaluation value of global torque drag. The local dynamic weighting coefficient for the wellbore instability risk of the i-th micro-element segment; The local dynamic weight coefficient The specific values ​​are determined by looking up a table using a pre-set engineering experience matrix, based on the lithology, permeability, and clay content of the formation unit traversed by the i-th trajectory micro-element, as well as the angle between the well trajectory at that micro-element and the formation bedding plane normal and the direction of the maximum horizontal principal stress; and for the entire well section, the weights of each item satisfy the normalization condition: .

5. The method according to claim 1, characterized in that, The process of modeling the spatial location uncertainty of underground obstacles and high-risk geological bodies as a three-dimensional probability density function, delineating the boundary of the probability exclusion body according to a preset risk tolerance, and maintaining a candidate drilling trajectory at a distance greater than a preset safety distance from the boundary as a second hard constraint condition includes: The second hard constraint condition includes modeling the spatial location uncertainty of each underground obstacle or high-risk geological body as a three-dimensional Gaussian distribution function using seismic interpretation data, setting a risk tolerance, defining the boundary of the integral region where the integral probability of the Gaussian distribution function in the spatial region reaches the tolerance value as the boundary of the probability exclusion body, and requiring that the minimum Euclidean distance from any point on the planned drilling trajectory to the boundary is greater than a preset safety distance value.

6. The method according to claim 1, characterized in that, The objective of minimizing the overall drilling cost objective function involves searching and determining an optimal set of optimization variables within the solution space that simultaneously satisfies the first and second hard constraints. This includes: The particle swarm optimization algorithm is adopted. The three-dimensional coordinates and the Cartesian component of the tangential vector of each key control point are used to form a high-dimensional optimization vector as a particle in the algorithm. In each iteration, for the candidate drilling trajectory represented by each particle, it is checked whether the candidate drilling trajectory satisfies the first and second hard constraints. If either condition is not satisfied, the cost objective function value of the particle is assigned a preset maximum value by the penalty function method. The velocity and position of the particle are updated iteratively until the maximum number of iterations is reached or the objective function value converges. The optimization variable corresponding to the globally optimal particle is output.

7. The method according to claim 1, characterized in that, The calculation of the torque drag comprehensive evaluation value includes discretizing the drilling trajectory into N micro-segments, using a soft rope model, and performing mechanical analysis segment by segment upwards from the drill bit. In each micro-segment, the weight, buoyancy, normal pressure generated by contact with the well wall, and friction force calculated based on the formation friction coefficient of the micro-segment drill string are comprehensively considered. The total tension and total torque at the wellhead are obtained through vector iteration calculation. The calculated total tension and total torque are then divided by preset reference tension and reference torque to achieve dimensionlessness. Finally, the torque drag comprehensive evaluation value is generated by combining the results according to a preset weighting method.

8. A controlled drilling trajectory optimization system for thick alluvial deposits, characterized in that, Includes the following modules: The generation module is used to set at least three key control points, and use the three-dimensional spatial coordinates of each key control point and the tangential vector of the well trajectory at the control point as optimization variables. A smooth candidate drilling trajectory is generated based on the optimization variables through a preset spline interpolation function that guarantees at least the continuity of the first derivative. A module is established to create a three-dimensional geomechanical model coupled with formation lithological and physical parameters. The first hard constraint condition is that the induced stress state around the well at any point on the candidate trajectory does not exceed the rock strength envelope defined based on a preset failure criterion. The spatial uncertainty of underground obstacles and high-risk geological bodies is modeled as a three-dimensional probability density function. The boundary of the probability exclusion body is delineated according to a preset risk tolerance. The second hard constraint condition is that the candidate drilling trajectory is kept at a distance greater than a preset safety distance from the boundary. The construction module is used to take the weighted sum of the total trajectory length, the comprehensive evaluation value of torque drag, and the comprehensive evaluation value of wellbore instability risk as the drilling comprehensive cost objective function. The weight coefficients of each item are determined according to the formation unit attributes and local drilling direction along the trajectory. With the goal of minimizing the drilling comprehensive cost objective function, the module searches and determines the optimal set of optimization variables in the solution space that simultaneously satisfies the first and second hard constraints. Based on the optimal optimization variables, the module generates an optimized drilling trajectory through the spline interpolation function.

9. The system according to claim 8, characterized in that, The step of generating a smooth candidate drilling trajectory based on the optimized variables using a preset spline interpolation function that guarantees at least the continuity of the first derivative includes: The spline interpolation function is a third-order Bessel-Hermitian spline interpolation function. It uses dimensionless parameters as parameterized variables, takes the three-dimensional spatial coordinates of each key control point as the interpolation point, and takes the tangential vector of the interpolation point as the first derivative of the dimensionless parameters at the interpolation point. It interpolates the well section between two adjacent key control points, and splices all the interpolated well sections to generate a globally continuous smooth drilling trajectory.

10. A computer-readable storage medium storing a computer program thereon, characterized in that, The computer program, when executed by a processor, implements the method as described in any one of claims 1-7.