A method for analytically calculating the stress around a wellbore of irregular shape

Abstract

The present application belongs to the technical field of wellbore stability analysis, and particularly relates to a wellbore stress analytical calculation method for irregular wellbore shape. The present application is based on radial projection iteration and fast Fourier transform, and combines the maximum change of the conformal mapping coefficient in the adjacent two iterations, the relative boundary fitting error and the high-order coefficient tail proportion factor to adaptively determine the truncation order of the conformal mapping function and the optimal conformal mapping coefficient, and use the optimal conformal mapping coefficient for the construction of the perturbation complex potential function and the calculation of the wellbore stress component. The present application can avoid invalid high-order term calculation, thereby improving the calculation efficiency and stability, and has both the advantages of fast calculation of analytical solution and the ability to accurately process real complex wellbore boundaries.

Description

Technical Field

[0001] This invention belongs to the field of wellbore stability analysis technology, specifically relating to an analytical calculation method for wellbore stress in irregular wellbore shapes. Background Technology

[0002] During oil and gas development, the actual wellbore shape changes due to factors such as geostress and drilling, often deviating from the ideal circular shape and exhibiting an irregular form. Currently, common methods for calculating wellbore stress are mainly divided into two categories: analytical methods and numerical simulation methods. Analytical methods mostly simplify the wellbore shape to an ideal circle or an approximate ellipse, failing to accurately calculate the wellbore stress of truly irregular wellbores. While numerical simulation methods can consider the actual wellbore shape, the calculation process is time-consuming and labor-intensive. Therefore, in the preliminary analysis of wellbore geometry parameterization, there is an urgent need for an analytical calculation method that can quickly and accurately calculate the wellbore stress of irregular wellbore shapes. Summary of the Invention

[0003] To overcome the shortcomings of existing technologies, this invention provides an analytical method for calculating wellbore stress in irregular wellbore shapes. This invention avoids the calculation of invalid higher-order terms, thereby improving computational efficiency and stability. It possesses both the advantages of rapid analytical solution calculation and the ability to accurately handle real, complex wellbore boundaries.

[0004] The present invention discloses an analytical calculation method for wellbore stress in irregular wellbore shapes. Based on radial projection iteration and fast Fourier transform, it combines the maximum change of conformal mapping coefficients in two adjacent iterations, the relative boundary fitting error, and the tail proportion factor of higher-order coefficients. Under the premise of satisfying the boundary fitting accuracy, it adaptively determines the truncation order of the conformal mapping function and the optimal conformal mapping coefficients, and uses the optimal conformal mapping coefficients for constructing the perturbation complex potential function and calculating the wellbore stress components.

[0005] To achieve the above objectives, the technical solution of the present invention includes the following core steps:

[0006] (1) The optimal conformal mapping coefficients are obtained by solving the conformal mapping and adaptive truncation order for the irregular wellbore boundary; specifically:

[0007] (1-1) Using conformal mapping function, construct discrete sampling points of unit circle boundary based on well wall boundary.

[0008] Mapping the infinitely large region outside the irregular well wall in the physical plane to the mathematical plane ( In a plane, the infinitely large region outside the unit circle corresponds to the boundary of the unit circle in the mathematical plane. For ease of numerical solution, the total circumference angle 2π of the unit circle in the mathematical plane is divided into equal parts. Share, received The discrete sampling angles of the unit circle boundary, and the corresponding discrete sampling points of the unit circle boundary are:

[0009] ;

[0010] ;

[0011] in, The total number of discrete sampling points at the boundary of the unit circle, dimensionless; The first unit circle boundary The number of each discrete sampling point is dimensionless. The first unit circle boundary Each discrete sampling angle, in rad; The first unit circle boundary A discrete sampling point, dimensionless, located in the mathematical plane ( On the unit circle boundary of the plane, it serves as the input to the conformal mapping function.

[0012] At the same time, in the Discrete sampling angle Above, calculate the polar radius of the corresponding true boundary point according to the wellbore boundary geometric equation. To obtain the physical plane ( The true boundary points of the well wall boundary in the plane:

[0013] ;

[0014] in, For the first Discrete sampling angle The corresponding polar radius of the true boundary point, m; Let m be the true boundary point, located on the wellbore boundary, and serve as the target that the conformal mapping function needs to approximate.

[0015] Conformal mapping functions are typically expressed in the form of Laurent series expansions of complex functions:

[0016] ;

[0017] in, For physical plane complex coordinates, , Let m be the physical plane coordinate. The imaginary unit; It is a mathematical plane, dimensionless; Let m be the conformal mapping function. For the first Conformal mapping coefficient, m; The truncation order is dimensionless; The order number of the conformal mapping coefficients is dimensionless.

[0018] ;

[0019] in, and They are respectively The real and imaginary parts of , m.

[0020] (1-2) Set the initial truncation order and initialize the conformal mapping coefficients.

[0021] Using the polar coordinate radius of the actual boundary point on the wellbore boundary relative to the origin as input, initialize the zeroth-order conformal mapping coefficient. Discrete sampling angle The polar radius of the corresponding true boundary point is denoted as Unit: m, initial value: The average polar radius of the true boundary points:

[0022] ;

[0023] The remaining higher-order conformal mapping coefficients are initialized as follows:

[0024] ;

[0025] Where: the superscript (1) indicates the first iteration; m represents the 0th order conformal mapping coefficient after the first iteration; After the first iteration Conformal mapping coefficient, m; The order numbering of the conformal mapping coefficients is dimensionless; The total number of discrete sampling points at the boundary of the unit circle, dimensionless; The first unit circle boundary Each discrete sampling point is numbered, dimensionless; For the first Discrete sampling angle The polar radius of the corresponding true boundary point, m.

[0026] (1-3) Substitute the discrete sampling points of the unit circle boundary into the conformal mapping function to obtain the corresponding current mapping boundary points.

[0027] Discrete sampling points are taken at the boundary of the unit circle: Substitute the current conformal mapping function, and based on the current mapping coefficients... Calculate the physical plane corresponding to each discrete sampling point on the boundary of the unit circle. The current mapping boundary point in the next iteration :

[0028] ;

[0029] in: The number of iterations is dimensionless. For the first During the nth iteration, the boundary of the unit circle... m are the current mapped boundary points corresponding to discrete sampling points; For the first During the nth iteration, the boundary of the unit circle... The conformal mapping function corresponding to each discrete sampling point, m; Polar angle, rad; The order numbering of the conformal mapping coefficients is dimensionless; The truncation order is dimensionless.

[0030] (1-4) Radial projection of the current mapped boundary point matches the target boundary point with the same argument.

[0031] For the current mapped boundary point Find its argument:

[0032] ;

[0033] in, Let be the argument of the current mapping boundary point in rad during the j-th iteration.

[0034] And find the target boundary points with the same argument on the well wall boundary:

[0035] ;

[0036] in, For the first During the nth iteration, the boundary of the unit circle... m are the current mapped boundary points corresponding to discrete sampling points; For the first The argument of the current mapped boundary point at the next iteration, in rad; For the target boundary point at the argument angle The polar radius function in the direction, m; Let m be the target boundary point with the same argument as the current mapped boundary point.

[0037] (1-5) Construct an auxiliary function and obtain the conformal mapping coefficients for the new round of iterations through fast Fourier transform.

[0038] To facilitate the extraction of updated conformal mapping coefficients from discrete sampling points, an auxiliary function is constructed:

[0039] ;

[0040] And because Therefore

[0041] ;

[0042] right Perform a Fast Fourier Transform (FFT) to obtain the conformal mapping coefficients for the next iteration.

[0043] ;

[0044] in, For the first During the nth iteration The auxiliary function values ​​corresponding to each discrete sampling point, m; For the first After the first iteration update Conformal mapping coefficient, m; The first unit circle boundary A discrete sampling point, dimensionless; The first unit circle boundary Each discrete sampling angle, in rad; Let m be the target boundary point with the same argument as the currently mapped boundary point; The order numbering of the conformal mapping coefficients is dimensionless; The total number of discrete sampling points at the boundary of the unit circle, dimensionless; The first unit circle boundary Each discrete sampling point is numbered, dimensionless; Let m be the function value of the auxiliary function at the unit circle boundary angle φ; Let m be the complex coordinates of the target boundary point with the same argument as the current mapped boundary point in the physical plane. Let be the complex variable sampling point corresponding to angle φ on the boundary of the unit circle, which is dimensionless.

[0045] (1-6) Determine whether the conformal mapping coefficients have converged iteratively; if they have converged, continue to step (1-7); otherwise, return to step (1-2) and iterate again.

[0046] The maximum change in conformal mapping coefficients between two consecutive iterations is:

[0047] ;

[0048] in, Let m be the conformal mapping coefficient of the kth order after the (j+1)th iteration update; Let m be the conformal mapping coefficient of the kth order after the j-th iteration update; Let m be the convergence error of the coefficients in the j-th iteration; To preset the convergence tolerance, take ;like ≤ If the iteration converges, then the iteration will be successful.

[0049] In the process of solving the conformal mapping coefficients, to balance the fitting accuracy, numerical stability, and computational efficiency of the wellbore boundary, this invention further introduces an adaptive truncation order determination strategy. Specifically, for any wellbore geometry, the boundary is first discretized into discrete sampling points, and then candidate truncation orders are used... The conformal mapping function is solved step by step, and the search is terminated based on a comprehensive assessment of the relative boundary fitting error index and the attenuation index of higher-order coefficients of the current mapping result. This automatically determines the optimal truncation order that meets the accuracy requirements. That is:

[0050] (1-7) Based on the convergence of the iteration in step (1-6), determine whether the adaptive truncation order strategy is satisfied; if satisfied, output the optimal conformal mapping coefficient; otherwise, return to step (1-3).

[0051] The adaptive truncation order strategy is as follows: <0.5% and <0.05;

[0052] The formula for calculating the relative boundary fitting error is:

[0053] ;

[0054] in, Let m be the polar radius of the physical plane obtained after conformal mapping of the nth discrete sampling point when the truncation order is N; Let m be the polar radius of the target boundary point with the same argument when the truncation order is N; The relative boundary fitting error of the nth discrete sampling point is given by the truncation order N, and is dimensionless; N is the truncation order, which is also dimensionless.

[0055] The tail proportion factor of higher-order coefficients under the current candidate truncation order is calculated using the following formula:

[0056] ;

[0057] This is a dimensionless tail proportion factor for higher-order coefficients when the truncation order is N. The statistical length of the tail coefficient is dimensionless. The starting order for tail statistics; For N, the k-th order conformal mapping coefficients are dimensionless. The corresponding rigid body translation term is not included in the tail proportion factor statistics. The smaller the tail proportion factor, the more likely the higher-order mapping terms at the current truncation order have reached decay stability, and further increasing the truncation order will bring limited boundary improvement.

[0058] The optimal conformal mapping coefficients obtained above are used for constructing the perturbation complex potential function and calculating the wellbore stress components.

[0059] This invention considers the anisotropic characteristics of formation materials and far-field stress, and uses the obtained optimal conformal mapping coefficients to construct the perturbation complex potential function in the mathematical plane variable. The perturbation complex potential functions described in this invention are expressed as follows:

[0060] ;

[0061] ;

[0062] in, Let be the complex potential function of the first perturbation, MPa·m; A1, A2, A3, A4 are the coefficients of the first perturbation complex potential function, MPa; A5, A6, A7, A8 are the coefficients of the second perturbation complex potential function, MPa. Let m be the real part of the optimal conformal mapping coefficient; m is the imaginary part of the optimal conformal mapping coefficients; m is the real part of the 0th-order conformal mapping coefficient; m is the imaginary part of the 0th-order conformal mapping coefficient; It is a mathematical plane, dimensionless.

[0063] After obtaining the perturbation complex potential function, the derivative of the perturbation complex potential function with respect to the generalized complex variable is obtained according to the strict chain rule. Then, the stress components in the Cartesian coordinate system on the wellbore boundary are calculated by combining the far-field stress and the internal pressure of the wellbore. The stress components are then transformed to the wellbore coordinate system to obtain the wellbore stress components, namely radial stress, circumferential stress and tangential shear stress.

[0064] The derivative of the perturbation complex potential function with respect to the generalized complex variable obtained in this invention is:

[0065] ;

[0066] ;

[0067] in, For the first perturbation complex potential function with respect to the generalized complex variable The derivative, MPa; For the second perturbation complex potential function with respect to the generalized complex variable The derivative of μ, MPa; μ1 and μ2 are characteristic roots, dimensionless complex numbers; Let m be the generalized complex variable corresponding to the first perturbation complex potential function; Let m be the generalized complex variable corresponding to the second perturbation complex potential function.

[0068] Stress components in Cartesian coordinates (with tensile stress as positive):

[0069] ;

[0070] ;

[0071] ;

[0072] in, The normal stress in the x-direction at the wellbore boundary is expressed in MPa. The normal stress in the y-direction at the wellbore boundary is given in MPa. The shear stress at the wellbore boundary is expressed in MPa.

[0073] The stress components in the wellbore coordinate system are obtained through coordinate transformation:

[0074] ;

[0075] ;

[0076] ;

[0077] Where δ is the directed angle between the unit outward normal vector and the positive x-axis of the global rectangular coordinate system, in rad; The radial stress in the wellbore coordinate system is expressed in MPa. Let be the circumferential stress in the wellbore coordinate system, in MPa; Let be the tangential shear stress in the wellbore coordinate system, in MPa.

[0078] When there is a fluid column pressure P inside the borehole w At that time, the effective stress superposition method based on linear elasticity theory is used to correct the stress field, and finally output the true stress distribution law around the well with arbitrary irregular shape.

[0079] The above methods are commonly used in the prior art and will not be elaborated further here.

[0080] Compared with the prior art, the present invention has the following beneficial effects:

[0081] (1) Introduce conformal mapping coefficients to parameterize the real irregular wellbore boundary into an analytical mapping function in the complex plane. Through conformal mapping coefficients, the outer domain of the unit circle in the mathematical plane can be mapped to the outer domain of the irregular wellbore in the physical plane, thereby transforming the complex boundary conditions into standard circular domain boundary conditions, which facilitates the construction of the subsequent perturbation complex potential function and the calculation of the wellbore stress components.

[0082] (2) Conformal mapping coefficients of different orders can characterize the geometric features of the wellbore boundary at different scales. In the process of solving the conformal mapping coefficients, this invention introduces an adaptive truncation order determination strategy. The relative boundary fitting error and the proportion of the tail of the higher-order coefficients are used to jointly determine whether to continue to increase the truncation order. The truncation order that matches its geometric complexity is automatically selected, avoiding the use of fixed higher-order expansion for all working conditions. This can reduce the calculation of invalid higher-order terms while meeting the boundary fitting accuracy, thereby improving computational efficiency and numerical stability.

[0083] (3) This invention utilizes conformal mapping functions and radial projection iterations to overcome the limitation that traditional analytical solutions are only applicable to circular or elliptical shapes. It can accurately characterize rough and extremely irregular wellbore walls, such as those detected by well logging dual-diameter instruments. This invention is not only applicable to conventional homogeneous formations, but can also accurately analyze irregular wellbore stress components in highly anisotropic formations. Attached Figure Description

[0084] Figure 1 A comparison diagram of the circumferential stress of the present invention and the Kirsch solution in isotropic formations with circular wellbores;

[0085] Figure 2 A comparison diagram of the circumferential stress in this invention and the Lekhnitskii theoretical solution in anisotropic formations with a circular wellbore;

[0086] Figure 3 A comparison diagram of the circumferential stress of the present invention and the analytical solution of the Savin elliptical borehole in isotropic formations with elliptical wellbores;

[0087] Figure 4 A comparison diagram of the circumferential stress in this invention and the Lekhnitskii theoretical solution in anisotropic formations with an elliptical wellbore;

[0088] Figure 5 A comparison diagram of the circumferential stress in the present invention and the finite element method in an isotropic formation of an explicit collapse wellbore;

[0089] Figure 6 A comparison diagram of circumferential stress in anisotropic formations with explicit collapse wellbore, using the methods of this invention and finite element analysis.

[0090] Figure 7 This is a diagram showing the circumferential stress distribution in an isotropic formation with an irregular, low-frequency, rough wellbore.

[0091] Figure 8 This is a diagram showing the circumferential stress distribution in anisotropic formations with irregular low-frequency rough wellbore. Detailed Implementation

[0092] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to specific basic parameters and four typical wellbore shapes.

[0093] In all the following embodiments, the far-field stress environment is defined as: the maximum horizontal principal stress σ H =50 MPa, minimum horizontal principal stress σ h =30 MPa.

[0094] If the directions of the maximum and minimum horizontal principal stresses are aligned with the global y-axis and x-axis, respectively, then the corresponding far-field stress components are:

[0095] =30MPa, =50MPa, =0Mpa.

[0096] Meanwhile, to verify the universality of the material properties, two types of formation media were examined:

[0097] Isotropic formation: elastic modulus E = 30 GPa, Poisson's ratio v = 0.25; corresponding compliance matrix coefficients are: a 11 = a 22 =1 / 30, a 12 =-0.25 / 30, a 66 =2(1+0.25) / 30,a 16 = a 26 =0;

[0098] Anisotropic strata: elastic modulus E1 = 1 GPa, E2 = 12 GPa, Poisson's ratio v 12 =0.03, shear modulus G 12 =0.5GPa; the corresponding compliance matrix coefficients are: a 11 =1, a 22 =1 / 12, a 12 =-0.03, a 66 =1 / 0.5, a 16 = a 26 =0.

[0099] Due to the circumferential stress in wellbore stability analysis The circumferential stress distribution is usually the main indicator for evaluating stress concentration, wellbore failure initiation location, and comparison with classical analytical solutions. Furthermore, the verification of Kirsch solutions for circular wells, Savin solutions for elliptical wells, and Lekhnitskii theoretical solutions for anisotropic media typically uses circumferential stress distribution as the primary comparison object. Therefore, the accompanying drawings of the invention focus on illustrating the circumferential stress distribution. The distribution of the inequality and its comparison with the theoretical solution.

[0100] Example 1: Calculation of an ideal circular wellbore

[0101] This embodiment is used to verify the accuracy and theoretical degradation capability of the present invention when dealing with the most basic undeformed circular wellbore.

[0102] Wellbore radius R = 1m, effective hydraulic column pressure P exists inside the wellbore. w =15 MPa.

[0103] Its boundary discrete expression is:

[0104] ;

[0105] (1) Perform conformal mapping and adaptive truncation order calculation on the well wall boundary to obtain the optimal conformal mapping coefficient;

[0106] (1-1) Using conformal mapping function, construct discrete sampling points of unit circle boundary based on well wall boundary.

[0107] Divide the circumference of the mathematical plane, 2π, into M = 4096 equal parts to construct discrete sampling points. Through the above processing, 4096 discrete sampling points are obtained on the boundary of the unit circle. These discrete sampling points serve as the input for subsequent conformal mapping solutions. Since R = 1 / m, the nth discrete sampling angle θ... n The corresponding polar radius r of the true boundary point n =1m, obtain the true boundary points of the circular well wall boundary: .

[0108] Conformal mapping functions are typically expressed in the form of Laurent series expansions of complex functions:

[0109] ;

[0110] .

[0111] (1-2) Initialize the conformal mapping coefficients and set the initial truncation order N0=2;

[0112] With the initial cutoff order N0=2, the adaptive mapping of the circular hole boundary satisfies the requirement at the minimum candidate cutoff order, and convergence can be achieved without further iteration.

[0113] ;

[0114] Since the radius of a circular wellbore is the same in all directions, the average polar radius of the actual boundary points is taken as the initial zeroth-order mapping coefficient m0, and r n Substituting =1m:

[0115] ;

[0116] have to:

[0117] m0=1.0000000000+0.0000000000i, m1≈0, m2≈0;

[0118] After this initialization step, the coefficients in the conformal mapping function have been assigned initial values, forming an initial conformal mapping function that can be used for the next calculation. The boundary of this initial conformal mapping function is now the unit circle boundary.

[0119] (1-3) Substitute the discrete sampling points of the unit circle boundary into the conformal mapping function to obtain the corresponding current mapped boundary points;

[0120] Discrete sampling points are taken at the boundary of the unit circle: Substituting into the conformal mapping function, we obtain the current mapping boundary point in the j-th iteration:

[0121] ;

[0122] In this embodiment, since the initial conformal mapping function is in circular mapping form, therefore: .

[0123] The current mapped boundary point coincides with the target boundary point of the wellbore boundary. The radius of the current mapped boundary point is r. map,n =1m, the radius of the target boundary point on the well wall boundary is: r target,n =1m.

[0124] (1-4) Radial projection of the current mapped boundary point matches the target boundary point with the same argument.

[0125] For the current mapped boundary point z n (j) Find its argument:

[0126] ;

[0127] Let be the argument of the current mapping boundary point in rad during the j-th iteration.

[0128] And find the target boundary points with the same argument on the well wall boundary:

[0129] ;

[0130] Since the radius of a circular wellbore is 1m in any direction, the current mapped boundary point is completely consistent with the target boundary point.

[0131] (1-5) Construct an auxiliary function and obtain the conformal mapping coefficients for the new round of iterations through fast Fourier transform;

[0132] Construct auxiliary function values ​​based on the target boundary points. And perform a fast Fourier transform on the M=4096 auxiliary function values ​​to update the conformal mapping coefficients.

[0133] The auxiliary function is:

[0134] ;

[0135] The conformal mapping coefficient is:

[0136] ;

[0137] After iteration, only the constant term m0=1 of the conformal mapping coefficients retains a valid value, while the remaining higher-order coefficients m k (k≥1) strictly approaches 0.

[0138] (1-6) Determine whether the conformal mapping coefficients have converged iteratively; if they have converged, continue to step (1-7); otherwise, return to step (1-2) and iterate again.

[0139] Calculate the maximum change in the mapping coefficients between two consecutive iterations:

[0140] ;

[0141] If Δm (j) If ≤ ε, then the iteration converges; Δm (j) Let be the coefficient convergence error of the j-th iteration, and ε be the preset convergence tolerance, which is taken as... .

[0142] (1-7) Based on the convergence of the iteration in step (1-6), determine whether the adaptive truncation order strategy is satisfied; if satisfied, output the optimal conformal mapping coefficient; otherwise, return to step (1-3).

[0143] The adaptive truncation order strategy is as follows: <0.5% and <0.05;

[0144] The relative boundary fitting error of the nth discrete sampling point on the unit circle boundary is calculated using the following formula:

[0145] ;

[0146] The tail proportion factor of higher-order coefficients under the current truncation order is calculated using the following formula:

[0147] ;

[0148] in, is the tail proportion factor of the higher-order coefficients under the current candidate truncation order, dimensionless; L is the statistical length of the tail coefficients, dimensionless; ks is the starting order of the tail statistics; m1 corresponds to the rigid body translation term and does not participate in the statistics of the tail proportion factor.

[0149] The circular wellbore meets the convergence requirement at the minimum candidate cutoff order. The cutoff order is adaptively determined as N=2, and the optimal conformal mapping coefficient is m0=1.0000000000+0.0000000000i.

[0150] (2) The optimal conformal mapping coefficients obtained above are used for constructing the perturbation complex potential function and calculating the wellbore stress components. The specific process is as follows:

[0151] (2-1) Construct a unified characteristic equation based on the compliance matrix coefficients of the formation material and solve for the characteristic roots

[0152] The unified characteristic equation is as follows:

[0153] ;

[0154] Among them, a 11 a 12 a 16 a 22 a 26 a 66 These are the compliance matrix coefficients, in GPa. −1 This equation applies to both isotropic and anisotropic strata.

[0155] Isotropic strata: elastic modulus E = 30 GPa, Poisson's ratio v = 0.25.

[0156] The corresponding flexibility matrix coefficients are: a 11 = a 22 =1 / 30, a 12 =-0.25 / 30, a 66 =2(1+0.25) / 30,a 16 = a 26 =0.

[0157] Anisotropic strata: elastic modulus E1 = 1 GPa, E2 = 12 GPa, Poisson's ratio v 12 =0.03, shear modulus G 12 =0.5GPa.

[0158] The corresponding flexibility matrix coefficients are: a 11 =1, a 22 =1 / 12, a 12 =-0.03, a 66 =1 / 0.5, a 16 = a26 =0.

[0159] Substituting the coefficients of the above compliance matrix into the unified characteristic equation, solving the fourth-degree polynomial, and selecting two eigenvalues ​​with imaginary parts greater than 0, we obtain:

[0160] .

[0161] (2-2) Define generalized complex variables

[0162] ;

[0163] ;

[0164] in, Let m be the generalized complex variable corresponding to the first perturbation complex potential function; Let m be the generalized complex variable corresponding to the second perturbation complex potential function; let m be the physical plane rectangular coordinates of x and y. μ1 and μ2 are the eigenvalues ​​obtained in step (2-1), which are dimensionless.

[0165] (2-3) Calculate the far-field stress coefficient and the coefficient of the perturbation complex potential function.

[0166] Far-field stress components are =30MPa, =50MPa, =0MPa and the obtained eigenvalues ​​are substituted into the formula for calculating the far-field stress coefficient.

[0167] ;

[0168] Calculate the far-field stress coefficients B, B', C':

[0169] .

[0170] Substituting the far-field stress coefficients B, B', and C' into the equations, we calculate the intermediate complex parameters K1, K2, … K8. Further calculation of the perturbation complex potential function coefficients yields:

[0171] ;

[0172] .

[0173] (2-4) Construct the perturbation complex potential function based on the optimal conformal mapping coefficients and the perturbation complex potential function coefficients;

[0174] ;

[0175] ;

[0176] (2-5) After obtaining the perturbation complex potential function, the derivative of the perturbation complex potential function with respect to the generalized complex variable is further obtained according to the strict chain rule. The derivative with θ as the intermediate variable is then transformed into the perturbation complex potential function with respect to the generalized complex variable. , The derivative of .

[0177] At the nth discrete sampling point on the well wall boundary, based on the boundary coordinates x n y n And the obtained eigenvalues, construct generalized complex variables:

[0178] ;

[0179] ;

[0180] For circular wellbores:

[0181] ;

[0182] ;

[0183] At the same time, differentiate the boundary coordinates with respect to θ:

[0184] ;

[0185] ;

[0186] again ,

[0187] ;

[0188] have to:

[0189] ;

[0190] ;

[0191] The obtained perturbation complex potential function coefficients A1~A8 and the obtained conformal mapping coefficients m k (Separate the real part) , and the virtual part , Substituting the expression for the perturbation complex potential function, we first obtain:

[0192] ;

[0193] ;

[0194] Then, using the strict chain rule, we obtain the derivative of the complex potential function:

[0195] ;

[0196] ;

[0197] To illustrate this chain-reaction differentiation process, four representative discrete sampling points on the wellbore at θ = 0°, 30°, 60°, and 90° are selected for illustration. Anisotropic formations are used as an example:

[0198] μ1 =-0.0000000000+1.3769709346i, μ2 =0.0000000000+0.2096450458i;

[0199] The generalized complex variable representing the discrete sampling points is:

[0200] .

[0201] Further differentiation yields:

[0202] .

[0203] The derivative of the complex potential function can be obtained through , Divide by the table respectively , The remaining discrete sampling points are then calculated point by point in the same manner.

[0204] (2-6) After obtaining the derivative of the perturbation complex potential function with respect to the generalized complex variable, the stress components in the Cartesian coordinate system on the well wall boundary are calculated by combining the far-field stress and the pressure inside the well. The stress components in the well wall coordinate system are obtained by coordinate transformation.

[0205] The stress components in Cartesian coordinates at the wellbore boundary are:

[0206] ;

[0207] ;

[0208] ;

[0209] in, The normal stress in the x-direction at the wellbore boundary is expressed in MPa. The normal stress in the y-direction at the wellbore boundary is given in MPa. The shear stress at the wellbore boundary is expressed in MPa.

[0210] The stress components in the wellbore coordinate system are obtained through coordinate transformation:

[0211] ;

[0212] ;

[0213] ;

[0214] Wherein, the direction angle δ is defined as the directed angle between the unit outward normal vector and the positive x-axis of the global rectangular coordinate system, in rad; The radial stress in the wellbore coordinate system is expressed in MPa. Let be the circumferential stress in the wellbore coordinate system, in MPa; Let be the tangential shear stress in the wellbore coordinate system, in MPa.

[0215] In this embodiment, there is an effective liquid column pressure inside the wellbore. By effectively correcting the stress on the wellbore wall, the stress distribution around the wellbore in isotropic and anisotropic formations can be obtained.

[0216] The wellbore circumferential stress distribution calculated using this method is compared with the classical analytical solution. In isotropic formations, the results completely coincide with the Kirsch exact solution, as shown below. Figure 1 In anisotropic strata, the results are basically consistent with the exact solution of Lekhnitskii theory, such as... Figure 2 .

[0217] Example 2: Calculation of Elliptical Wellbore

[0218] The major semi-axis of the wellbore is a = 3 m, the minor semi-axis is b = 1 m, and P w =0 MPa, the boundary discrete expression is:

[0219] ;

[0220] (1) Perform conformal mapping and adaptive truncation order calculation on the well wall boundary to obtain the optimal conformal mapping coefficient.

[0221] (1-1) Using conformal mapping function, construct discrete sampling points of unit circle boundary based on well wall boundary.

[0222] Divide the circumference of the mathematical plane, 2π, into M=4096 equal parts to construct discrete sampling points for subsequent solution of conformal mapping coefficients.

[0223] At each discrete sampling angle θ n Above, calculate the actual boundary points of the wellbore boundary according to the elliptical wellbore boundary equation:

[0224] ;

[0225] The corresponding complex coordinates are:

[0226] ;

[0227] (1-2) Initialize the conformal mapping coefficients and set the initial truncation order N0=2.

[0228] The conformal mapping function is in Laurent series form:

[0229] ;

[0230] First, initialize the zeroth-order mapping coefficients m0 based on the average polar radius of the true boundary points. (1) The coefficient is approximately 1.2163, with the remaining coefficients initialized to 0. For elliptical wellbores, the radii of the boundary vary along different directions, and the zeroth-order term alone cannot accurately describe the elliptical boundary. Therefore, higher-order mapping coefficients need to be updated through subsequent radial projection iterations and FFT.

[0231] (1-3) Substitute the discrete sampling points of the unit circle boundary into the conformal mapping function to obtain the corresponding current mapped boundary points;

[0232] Discrete sampling points are taken at the boundary of the unit circle: =𝑒 𝑖𝜃n Substituting into the current conformal mapping function, we obtain the current mapping boundary point in the j-th iteration:

[0233] ;

[0234] The current mapped boundary points usually cannot completely coincide with the elliptical well wall boundary in one go, so radial projection correction is required.

[0235] (1-4) Radial projection of the current mapped boundary point matches the target boundary point with the same argument.

[0236] For the current mapped boundary point z n (j) Find its argument:

[0237] ;

[0238] And find the target boundary points with the same argument on the well wall boundary:

[0239] ;

[0240] For elliptical boundaries:

[0241] ;

[0242] When the argument is When, the polar radius function of the target boundary point in this direction can be written as:

[0243] ;

[0244] (1-5) Construct an auxiliary function and obtain the conformal mapping coefficients for the new round of iterations through fast Fourier transform;

[0245] Construct auxiliary function values ​​based on the target boundary points. And perform a fast Fourier transform on the M=4096 auxiliary function values ​​to update the conformal mapping coefficients.

[0246] The auxiliary function is:

[0247] ;

[0248] The conformal mapping coefficient is:

[0249] ;

[0250] In this embodiment, the elliptical boundary has obvious second-order geometric features compared to the circular boundary. Therefore, after the FFT update, except for the zero-order scale term, the second-order mapping coefficient m2 is the main higher-order term.

[0251] The obtained conformal mapping coefficients for the first few orders are: m0 = 2.0000018897 + 0.0000000000i, m1 ≈ 0, m2 = 0.9999860994 - 0.0000000000i, m3 ≈ 0, m4 = 0.0000056552 - 0.0000000000i, m5 ≈ 0, ..., m 21 ≈0,m 22 =0.0000000268+0.0000000000i;

[0252] All other higher-order terms show significant decay. Among them:

[0253]

[0254] (1-6) Determine whether the conformal mapping coefficients have converged iteratively; if they have converged, continue to step (1-7); otherwise, return to step (1-2) and iterate again.

[0255] Calculate the maximum change in the mapping coefficients between two consecutive iterations:

[0256] ;

[0257] like ≤ If the iteration converges, then the iteration will be successful. Let be the convergence error of the coefficients in the j-th iteration. To preset the convergence tolerance, take .

[0258] (1-7) Based on the convergence of the iteration in step (1-6), determine whether the adaptive truncation order strategy is satisfied; if satisfied, output the optimal conformal mapping coefficient; otherwise, return to step (1-3).

[0259] The adaptive truncation order strategy is as follows: <0.5% and <0.05;

[0260] The relative boundary fitting error of the nth discrete sampling point on the unit circle boundary is calculated using the following formula:

[0261] .

[0262] The tail proportion factor of higher-order coefficients under the current truncation order is calculated using the following formula:

[0263] .

[0264] In this embodiment, when the truncation order is N=22, the boundary fitting error meets the preset requirements, and the higher-order mapping coefficients have significantly decreased. The iteration count is 98, and the optimal conformal mapping coefficients are: m0=2.0000018897+0.0000000000i, m1≈0, m2= 0.9999860994-0.0000000000i, m3≈0, m4=0.0000056552-0.0000000000i, m5≈0, ..., m 21 ≈0,m 22 =0.0000000268+0.0000000000i.

[0265] (2) The optimal conformal mapping coefficients obtained above are used for the construction of the perturbation complex potential function and the calculation of the well perimeter stress components. The specific process is the same as in Example 1, and will not be repeated here.

[0266] To illustrate this chain-reaction differentiation process, four representative wellbore points with θ = 0°, 30°, 60°, and 90° are selected for illustration. Anisotropic formations are used as an example:

[0267] μ1 =-0.0000000000+1.3769709346i, μ2 =0.0000000000+0.2096450458i;

[0268] The generalized complex variable representing the discrete sampling points is:

[0269] .

[0270] Further differentiation yields:

[0271] .

[0272] The derivative of the complex potential function can be obtained through , Divide by the table respectively , The remaining discrete sampling points are then calculated point by point in the same manner.

[0273] In this embodiment, there is no effective liquid column pressure inside the wellbore, so there is no need to perform effective stress correction on the wellbore stress. The wellbore stress distribution in isotropic and anisotropic formations can be obtained for elliptical wellbores.

[0274] For elliptical wellbores, the method can automatically determine the optimal order as 22 through adaptive conformal mapping, with m0 and m2 being the dominant terms in the mapping coefficients, accurately reflecting the geometric characteristics of the elliptical boundary. In isotropic formations, the results are highly consistent with the exact solution for the Savin ellipse, such as... Figure 3 In anisotropic strata, the results are in high agreement with the Lekhnitskii theoretical solution, such as... Figure 4 .

[0275] Example 3: Calculation of complex wellbore with explicit collapse features

[0276] This embodiment is used to verify the adaptability of the present invention to complex geometries with locally extremely high curvature (such as the tip of a wellbore collapse). An explicit nonlinear boundary model is introduced, and the extreme radius function is set as... The depth control factor is 0.5, and the sharpness factor is 8. An effective fluid column pressure P exists inside the wellbore. w =15 MPa, and its boundary discrete expression is:

[0277] ;

[0278] (1) Perform conformal mapping and adaptive truncation order calculation on the well wall boundary to obtain the optimal conformal mapping coefficient;

[0279] (1-1) Based on the well wall boundary, construct discrete sampling points of the unit circle boundary to obtain the conformal mapping function.

[0280] Divide the circumference of the mathematical plane, 2π, into M=4096 equal parts to construct discrete sampling points for subsequent solution of conformal mapping coefficients.

[0281] At each discrete sampling angle θ n Above, calculate the actual boundary points of the wellbore according to the above wellbore boundary equation:

[0282] ;

[0283] The corresponding complex coordinates are:

[0284] ;

[0285] (1-2) Initialize the conformal mapping coefficients and set the initial truncation order N0=2.

[0286] The initial truncation order N0=2, and the conformal mapping function adopts the Laurent series form:

[0287] ;

[0288] First, initialize the zeroth-order mapping coefficients m0 based on the average polar radius of the true boundary points. (1) The value is approximately 1.1871, and the remaining higher-order coefficients are initialized to 0. For the aforementioned boundary, the radius varies along different directions, and the zeroth-order term alone cannot accurately describe this boundary. Therefore, subsequent radial projection iterations and FFT are needed to update the higher-order mapping coefficients.

[0289] (1-3) Substitute the discrete sampling points of the unit circle boundary into the conformal mapping function to obtain the corresponding current mapped boundary points;

[0290] Discrete sampling points are taken at the boundary of the unit circle: Substituting into the current conformal mapping function, we obtain the current mapping boundary point in the j-th iteration:

[0291] ;

[0292] The current mapped boundary points usually cannot completely coincide with the real boundary points in one go, so radial projection correction is still required.

[0293] (1-4) Radial projection of the current mapped boundary point matches the target boundary point with the same argument.

[0294] For the current mapped boundary point z n (j) Find its argument:

[0295] ;

[0296] And find the target boundary points with the same argument on the well wall boundary:

[0297] ;

[0298] For the above target boundary points:

[0299] ;

[0300] When the argument is When the polar radius of the target boundary point in this direction is:

[0301] ;

[0302] The corresponding target boundary points are:

[0303] ;

[0304] (1-5) Construct an auxiliary function and obtain the conformal mapping coefficients for the new round of iterations through fast Fourier transform;

[0305] Construct auxiliary function values ​​based on the target boundary points. And perform a fast Fourier transform on the M=4096 auxiliary function values ​​to update the conformal mapping coefficients.

[0306] The auxiliary function is:

[0307] ;

[0308] The conformal mapping coefficient is:

[0309] ;

[0310] In this embodiment, since the wellbore boundary has local collapse characteristics, in addition to the 0th principal scale coefficient, there will also be obvious even-order higher-order mapping coefficients, which are used to characterize the local enlargement and high curvature boundary characteristics.

[0311] The obtained conformal mapping coefficients for the first few orders are: m0 = 1.1871231331 + 0.0000000000i, m1 ≈ 0, m2 = 0.2652220407 - 0.0000000000i, m3 ≈ 0, m4 = 0.0679830845 - 0.0000000000i, m5 ≈ 0, ..., m 25 ≈0,m 26 =-0.0000747019+0.0000000000i;

[0312] (1-6) Determine whether the conformal mapping coefficients have converged iteratively; if they have converged, continue to step (1-7); otherwise, return to step (1-2) and iterate again.

[0313] Calculate the maximum change in the mapping coefficients between two consecutive iterations:

[0314] ;

[0315] If Δm (j) If Δm ≤ ε, then the iteration converges. (j) Let be the coefficient convergence error of the j-th iteration, and ε be the preset convergence tolerance, which is taken as... .

[0316] (1-7) Based on the convergence of the iteration in step (1-6), determine whether the adaptive truncation order strategy is satisfied; if satisfied, output the optimal conformal mapping coefficient; otherwise, return to step (1-3).

[0317] The adaptive truncation order strategy is as follows: <0.5% and <0.05;

[0318] The relative boundary fitting error of the nth discrete sampling point on the unit circle boundary is calculated using the following formula:

[0319] .

[0320] The tail proportion factor of higher-order coefficients under the current candidate truncation order is calculated using the following formula:

[0321] .

[0322] In this embodiment, when the truncation order is N=26, the boundary fitting error meets the preset requirements, and the higher-order mapping coefficients have significantly decreased. The optimal conformal mapping coefficients are: m0=1.1871231331+0.0000000000i, m1≈0, m2=0.2652220407-0.0000000000i, m3≈0, m4=0.0679830845-0.0000000000i, m5≈0, …, m 25 ≈0,m 26 =-0.0000747019+0.0000000000i.

[0323] (2) The optimal conformal mapping coefficients obtained above are used for the construction of the perturbation complex potential function and the calculation of the well perimeter stress components. The specific process is the same as in Example 1, and will not be repeated here.

[0324] To illustrate this chain-reaction differentiation process, four representative wellbore points with θ = 0°, 30°, 60°, and 90° are selected for illustration. Anisotropic formations are used as an example:

[0325] μ1 =-0.0000000000+1.3769709346i, μ2 =0.0000000000+0.2096450458i;

[0326] The generalized complex variable representing the discrete sampling points is:

[0327] .

[0328] Further differentiation yields:

[0329] .

[0330] The derivative of the complex potential function can be obtained through , Divide by the table respectively , The remaining discrete sampling points are then calculated point by point in the same manner.

[0331] In this embodiment, there is an effective liquid column pressure inside the wellbore. By performing effective stress correction on the wellbore wall stress, the wellbore stress distribution in isotropic and anisotropic formations can be obtained.

[0332] For irregular wellbore boundaries with high curvature notches and collapse-expanded diameter characteristics, the method of this invention can still stably obtain high-precision boundary representation through adaptive conformal mapping, with an optimal truncation order of 26. The stress calculation results are highly consistent with the finite element analysis. Error statistical analysis shows that under complex explicit collapse shapes, the global average relative error for isotropic methods is 9.31%, as shown below. Figure 5 The global average relative error for anisotropic conditions is 4.47%, such as... Figure 6 .

[0333] Example 4: Real Calculation of Irregular Low-Frequency Rough Borehole

[0334] This embodiment is used to simulate the irregular, asymmetrical, and rough wellbore morphology detected by dual-diameter logging in real engineering projects. The wellbore morphology function is set as an irregular profile with multiple superimposed frequencies: P w =0 MPa.

[0335] Its boundary discrete expression is:

[0336] ;

[0337] ;

[0338] (1) Perform conformal mapping and adaptive truncation order calculation on the well wall boundary to obtain the optimal conformal mapping coefficient;

[0339] (1-1) Using conformal mapping function, construct discrete sampling points of unit circle boundary based on well wall boundary.

[0340] Divide the circumference of the mathematical plane, 2π, into M=4096 equal parts to construct discrete sampling points for subsequent solution of conformal mapping coefficients.

[0341] At each discrete sampling angle θ n Above, calculate the actual boundary points of the wellbore according to the above wellbore boundary equation:

[0342] ;

[0343] ;

[0344] The corresponding complex coordinates are:

[0345] ;

[0346] (1-2) Initialize the conformal mapping coefficients and set the initial truncation order N0=2.

[0347] The initial truncation order N0=2, and the conformal mapping function adopts the Laurent series form.

[0348] ;

[0349] First, initialize the zeroth-order mapping coefficients m0 based on the average polar radius of the true boundary points. (1) The remaining higher-order coefficients are initialized to 0. For the aforementioned boundary, the radius varies along different directions, and the zeroth-order term alone cannot accurately describe this boundary. Therefore, the higher-order mapping coefficients need to be updated through subsequent radial projection iterations and FFT.

[0350] (1-3) Substitute the discrete sampling points of the unit circle boundary into the conformal mapping function to obtain the corresponding current mapped boundary points;

[0351] Discrete sampling points are taken at the boundary of the unit circle: Substituting into the current conformal mapping function, we obtain the current mapping boundary point in the j-th iteration:

[0352] ;

[0353] The current mapped boundary points usually cannot completely coincide with the real boundary points in one go, so radial projection correction is still required.

[0354] (1-4) Radial projection of the current mapped boundary point matches the target boundary point with the same argument.

[0355] For the current mapped boundary point z n (j) Find its argument:

[0356] ;

[0357] And find the target boundary points with the same argument on the well wall boundary:

[0358] ;

[0359] Regarding the above boundary:

[0360] When the argument is When the polar radius of the target boundary point in this direction is:

[0361] ;

[0362] The corresponding target boundary points are:

[0363] ;

[0364] (1-5) Construct an auxiliary function and obtain the conformal mapping coefficients for the new round of iterations through fast Fourier transform;

[0365] Construct auxiliary function values ​​based on the target boundary points. And perform a fast Fourier transform on the M=4096 auxiliary function values ​​to update the conformal mapping coefficients.

[0366] The auxiliary function is:

[0367] ;

[0368] The conformal mapping coefficient is:

[0369] ;

[0370] In this embodiment, due to the roughness of the wellbore boundary, in addition to the 0th order principal scale coefficient, the third and fourth order terms and several higher order terms generated by their interaction have relatively significant contributions.

[0371] The obtained conformal mapping coefficients for the first few orders are: m0=1.0071517427-0.0000000265i, m1=-0.0045014835+0.0053460767i, m2=0.0000012857+0.0000576648i, m3=0.0517969744-0.0282957906i, m4=-0.0075769439+0.0373761491i, m5=0.0005079195+0.0006666806i, ..., m 27 =0.0000050796-0.0000024458i,m 28 =-0.0000023657+0.0000051742i.

[0372] (1-6) Determine whether the conformal mapping coefficients have converged iteratively; if they have converged, continue to step (1-7); otherwise, return to step (1-2) and iterate again.

[0373] Calculate the maximum change in the mapping coefficients between two consecutive iterations:

[0374] ;

[0375] If Δm(j) If ≤ ε, then the iteration converges; Δm (j) Let be the coefficient convergence error of the j-th iteration, and ε be the preset convergence tolerance, which is taken as... .

[0376] (1-7) Based on the convergence of the iteration in step (1-6), determine whether the adaptive truncation order strategy is satisfied; if satisfied, output the optimal conformal mapping coefficient; otherwise, return to step (1-3).

[0377] The adaptive truncation order strategy is as follows: <0.5% and <0.05;

[0378] The relative boundary fitting error of the nth discrete sampling point on the unit circle boundary is calculated using the following formula:

[0379] ;

[0380] The tail proportion factor of higher-order coefficients under the current candidate truncation order is calculated using the following formula:

[0381] ;

[0382] In this embodiment, when the truncation order is N=28, the boundary fitting error meets the preset requirements, and the higher-order mapping coefficients have significantly decreased. The optimal conformal mapping coefficients are: m0=1.0071517427-0.0000000265i, m1=-0.0045014835+0.0053460767i, m2=0.0000012857+0.0000576648i, m3= 0.0517969744-0.0282957906i, m4=-0.0075769439+0.0373761491i, m5=0.0005079195+0.0006666806i, ..., m 27 =0.0000050796-0.0000024458i,m 28 =-0.0000023657+0.0000051742i.

[0383] (2) The optimal conformal mapping coefficients obtained above are used for the construction of the perturbation complex potential function and the calculation of the well perimeter stress components. The specific process is the same as in Example 1, and will not be repeated here.

[0384] To illustrate this chain-reaction differentiation process, four representative wellbore points with θ = 0°, 30°, 60°, and 90° are selected for illustration. Anisotropic formations are used as an example:

[0385] μ1 =-0.0000000000+1.3769709346i, μ2 =0.0000000000+0.2096450458i;

[0386] The generalized complex variable representing the discrete sampling points is:

[0387] .

[0388] Further differentiation yields:

[0389] .

[0390] The derivative of the complex potential function can be obtained through , Divide by the table respectively , The remaining discrete sampling points are then calculated point by point in the same manner.

[0391] In this embodiment, there is no effective fluid column pressure inside the wellbore, so there is no need for wellbore fluid column pressure correction. The wellbore stress distribution of irregular low-frequency rough wellbore in isotropic and anisotropic formations can be directly output as follows: Figure 7 and Figure 8 As shown.

Claims

1. A method for analytical calculation of wellbore stress in irregular wellbore shapes, characterized in that, Based on radial projection iteration and fast Fourier transform, and combining the maximum change of conformal mapping coefficients in two adjacent iterations, relative boundary fitting error, and tail proportion factor of higher-order coefficients, the cutoff order of the conformal mapping function and the optimal conformal mapping coefficients are adaptively determined, and the optimal conformal mapping coefficients are used for the construction of perturbation complex potential function and the calculation of wellbore stress components.

2. The analytical calculation method for wellbore stress of irregular wellbore shape according to claim 1, characterized in that, The specific steps to determine the optimal conformal mapping coefficients are as follows: (1-1) Using conformal mapping functions, construct discrete sampling points of the unit circle boundary based on the well wall boundary; (1-2) Set the initial truncation order and initialize the conformal mapping coefficients; (1-3) Substitute the discrete sampling points of the unit circle boundary into the conformal mapping function to obtain the corresponding current mapped boundary points; (1-4) Radial projection of the current mapped boundary point matches the target boundary point with the same argument; (1-5) Construct an auxiliary function and obtain the conformal mapping coefficients for the new round of iterations through fast Fourier transform; (1-6) Determine whether the conformal mapping coefficients have converged in iteration by combining the maximum change in the conformal mapping coefficients in two adjacent iterations; if they have converged, continue to step (1-7); otherwise, return to step (1-2) and iterate again. (1-7) Based on the convergence of the iteration in step (1-6), determine whether the adaptive truncation order strategy is satisfied according to the relative boundary fitting error and the tail proportion factor of the higher order coefficients; if satisfied, output the optimal conformal mapping coefficients. Otherwise, return to steps (1-3).

3. The analytical calculation method for wellbore stress of irregular wellbore shape according to claim 2, characterized in that, In steps (1-5), the auxiliary function is: ； The conformal mapping coefficient is: ； in, Let m be the auxiliary function value corresponding to the nth discrete sampling point in the j-th iteration; Let m be the conformal mapping coefficient of the kth order after the (j+1)th iteration update; θ is the nth discrete sampling point on the boundary of the unit circle, dimensionless; n The nth discrete sampling angle of the unit circle boundary, in rad; The target boundary point is at the same argument as the current mapped boundary point, m; k is the order number of the conformal mapping coefficient, dimensionless; M is the total number of discrete sampling points on the unit circle boundary, dimensionless; n is the number of the nth discrete sampling point on the unit circle boundary, dimensionless.

4. The analytical calculation method for wellbore stress of irregular wellbore shape according to claim 1, characterized in that, The maximum change in conformal mapping coefficients between two consecutive iterations is: ； in, Let m be the conformal mapping coefficient of the kth order after the (j+1)th iteration update; Let m be the conformal mapping coefficient of the kth order after the j-th iteration update; Let m be the convergence error of the coefficients in the j-th iteration; To preset the convergence tolerance, we take 10. -10 ;like ≤ If the iteration converges, then the iteration will be successful.

5. The analytical calculation method for wellbore stress of irregular wellbore shape according to claim 1, characterized in that, The relative boundary fitting error is: ； in, Let m be the polar radius of the physical plane obtained after conformal mapping of the nth discrete sampling point when the truncation order is N; Let m be the polar radius of the target boundary point with the same argument when the truncation order is N; The relative boundary fitting error of the nth discrete sampling point is given by the truncation order N, and is dimensionless; N is the truncation order, which is also dimensionless.

6. The analytical calculation method for wellbore stress of irregular wellbore shape according to claim 1, characterized in that, The tail proportion factor of higher-order coefficients is: ； in, is the tail proportion factor of higher-order coefficients when the truncation order is N, dimensionless; L is the statistical length of the tail coefficients, dimensionless; ks is the starting order of the tail statistics, dimensionless. Let m be the conformal mapping coefficient of the kth order when the truncation order is N.

7. The analytical calculation method for wellbore stress of irregular wellbore shape according to claim 2, characterized in that, In steps (1-7), the adaptive truncation order strategy is as follows: <0.5% and <0.

05.

8. The analytical calculation method for wellbore stress of irregular wellbore shape according to claim 1, characterized in that, The perturbation complex potential function is: ； ； in, Let MPa be the complex potential function of the first disturbance. A1, A2, A3, A4 are the coefficients of the first perturbation complex potential function, MPa; A5, A6, A7, A8 are the coefficients of the second perturbation complex potential function, MPa. Let m be the real part of the optimal conformal mapping coefficient; m is the imaginary part of the optimal conformal mapping coefficients; m is the real part of the 0th-order conformal mapping coefficient; m is the imaginary part of the 0th-order conformal mapping coefficient; It is a mathematical plane, dimensionless.

9. The analytical calculation method for wellbore stress of irregular wellbore shape according to claim 1, characterized in that, The formula for calculating the wellbore stress components is: ； ； ； Wherein, the direction angle δ is the directed angle between the unit outward normal vector and the positive x-axis of the global rectangular coordinate system, in rad; The radial stress in the wellbore coordinate system is expressed in MPa. Let be the circumferential stress in the wellbore coordinate system, in MPa; The stress is the tangential shear stress in the wellbore coordinate system, in MPa; The normal stress in the x-direction at the wellbore boundary is expressed in MPa. The normal stress in the y-direction at the wellbore boundary is given in MPa. The shear stress at the wellbore boundary is expressed in MPa.

Images